Bug Summary

File:src/pull/binningPull.c
Location:line 460, column 9
Description:Value stored to 'point' is never read

Annotated Source Code

1/*
2 Teem: Tools to process and visualize scientific data and images .
3 Copyright (C) 2013, 2012, 2011, 2010, 2009 University of Chicago
4 Copyright (C) 2008, 2007, 2006, 2005 Gordon Kindlmann
5 Copyright (C) 2004, 2003, 2002, 2001, 2000, 1999, 1998 University of Utah
6
7 This library is free software; you can redistribute it and/or
8 modify it under the terms of the GNU Lesser General Public License
9 (LGPL) as published by the Free Software Foundation; either
10 version 2.1 of the License, or (at your option) any later version.
11 The terms of redistributing and/or modifying this software also
12 include exceptions to the LGPL that facilitate static linking.
13
14 This library is distributed in the hope that it will be useful,
15 but WITHOUT ANY WARRANTY; without even the implied warranty of
16 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
17 Lesser General Public License for more details.
18
19 You should have received a copy of the GNU Lesser General Public License
20 along with this library; if not, write to Free Software Foundation, Inc.,
21 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
22*/
23
24#include "pull.h"
25#include "privatePull.h"
26
27/*
28** because the pullContext keeps an array of bins (not pointers to them)
29** we have Init and Done functions (not New and Nix)
30*/
31void
32_pullBinInit(pullBin *bin) {
33
34 bin->point = NULL((void*)0);
35 bin->pointNum = 0;
36 bin->pointArr = NULL((void*)0);
37 bin->neighBin = NULL((void*)0);
38 return;
39}
40
41/*
42** bins own the points they contain- so this frees them
43*/
44void
45_pullBinDone(pullBin *bin) {
46 unsigned int idx;
47
48 for (idx=0; idx<bin->pointNum; idx++) {
49 bin->point[idx] = pullPointNix(bin->point[idx]);
50 }
51 bin->pointArr = airArrayNuke(bin->pointArr);
52 bin->neighBin = (pullBin **)airFree(bin->neighBin);
53 return;
54}
55
56int
57_pullBinNeighborSet(pullContext *pctx, pullBin *bin) {
58 static const char me[]="_pullBinNeighborSet";
59 unsigned int neiIdx, neiNum, be[4], binIdx;
60 unsigned int xi, yi, zi, si, xx, yy, zz, ss,
61 xmax, ymax, zmax, smax;
62 pullBin *nei[3*3*3*3];
63 int xmin, ymin, zmin, smin;
64
65 binIdx = AIR_UINT(bin - pctx->bin)((unsigned int)(bin - pctx->bin));
66 /* annoyingly, have to recover the bin coordinates */
67 ELL_4V_COPY(be, pctx->binsEdge)((be)[0] = (pctx->binsEdge)[0], (be)[1] = (pctx->binsEdge
)[1], (be)[2] = (pctx->binsEdge)[2], (be)[3] = (pctx->binsEdge
)[3]);
68 xi = binIdx % be[0];
69 binIdx = (binIdx - xi)/be[0];
70 yi = binIdx % be[1];
71 binIdx = (binIdx - yi)/be[1];
72 zi = binIdx % be[2];
73 si = (binIdx - zi)/be[2];
74 neiNum = 0;
75 bin->neighBin = (pullBin **)airFree(bin->neighBin);
76 smin = AIR_MAX(0, (int)si-1)((0) > ((int)si-1) ? (0) : ((int)si-1));
77 smax = AIR_MIN(si+1, be[3]-1)((si+1) < (be[3]-1) ? (si+1) : (be[3]-1));
78 zmin = AIR_MAX(0, (int)zi-1)((0) > ((int)zi-1) ? (0) : ((int)zi-1));
79 zmax = AIR_MIN(zi+1, be[2]-1)((zi+1) < (be[2]-1) ? (zi+1) : (be[2]-1));
80 ymin = AIR_MAX(0, (int)yi-1)((0) > ((int)yi-1) ? (0) : ((int)yi-1));
81 ymax = AIR_MIN(yi+1, be[1]-1)((yi+1) < (be[1]-1) ? (yi+1) : (be[1]-1));
82 xmin = AIR_MAX(0, (int)xi-1)((0) > ((int)xi-1) ? (0) : ((int)xi-1));
83 xmax = AIR_MIN(xi+1, be[0]-1)((xi+1) < (be[0]-1) ? (xi+1) : (be[0]-1));
84 for (ss=smin; ss<=smax; ss++) {
85 for (zz=zmin; zz<=zmax; zz++) {
86 for (yy=ymin; yy<=ymax; yy++) {
87 for (xx=xmin; xx<=xmax; xx++) {
88 nei[neiNum++] = pctx->bin + xx + be[0]*(yy + be[1]*(zz + be[2]*ss));
89 }
90 }
91 }
92 }
93 if (!( bin->neighBin = AIR_CALLOC(1+neiNum, pullBin*)(pullBin**)(calloc((1+neiNum), sizeof(pullBin*))) )) {
94 biffAddf(PULLpullBiffKey, "%s: couldn't calloc array of %u neighbor pointers", me, 1+neiNum);
95 return 1;
96 }
97 for (neiIdx=0; neiIdx<neiNum; neiIdx++) {
98 bin->neighBin[neiIdx] = nei[neiIdx];
99 }
100 /* NULL-terminate the bin array */
101 bin->neighBin[neiIdx] = NULL((void*)0);
102 return 0;
103}
104
105/*
106** bins on boundary now extend to infinity; so the only time this
107** returns NULL (indicating error) is for non-existent positions
108*/
109pullBin *
110_pullBinLocate(pullContext *pctx, double *posWorld) {
111 static const char me[]="_pullBinLocate";
112 unsigned int axi, eidx[4], binIdx;
113
114 if (!ELL_4V_EXISTS(posWorld)((((int)(!(((posWorld)[0]) - ((posWorld)[0]))))) && (
((int)(!(((posWorld)[1]) - ((posWorld)[1]))))) && (((
int)(!(((posWorld)[2]) - ((posWorld)[2]))))) && (((int
)(!(((posWorld)[3]) - ((posWorld)[3]))))))) {
115 biffAddf(PULLpullBiffKey, "%s: non-existent position (%g,%g,%g,%g)", me,
116 posWorld[0], posWorld[1], posWorld[2], posWorld[3]);
117 return NULL((void*)0);
118 }
119
120 if (pctx->flag.binSingle) {
121 binIdx = 0;
122 } else {
123 for (axi=0; axi<4; axi++) {
124 eidx[axi] = airIndexClamp(pctx->bboxMin[axi],
125 posWorld[axi],
126 pctx->bboxMax[axi],
127 pctx->binsEdge[axi]);
128 }
129 binIdx = (eidx[0]
130 + pctx->binsEdge[0]*(
131 eidx[1] + pctx->binsEdge[1]*(
132 eidx[2] + pctx->binsEdge[2] * eidx[3])));
133 }
134
135 return pctx->bin + binIdx;
136}
137
138/*
139** this makes the bin the owner of the point
140*/
141int
142_pullBinPointAdd(pullContext *pctx, pullBin *bin, pullPoint *point) {
143 static const char me[]="_pullBinPointAdd";
144 int pntI;
145 pullPtrPtrUnion pppu;
146
147 AIR_UNUSED(pctx)(void)(pctx);
148 if (!(bin->pointArr)) {
149 pppu.points = &(bin->point);
150 bin->pointArr = airArrayNew(pppu.v, &(bin->pointNum),
151 sizeof(pullPoint *), _PULL_BIN_INCR32);
152 if (!( bin->pointArr )) {
153 biffAddf(PULLpullBiffKey, "%s: couldn't create point array", me);
154 return 1;
155 }
156 }
157 if (!( bin->neighBin )) {
158 /* set up neighbor bin vector if not done so already */
159 if (_pullBinNeighborSet(pctx, bin)) {
160 biffAddf(PULLpullBiffKey, "%s: couldn't initialize neighbor bins", me);
161 return 1;
162 }
163 }
164 pntI = airArrayLenIncr(bin->pointArr, 1);
165 bin->point[pntI] = point;
166 return 0;
167}
168
169/*
170** the bin loses track of the point, caller responsible for ownership,
171** even though caller never identifies it by pointer, which is weird
172*/
173void
174_pullBinPointRemove(pullContext *pctx, pullBin *bin, int loseIdx) {
175
176 AIR_UNUSED(pctx)(void)(pctx);
177 bin->point[loseIdx] = bin->point[bin->pointNum-1];
178 airArrayLenIncr(bin->pointArr, -1);
179 return;
180}
181
182/*
183** adds point to context
184*/
185int
186pullBinsPointAdd(pullContext *pctx, pullPoint *point, pullBin **binP) {
187 static const char me[]="pullBinsPointAdd";
188 pullBin *bin;
189
190 if (binP) {
191 *binP = NULL((void*)0);
192 }
193 if (!( bin = _pullBinLocate(pctx, point->pos) )) {
194 biffAddf(PULLpullBiffKey, "%s: can't locate point %p %u",
195 me, AIR_CAST(void*, point)((void*)(point)), point->idtag);
196 return 1;
197 }
198 if (binP) {
199 *binP = bin;
200 }
201 if (_pullBinPointAdd(pctx, bin, point)) {
202 biffAddf(PULLpullBiffKey, "%s: trouble adding point %p %u",
203 me, AIR_CAST(void*, point)((void*)(point)), point->idtag);
204 return 1;
205 }
206 return 0;
207}
208
209int
210pullBinsPointMaybeAdd(pullContext *pctx, pullPoint *point,
211 /* output */
212 pullBin **binP, int *added) {
213 static const char me[]="pullBinsPointMaybeAdd";
214 pullBin *bin;
215 unsigned int idx;
216 int okay;
217
218 if (binP) {
219 *binP = NULL((void*)0);
220 }
221 if (!(pctx && point && added)) {
222 biffAddf(PULLpullBiffKey, "%s: got NULL pointer", me);
223 return 1;
224 }
225 if (!( bin = _pullBinLocate(pctx, point->pos) )) {
226 biffAddf(PULLpullBiffKey, "%s: can't locate point %p %u",
227 me, AIR_CAST(void*, point)((void*)(point)), point->idtag);
228 return 1;
229 }
230 if (binP) {
231 *binP = bin;
232 }
233 if (pctx->flag.restrictiveAddToBins) {
234 okay = AIR_TRUE1;
235 for (idx=0; idx<bin->pointNum; idx++) {
236 double diff[4], len;
237 ELL_4V_SUB(diff, point->pos, bin->point[idx]->pos)((diff)[0] = (point->pos)[0] - (bin->point[idx]->pos
)[0], (diff)[1] = (point->pos)[1] - (bin->point[idx]->
pos)[1], (diff)[2] = (point->pos)[2] - (bin->point[idx]
->pos)[2], (diff)[3] = (point->pos)[3] - (bin->point
[idx]->pos)[3]);
238 ELL_3V_SCALE(diff, 1/pctx->sysParm.radiusSpace, diff)((diff)[0] = (1/pctx->sysParm.radiusSpace)*(diff)[0], (diff
)[1] = (1/pctx->sysParm.radiusSpace)*(diff)[1], (diff)[2] =
(1/pctx->sysParm.radiusSpace)*(diff)[2]);
239 diff[3] /= pctx->sysParm.radiusScale;
240 len = ELL_4V_LEN(diff)(sqrt((((diff))[0]*((diff))[0] + ((diff))[1]*((diff))[1] + ((
diff))[2]*((diff))[2] + ((diff))[3]*((diff))[3])));
241 if (len < _PULL_BINNING_MAYBE_ADD_THRESH0.1) {
242 okay = AIR_FALSE0;
243 break;
244 }
245 }
246 if (okay) {
247 if (_pullBinPointAdd(pctx, bin, point)) {
248 biffAddf(PULLpullBiffKey, "%s: trouble adding point %p %u",
249 me, AIR_CAST(void*, point)((void*)(point)), point->idtag);
250 return 1;
251 }
252 *added = AIR_TRUE1;
253 } else {
254 *added = AIR_FALSE0;
255 }
256 } else {
257 if (_pullBinPointAdd(pctx, bin, point)) {
258 biffAddf(PULLpullBiffKey, "%s: trouble adding point %p %u",
259 me, AIR_CAST(void*, point)((void*)(point)), point->idtag);
260 return 1;
261 }
262 *added = AIR_TRUE1;
263 }
264 return 0;
265}
266
267int
268_pullBinSetup(pullContext *pctx) {
269 static const char me[]="_pullBinSetup";
270 unsigned ii;
271 double volEdge[4], width;
272
273 /* the maximum distance of interaction is when one particle is sitting
274 on the edge of another particle's sphere of influence, *NOT*, when
275 the spheres of influence of two particles are tangent: particles
276 interact with potential fields of other particles, but there is
277 no interaction between potential fields. */
278 width = (pctx->sysParm.radiusSpace ? pctx->sysParm.radiusSpace : 0.1);
279 pctx->maxDistSpace = pctx->sysParm.binWidthSpace*width;
280 width = (pctx->sysParm.radiusScale ? pctx->sysParm.radiusScale : 0.1);
281 pctx->maxDistScale = 1*width;
282
283 if (pctx->verbose) {
284 printf("%s: radiusSpace = %g -(%g)-> maxDistSpace = %g\n", me,
285 pctx->sysParm.radiusSpace, pctx->sysParm.binWidthSpace,
286 pctx->maxDistSpace);
287 printf("%s: radiusScale = %g ----> maxDistScale = %g\n", me,
288 pctx->sysParm.radiusScale, pctx->maxDistScale);
289 }
290
291 if (pctx->flag.binSingle) {
292 pctx->binsEdge[0] = 1;
293 pctx->binsEdge[1] = 1;
294 pctx->binsEdge[2] = 1;
295 pctx->binsEdge[3] = 1;
296 pctx->binNum = 1;
297 } else {
298 volEdge[0] = pctx->bboxMax[0] - pctx->bboxMin[0];
299 volEdge[1] = pctx->bboxMax[1] - pctx->bboxMin[1];
300 volEdge[2] = pctx->bboxMax[2] - pctx->bboxMin[2];
301 volEdge[3] = pctx->bboxMax[3] - pctx->bboxMin[3];
302 if (pctx->verbose) {
303 printf("%s: volEdge = %g %g %g %g\n", me,
304 volEdge[0], volEdge[1], volEdge[2], volEdge[3]);
305 }
306 pctx->binsEdge[0] = AIR_CAST(unsigned int,((unsigned int)(floor(volEdge[0]/pctx->maxDistSpace)))
307 floor(volEdge[0]/pctx->maxDistSpace))((unsigned int)(floor(volEdge[0]/pctx->maxDistSpace)));
308 pctx->binsEdge[0] = pctx->binsEdge[0] ? pctx->binsEdge[0] : 1;
309 pctx->binsEdge[1] = AIR_CAST(unsigned int,((unsigned int)(floor(volEdge[1]/pctx->maxDistSpace)))
310 floor(volEdge[1]/pctx->maxDistSpace))((unsigned int)(floor(volEdge[1]/pctx->maxDistSpace)));
311 pctx->binsEdge[1] = pctx->binsEdge[1] ? pctx->binsEdge[1] : 1;
312 pctx->binsEdge[2] = AIR_CAST(unsigned int,((unsigned int)(floor(volEdge[2]/pctx->maxDistSpace)))
313 floor(volEdge[2]/pctx->maxDistSpace))((unsigned int)(floor(volEdge[2]/pctx->maxDistSpace)));
314 pctx->binsEdge[2] = pctx->binsEdge[2] ? pctx->binsEdge[2] : 1;
315 pctx->binsEdge[3] = AIR_CAST(unsigned int,((unsigned int)(floor(volEdge[3]/pctx->maxDistScale)))
316 floor(volEdge[3]/pctx->maxDistScale))((unsigned int)(floor(volEdge[3]/pctx->maxDistScale)));
317 pctx->binsEdge[3] = pctx->binsEdge[3] ? pctx->binsEdge[3] : 1;
318 /* hack to observe things at bin boundaries
319 ELL_3V_SET(pctx->binsEdge, 3, 3, 3);
320 */
321 if (pctx->verbose) {
322 printf("%s: binsEdge=(%u,%u,%u,%u)\n", me,
323 pctx->binsEdge[0], pctx->binsEdge[1],
324 pctx->binsEdge[2], pctx->binsEdge[3]);
325 }
326 pctx->binNum = (pctx->binsEdge[0]*pctx->binsEdge[1]
327 *pctx->binsEdge[2]*pctx->binsEdge[3]);
328 }
329 if (pctx->binNum > PULL_BIN_MAXNUM40000000) {
330 biffAddf(PULLpullBiffKey, "%s: sorry, #bins %u > PULL_BIN_MAXNUM %u. Try "
331 "increasing pctx->sysParm.binWidthSpace (%g)", me,
332 pctx->binNum, PULL_BIN_MAXNUM40000000, pctx->sysParm.binWidthSpace);
333 return 1;
334 }
335 if (pctx->verbose) {
336 printf("%s: trying to allocate %u bins . . . \n", me, pctx->binNum);
337 }
338 pctx->bin = AIR_CALLOC(pctx->binNum, pullBin)(pullBin*)(calloc((pctx->binNum), sizeof(pullBin)));
339 if (!( pctx->bin )) {
340 biffAddf(PULLpullBiffKey, "%s: couln't allocate %u bins", me, pctx->binNum);
341 return 1;
342 }
343 if (pctx->verbose) {
344 printf("%s: done allocating. Initializing . . . \n", me);
345 }
346 for (ii=0; ii<pctx->binNum; ii++) {
347 _pullBinInit(pctx->bin + ii);
348 }
349 if (pctx->verbose) {
350 printf("%s: done initializing.\n", me);
351 }
352 if (pctx->flag.binSingle) {
353 if (!( pctx->bin[0].neighBin = AIR_CALLOC(2, pullBin*)(pullBin**)(calloc((2), sizeof(pullBin*))) )) {
354 biffAddf(PULLpullBiffKey, "%s: trouble allocating for single bin?", me);
355 return 1;
356 }
357 pctx->bin[0].neighBin[0] = pctx->bin + 0;
358 pctx->bin[0].neighBin[1] = NULL((void*)0);
359 }
360 return 0;
361}
362
363void
364_pullBinFinish(pullContext *pctx) {
365 unsigned int ii;
366
367 for (ii=0; ii<pctx->binNum; ii++) {
368 _pullBinDone(pctx->bin + ii);
369 }
370 pctx->bin = (pullBin *)airFree(pctx->bin);
371 ELL_4V_SET(pctx->binsEdge, 0, 0, 0, 0)((pctx->binsEdge)[0] = (0), (pctx->binsEdge)[1] = (0), (
pctx->binsEdge)[2] = (0), (pctx->binsEdge)[3] = (0));
372 pctx->binNum = 0;
373}
374
375/*
376** sets pctx->stuckNum
377** resets all task[]->stuckNum
378** reallocates pctx->tmpPointPerm and pctx->tmpPointPtr
379** the point of this is to do rebinning
380**
381** This function is only called by the master thread, this
382** does *not* have to be thread-safe in any way
383*/
384int
385_pullIterFinishDescent(pullContext *pctx) {
386 static const char me[]="_pullIterFinishDescent";
387 unsigned int oldBinIdx, pointIdx, taskIdx, runIdx, pointNum;
388 pullBin *oldBin;
389 pullPoint *point;
390 int added;
391
392 _pullNixTheNixed(pctx);
393
394 pctx->stuckNum = 0;
395 for (taskIdx=0; taskIdx<pctx->threadNum; taskIdx++) {
396 pctx->stuckNum += pctx->task[taskIdx]->stuckNum;
397 pctx->task[taskIdx]->stuckNum = 0;
398 }
399
400 /* even w/ a single bin, we may still have to permute the points */
401 pointNum = pullPointNumber(pctx);
402 if (pointNum != pctx->tmpPointNum) {
403 if (pctx->verbose) {
404 printf("!%s: changing total point # %u --> %u\n", me,
405 pctx->tmpPointNum, pointNum);
406 }
407 airFree(pctx->tmpPointPerm);
408 airFree(pctx->tmpPointPtr);
409 pctx->tmpPointPtr = AIR_CAST(pullPoint **,((pullPoint **)(calloc(pointNum, sizeof(pullPoint*))))
410 calloc(pointNum, sizeof(pullPoint*)))((pullPoint **)(calloc(pointNum, sizeof(pullPoint*))));
411 pctx->tmpPointPerm = AIR_CAST(unsigned int *,((unsigned int *)(calloc(pointNum, sizeof(unsigned int))))
412 calloc(pointNum, sizeof(unsigned int)))((unsigned int *)(calloc(pointNum, sizeof(unsigned int))));
413 if (!( pctx->tmpPointPtr && pctx->tmpPointPerm )) {
414 biffAddf(PULLpullBiffKey, "%s: couldn't allocate tmp buffers %p %p", me,
415 AIR_VOIDP(pctx->tmpPointPtr)((void *)(pctx->tmpPointPtr)), AIR_VOIDP(pctx->tmpPointPerm)((void *)(pctx->tmpPointPerm)));
416 return 1;
417 }
418 pctx->tmpPointNum = pointNum;
419 }
420 runIdx = 0;
421 for (oldBinIdx=0; oldBinIdx<pctx->binNum; oldBinIdx++) {
422 oldBin = pctx->bin + oldBinIdx;
423 while (oldBin->pointNum) {
424 /* tricky: we can't traverse bin->point[], because of how it is
425 re-ordered on point removal, so we always grab point[0] */
426 pctx->tmpPointPtr[runIdx++] = oldBin->point[0];
427 _pullBinPointRemove(pctx, oldBin, 0);
428 }
429 }
430 airShuffle_r(pctx->task[0]->rng,
431 pctx->tmpPointPerm, pointNum, pctx->flag.permuteOnRebin);
432 if (pctx->flag.permuteOnRebin && pctx->verbose) {
433 printf("%s: permuting %u points\n", me, pointNum);
434 }
435 for (pointIdx=0; pointIdx<pointNum; pointIdx++) {
436 point = pctx->tmpPointPtr[pctx->tmpPointPerm[pointIdx]];
437 /*
438 if (1 == pctx->iter && 102 == point->idtag) {
439 _pullDebugSanity(pctx->task[0], point);
440 }
441 */
442 /* Sun Jul 14 01:30:41 CDT 2013: the previous code was basically
443 just pullBinsPointAdd(). But when working with codim-3
444 constraints (like finding spatial maxima with maximal strength
445 along scale), its easy for many points to start piling on top
446 of each other; that is the problem that
447 pullFlagRestrictiveAddToBins was designed to solve. */
448 if (pctx->constraint && 0 == pctx->constraintDim) {
449 if (pullBinsPointMaybeAdd(pctx, point, NULL((void*)0), &added)) {
450 biffAddf(PULLpullBiffKey, "%s: trouble binning? point %u", me, point->idtag);
451 return 1;
452 }
453 if (!added) {
454 /* the point wasn't owned by any bin, and now it turns out
455 no bin wanted to own it, so we have to remove it */
456 /* in the case of point maxima searching; this mainly happened
457 because two points at the same spatial positition slowly
458 descended along scale towards each other at the scale of
459 maximal strength */
460 point = pullPointNix(point);
Value stored to 'point' is never read
461 }
462 } else {
463 if (pullBinsPointAdd(pctx, point, NULL((void*)0))) {
464 biffAddf(PULLpullBiffKey, "%s: trouble binning point %u", me, point->idtag);
465 return 1;
466 }
467 point = NULL((void*)0);
468 }
469 pctx->tmpPointPtr[pctx->tmpPointPerm[pointIdx]] = NULL((void*)0);
470 }
471
472 return 0;
473}
474