Bug Summary

File:src/nrrd/cc.c
Location:line 339, column 3
Description:Value stored to 'ret' 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 "nrrd.h"
25#include "privateNrrd.h"
26
27/*
28** learned: if you have globals, such as _nrrdCC_verb, which are
29** defined and declared here, but which are NOT initialized, then
30** C++ apps which are linking against Teem will have problems!!!
31** This was first seen on the mac.
32*/
33int _nrrdCC_verb = 0;
34int _nrrdCC_EqvIncr = 10000; /* HEY: this has to be big so that ccfind is not
35 stuck constantly re-allocating the eqv array.
36 This is will be one of the best places to
37 test the new multiplicative reallocation
38 strategy, planned for Teem 2.0 */
39
40int
41_nrrdCCFind_1(Nrrd *nout, unsigned int *numid, const Nrrd *nin) {
42 /* static const char me[]="_nrrdCCFind_1"; */
43 unsigned int sx, I, id, lval, val, *out, (*lup)(const void *, size_t);
44
45 lup = nrrdUILookup[nin->type];
46 out = AIR_CAST(unsigned int*, nout->data)((unsigned int*)(nout->data));
47 out[0] = id = 0;
48 *numid = 1;
49
50 sx = AIR_CAST(unsigned int, nin->axis[0].size)((unsigned int)(nin->axis[0].size));
51 lval = lup(nin->data, 0);
52 for (I=1; I<sx; I++) {
53 val = lup(nin->data, I);
54 if (lval != val) {
55 id++;
56 (*numid)++;
57 }
58 out[I] = id;
59 lval = val;
60 }
61
62 return 0;
63}
64
65/*
66** layout of value (pvl) and index (pid) cache:
67**
68** 2 3 4 --> X
69** 1 . . oddly, index 0 is never used
70** . . .
71** |
72** v Y
73*/
74int
75_nrrdCCFind_2(Nrrd *nout, unsigned int *numid, airArray *eqvArr,
76 const Nrrd *nin, unsigned int conny) {
77 static const char me[]="_nrrdCCFind_2";
78 double vl=0, pvl[5]={0,0,0,0,0};
79 unsigned int id, pid[5]={0,0,0,0,0}, (*lup)(const void *, size_t), *out;
80 unsigned int p, x, y, sx, sy;
81
82 id = 0; /* sssh! compiler warnings */
83 lup = nrrdUILookup[nin->type];
84 out = AIR_CAST(unsigned int*, nout->data)((unsigned int*)(nout->data));
85 sx = AIR_CAST(unsigned int, nin->axis[0].size)((unsigned int)(nin->axis[0].size));
86 sy = AIR_CAST(unsigned int, nin->axis[1].size)((unsigned int)(nin->axis[1].size));
87#define GETV_2(x,y)((((0) <= (((int)(x))) && (((int)(x))) <= (((int
)(sx-1)))) && ((0) <= (((int)(y))) && (((int
)(y))) <= (((int)(sy-1))))) ? lup(nin->data, (x) + sx*(
y)) : 0.5) ((AIR_IN_CL(0, AIR_CAST(int, x), AIR_CAST(int, sx-1))((0) <= (((int)(x))) && (((int)(x))) <= (((int)
(sx-1)))) \
88 && AIR_IN_CL(0, AIR_CAST(int, y), AIR_CAST(int, sy-1))((0) <= (((int)(y))) && (((int)(y))) <= (((int)
(sy-1))))) \
89 ? lup(nin->data, (x) + sx*(y)) \
90 : 0.5) /* value that can't come from an array of uints */
91#define GETI_2(x,y)((((0) <= (((int)(x))) && (((int)(x))) <= (((int
)(sx-1)))) && ((0) <= (((int)(y))) && (((int
)(y))) <= (((int)(sy-1))))) ? out[(x) + sx*(y)] : ((unsigned
int)(-1))) ((AIR_IN_CL(0, AIR_CAST(int, x), AIR_CAST(int, sx-1))((0) <= (((int)(x))) && (((int)(x))) <= (((int)
(sx-1)))) \
92 && AIR_IN_CL(0, AIR_CAST(int, y), AIR_CAST(int, sy-1))((0) <= (((int)(y))) && (((int)(y))) <= (((int)
(sy-1))))) \
93 ? out[(x) + sx*(y)] \
94 : AIR_CAST(unsigned int, -1)((unsigned int)(-1))) /* CC index (probably!)
95 never assigned */
96
97 *numid = 0;
98 for (y=0; y<sy; y++) {
99 for (x=0; x<sx; x++) {
100 if (_nrrdCC_verb) {
101 fprintf(stderr__stderrp, "%s(%d,%d) -----------\n", me, x, y);
102 }
103 if (!x) {
104 pvl[1] = GETV_2(-1, y)((((0) <= (((int)(-1))) && (((int)(-1))) <= (((
int)(sx-1)))) && ((0) <= (((int)(y))) && (
((int)(y))) <= (((int)(sy-1))))) ? lup(nin->data, (-1) +
sx*(y)) : 0.5); pid[1] = GETI_2(-1, y)((((0) <= (((int)(-1))) && (((int)(-1))) <= (((
int)(sx-1)))) && ((0) <= (((int)(y))) && (
((int)(y))) <= (((int)(sy-1))))) ? out[(-1) + sx*(y)] : ((
unsigned int)(-1)));
105 pvl[2] = GETV_2(-1, y-1)((((0) <= (((int)(-1))) && (((int)(-1))) <= (((
int)(sx-1)))) && ((0) <= (((int)(y-1))) &&
(((int)(y-1))) <= (((int)(sy-1))))) ? lup(nin->data, (
-1) + sx*(y-1)) : 0.5); pid[2] = GETI_2(-1, y-1)((((0) <= (((int)(-1))) && (((int)(-1))) <= (((
int)(sx-1)))) && ((0) <= (((int)(y-1))) &&
(((int)(y-1))) <= (((int)(sy-1))))) ? out[(-1) + sx*(y-1)
] : ((unsigned int)(-1)));
106 pvl[3] = GETV_2(0, y-1)((((0) <= (((int)(0))) && (((int)(0))) <= (((int
)(sx-1)))) && ((0) <= (((int)(y-1))) && ((
(int)(y-1))) <= (((int)(sy-1))))) ? lup(nin->data, (0) +
sx*(y-1)) : 0.5); pid[3] = GETI_2(0, y-1)((((0) <= (((int)(0))) && (((int)(0))) <= (((int
)(sx-1)))) && ((0) <= (((int)(y-1))) && ((
(int)(y-1))) <= (((int)(sy-1))))) ? out[(0) + sx*(y-1)] : (
(unsigned int)(-1)));
107
108 } else {
109 pvl[1] = vl; pid[1] = id;
110 pvl[2] = pvl[3]; pid[2] = pid[3];
111 pvl[3] = pvl[4]; pid[3] = pid[4];
112 }
113 pvl[4] = GETV_2(x+1, y-1)((((0) <= (((int)(x+1))) && (((int)(x+1))) <= (
((int)(sx-1)))) && ((0) <= (((int)(y-1))) &&
(((int)(y-1))) <= (((int)(sy-1))))) ? lup(nin->data, (
x+1) + sx*(y-1)) : 0.5); pid[4] = GETI_2(x+1, y-1)((((0) <= (((int)(x+1))) && (((int)(x+1))) <= (
((int)(sx-1)))) && ((0) <= (((int)(y-1))) &&
(((int)(y-1))) <= (((int)(sy-1))))) ? out[(x+1) + sx*(y-1
)] : ((unsigned int)(-1)));
114 vl = GETV_2(x, y)((((0) <= (((int)(x))) && (((int)(x))) <= (((int
)(sx-1)))) && ((0) <= (((int)(y))) && (((int
)(y))) <= (((int)(sy-1))))) ? lup(nin->data, (x) + sx*(
y)) : 0.5);
115 p = 0;
116 if (vl == pvl[1]) {
117 id = pid[p=1];
118 }
119#define TEST(P)if (vl == pvl[(P)]) { if (p) { if (id != pid[(P)]) { airEqvAdd
(eqvArr, pid[(P)], id); } } else { id = pid[p=(P)]; } } \
120 if (vl == pvl[(P)]) { \
121 if (p) { /* we already had a value match */ \
122 if (id != pid[(P)]) { \
123 airEqvAdd(eqvArr, pid[(P)], id); \
124 } \
125 } else { \
126 id = pid[p=(P)]; \
127 } \
128 }
129 TEST(3)if (vl == pvl[(3)]) { if (p) { if (id != pid[(3)]) { airEqvAdd
(eqvArr, pid[(3)], id); } } else { id = pid[p=(3)]; } };
130 if (2 == conny) {
131 TEST(2)if (vl == pvl[(2)]) { if (p) { if (id != pid[(2)]) { airEqvAdd
(eqvArr, pid[(2)], id); } } else { id = pid[p=(2)]; } };
132 TEST(4)if (vl == pvl[(4)]) { if (p) { if (id != pid[(4)]) { airEqvAdd
(eqvArr, pid[(4)], id); } } else { id = pid[p=(4)]; } };
133 }
134 if (!p) {
135 /* didn't match anything previous */
136 id = *numid;
137 (*numid)++;
138 }
139 if (_nrrdCC_verb) {
140 fprintf(stderr__stderrp, "%s: pvl: %g %g %g %g (vl = %g)\n", me,
141 pvl[1], pvl[2], pvl[3], pvl[4], vl);
142 fprintf(stderr__stderrp, " pid: %d %d %d %d\n",
143 pid[1], pid[2], pid[3], pid[4]);
144 fprintf(stderr__stderrp, " --> p = %d, id = %d, *numid = %d\n",
145 p, id, *numid);
146 }
147 out[x + sx*y] = id;
148 }
149 }
150
151 return 0;
152}
153
154/*
155**
156** 5 6 7 --> X
157** 8 9 10
158** 11 12 13
159** |
160** v Y
161** 2 3 4
162** / 1 . . again, 0 index never used, for reasons forgotten
163** Z . . .
164*/
165int
166_nrrdCCFind_3(Nrrd *nout, unsigned int *numid, airArray *eqvArr,
167 const Nrrd *nin, unsigned int conny) {
168 /* static const char me[]="_nrrdCCFind_3" ; */
169 double pvl[14]={0,0,0,0,0,0,0,0,0,0,0,0,0,0}, vl=0;
170 unsigned int id, *out, (*lup)(const void *, size_t),
171 pid[14]={0,0,0,0,0,0,0,0,0,0,0,0,0,0};
172 unsigned int p, x, y, z, sx, sy, sz;
173
174 id = 0; /* sssh! compiler warnings */
175 lup = nrrdUILookup[nin->type];
176 out = AIR_CAST(unsigned int*, nout->data)((unsigned int*)(nout->data));
177 sx = AIR_CAST(unsigned int, nin->axis[0].size)((unsigned int)(nin->axis[0].size));
178 sy = AIR_CAST(unsigned int, nin->axis[1].size)((unsigned int)(nin->axis[1].size));
179 sz = AIR_CAST(unsigned int, nin->axis[2].size)((unsigned int)(nin->axis[2].size));
180#define GETV_3(x,y,z)((((0) <= (((int)(x))) && (((int)(x))) <= (((int
)(sx-1)))) && ((0) <= (((int)(y))) && (((int
)(y))) <= (((int)(sy-1)))) && ((0) <= (((int)(z
))) && (((int)(z))) <= (((int)(sz-1))))) ? lup(nin
->data, (x) + sx*((y) + sy*(z))) : 0.5) ((AIR_IN_CL(0, AIR_CAST(int, x), AIR_CAST(int, sx-1))((0) <= (((int)(x))) && (((int)(x))) <= (((int)
(sx-1)))) \
181 && AIR_IN_CL(0, AIR_CAST(int, y), AIR_CAST(int, sy-1))((0) <= (((int)(y))) && (((int)(y))) <= (((int)
(sy-1)))) \
182 && AIR_IN_CL(0, AIR_CAST(int, z), AIR_CAST(int, sz-1))((0) <= (((int)(z))) && (((int)(z))) <= (((int)
(sz-1)))))\
183 ? lup(nin->data, (x) + sx*((y) + sy*(z))) \
184 : 0.5)
185#define GETI_3(x,y,z)((((0) <= (((int)(x))) && (((int)(x))) <= (((int
)(sx-1)))) && ((0) <= (((int)(y))) && (((int
)(y))) <= (((int)(sy-1)))) && ((0) <= (((int)(z
))) && (((int)(z))) <= (((int)(sz-1))))) ? out[(x)
+ sx*((y) + sy*(z))] : ((unsigned int)(-1))) ((AIR_IN_CL(0, AIR_CAST(int, x), AIR_CAST(int, sx-1))((0) <= (((int)(x))) && (((int)(x))) <= (((int)
(sx-1)))) \
186 && AIR_IN_CL(0, AIR_CAST(int, y), AIR_CAST(int, sy-1))((0) <= (((int)(y))) && (((int)(y))) <= (((int)
(sy-1)))) \
187 && AIR_IN_CL(0, AIR_CAST(int, z), AIR_CAST(int, sz-1))((0) <= (((int)(z))) && (((int)(z))) <= (((int)
(sz-1)))))\
188 ? out[(x) + sx*((y) + sy*(z))] \
189 : AIR_CAST(unsigned int, -1)((unsigned int)(-1)))
190
191 *numid = 0;
192 for (z=0; z<sz; z++) {
193 for (y=0; y<sy; y++) {
194 for (x=0; x<sx; x++) {
195 if (!x) {
196 pvl[ 1] = GETV_3( -1, y, z)((((0) <= (((int)(-1))) && (((int)(-1))) <= (((
int)(sx-1)))) && ((0) <= (((int)(y))) && (
((int)(y))) <= (((int)(sy-1)))) && ((0) <= (((int
)(z))) && (((int)(z))) <= (((int)(sz-1))))) ? lup(
nin->data, (-1) + sx*((y) + sy*(z))) : 0.5); pid[ 1] = GETI_3( -1, y, z)((((0) <= (((int)(-1))) && (((int)(-1))) <= (((
int)(sx-1)))) && ((0) <= (((int)(y))) && (
((int)(y))) <= (((int)(sy-1)))) && ((0) <= (((int
)(z))) && (((int)(z))) <= (((int)(sz-1))))) ? out[
(-1) + sx*((y) + sy*(z))] : ((unsigned int)(-1)));
197 pvl[ 2] = GETV_3( -1, y-1, z)((((0) <= (((int)(-1))) && (((int)(-1))) <= (((
int)(sx-1)))) && ((0) <= (((int)(y-1))) &&
(((int)(y-1))) <= (((int)(sy-1)))) && ((0) <= (
((int)(z))) && (((int)(z))) <= (((int)(sz-1))))) ?
lup(nin->data, (-1) + sx*((y-1) + sy*(z))) : 0.5); pid[ 2] = GETI_3( -1, y-1, z)((((0) <= (((int)(-1))) && (((int)(-1))) <= (((
int)(sx-1)))) && ((0) <= (((int)(y-1))) &&
(((int)(y-1))) <= (((int)(sy-1)))) && ((0) <= (
((int)(z))) && (((int)(z))) <= (((int)(sz-1))))) ?
out[(-1) + sx*((y-1) + sy*(z))] : ((unsigned int)(-1)));
198 pvl[ 3] = GETV_3( 0, y-1, z)((((0) <= (((int)(0))) && (((int)(0))) <= (((int
)(sx-1)))) && ((0) <= (((int)(y-1))) && ((
(int)(y-1))) <= (((int)(sy-1)))) && ((0) <= (((
int)(z))) && (((int)(z))) <= (((int)(sz-1))))) ? lup
(nin->data, (0) + sx*((y-1) + sy*(z))) : 0.5); pid[ 3] = GETI_3( 0, y-1, z)((((0) <= (((int)(0))) && (((int)(0))) <= (((int
)(sx-1)))) && ((0) <= (((int)(y-1))) && ((
(int)(y-1))) <= (((int)(sy-1)))) && ((0) <= (((
int)(z))) && (((int)(z))) <= (((int)(sz-1))))) ? out
[(0) + sx*((y-1) + sy*(z))] : ((unsigned int)(-1)));
199 pvl[ 5] = GETV_3( -1, y-1, z-1)((((0) <= (((int)(-1))) && (((int)(-1))) <= (((
int)(sx-1)))) && ((0) <= (((int)(y-1))) &&
(((int)(y-1))) <= (((int)(sy-1)))) && ((0) <= (
((int)(z-1))) && (((int)(z-1))) <= (((int)(sz-1)))
)) ? lup(nin->data, (-1) + sx*((y-1) + sy*(z-1))) : 0.5); pid[ 5] = GETI_3( -1, y-1, z-1)((((0) <= (((int)(-1))) && (((int)(-1))) <= (((
int)(sx-1)))) && ((0) <= (((int)(y-1))) &&
(((int)(y-1))) <= (((int)(sy-1)))) && ((0) <= (
((int)(z-1))) && (((int)(z-1))) <= (((int)(sz-1)))
)) ? out[(-1) + sx*((y-1) + sy*(z-1))] : ((unsigned int)(-1))
);
200 pvl[ 8] = GETV_3( -1, y, z-1)((((0) <= (((int)(-1))) && (((int)(-1))) <= (((
int)(sx-1)))) && ((0) <= (((int)(y))) && (
((int)(y))) <= (((int)(sy-1)))) && ((0) <= (((int
)(z-1))) && (((int)(z-1))) <= (((int)(sz-1))))) ? lup
(nin->data, (-1) + sx*((y) + sy*(z-1))) : 0.5); pid[ 8] = GETI_3( -1, y, z-1)((((0) <= (((int)(-1))) && (((int)(-1))) <= (((
int)(sx-1)))) && ((0) <= (((int)(y))) && (
((int)(y))) <= (((int)(sy-1)))) && ((0) <= (((int
)(z-1))) && (((int)(z-1))) <= (((int)(sz-1))))) ? out
[(-1) + sx*((y) + sy*(z-1))] : ((unsigned int)(-1)));
201 pvl[11] = GETV_3( -1, y+1, z-1)((((0) <= (((int)(-1))) && (((int)(-1))) <= (((
int)(sx-1)))) && ((0) <= (((int)(y+1))) &&
(((int)(y+1))) <= (((int)(sy-1)))) && ((0) <= (
((int)(z-1))) && (((int)(z-1))) <= (((int)(sz-1)))
)) ? lup(nin->data, (-1) + sx*((y+1) + sy*(z-1))) : 0.5); pid[11] = GETI_3( -1, y+1, z-1)((((0) <= (((int)(-1))) && (((int)(-1))) <= (((
int)(sx-1)))) && ((0) <= (((int)(y+1))) &&
(((int)(y+1))) <= (((int)(sy-1)))) && ((0) <= (
((int)(z-1))) && (((int)(z-1))) <= (((int)(sz-1)))
)) ? out[(-1) + sx*((y+1) + sy*(z-1))] : ((unsigned int)(-1))
);
202 pvl[ 6] = GETV_3( 0, y-1, z-1)((((0) <= (((int)(0))) && (((int)(0))) <= (((int
)(sx-1)))) && ((0) <= (((int)(y-1))) && ((
(int)(y-1))) <= (((int)(sy-1)))) && ((0) <= (((
int)(z-1))) && (((int)(z-1))) <= (((int)(sz-1)))))
? lup(nin->data, (0) + sx*((y-1) + sy*(z-1))) : 0.5); pid[ 6] = GETI_3( 0, y-1, z-1)((((0) <= (((int)(0))) && (((int)(0))) <= (((int
)(sx-1)))) && ((0) <= (((int)(y-1))) && ((
(int)(y-1))) <= (((int)(sy-1)))) && ((0) <= (((
int)(z-1))) && (((int)(z-1))) <= (((int)(sz-1)))))
? out[(0) + sx*((y-1) + sy*(z-1))] : ((unsigned int)(-1)));
203 pvl[ 9] = GETV_3( 0, y, z-1)((((0) <= (((int)(0))) && (((int)(0))) <= (((int
)(sx-1)))) && ((0) <= (((int)(y))) && (((int
)(y))) <= (((int)(sy-1)))) && ((0) <= (((int)(z
-1))) && (((int)(z-1))) <= (((int)(sz-1))))) ? lup
(nin->data, (0) + sx*((y) + sy*(z-1))) : 0.5); pid[ 9] = GETI_3( 0, y, z-1)((((0) <= (((int)(0))) && (((int)(0))) <= (((int
)(sx-1)))) && ((0) <= (((int)(y))) && (((int
)(y))) <= (((int)(sy-1)))) && ((0) <= (((int)(z
-1))) && (((int)(z-1))) <= (((int)(sz-1))))) ? out
[(0) + sx*((y) + sy*(z-1))] : ((unsigned int)(-1)));
204 pvl[12] = GETV_3( 0, y+1, z-1)((((0) <= (((int)(0))) && (((int)(0))) <= (((int
)(sx-1)))) && ((0) <= (((int)(y+1))) && ((
(int)(y+1))) <= (((int)(sy-1)))) && ((0) <= (((
int)(z-1))) && (((int)(z-1))) <= (((int)(sz-1)))))
? lup(nin->data, (0) + sx*((y+1) + sy*(z-1))) : 0.5); pid[12] = GETI_3( 0, y+1, z-1)((((0) <= (((int)(0))) && (((int)(0))) <= (((int
)(sx-1)))) && ((0) <= (((int)(y+1))) && ((
(int)(y+1))) <= (((int)(sy-1)))) && ((0) <= (((
int)(z-1))) && (((int)(z-1))) <= (((int)(sz-1)))))
? out[(0) + sx*((y+1) + sy*(z-1))] : ((unsigned int)(-1)));
205 } else {
206 pvl[ 1] = vl; pid[ 1] = id;
207 pvl[ 2] = pvl[ 3]; pid[ 2] = pid[ 3];
208 pvl[ 3] = pvl[ 4]; pid[ 3] = pid[ 4];
209 pvl[ 5] = pvl[ 6]; pid[ 5] = pid[ 6];
210 pvl[ 8] = pvl[ 9]; pid[ 8] = pid[ 9];
211 pvl[11] = pvl[12]; pid[11] = pid[12];
212 pvl[ 6] = pvl[ 7]; pid[ 6] = pid[ 7];
213 pvl[ 9] = pvl[10]; pid[ 9] = pid[10];
214 pvl[12] = pvl[13]; pid[12] = pid[13];
215 }
216 pvl[ 4] = GETV_3(x+1, y-1, z)((((0) <= (((int)(x+1))) && (((int)(x+1))) <= (
((int)(sx-1)))) && ((0) <= (((int)(y-1))) &&
(((int)(y-1))) <= (((int)(sy-1)))) && ((0) <= (
((int)(z))) && (((int)(z))) <= (((int)(sz-1))))) ?
lup(nin->data, (x+1) + sx*((y-1) + sy*(z))) : 0.5); pid[ 4] = GETI_3(x+1, y-1, z)((((0) <= (((int)(x+1))) && (((int)(x+1))) <= (
((int)(sx-1)))) && ((0) <= (((int)(y-1))) &&
(((int)(y-1))) <= (((int)(sy-1)))) && ((0) <= (
((int)(z))) && (((int)(z))) <= (((int)(sz-1))))) ?
out[(x+1) + sx*((y-1) + sy*(z))] : ((unsigned int)(-1)));
217 pvl[ 7] = GETV_3(x+1, y-1, z-1)((((0) <= (((int)(x+1))) && (((int)(x+1))) <= (
((int)(sx-1)))) && ((0) <= (((int)(y-1))) &&
(((int)(y-1))) <= (((int)(sy-1)))) && ((0) <= (
((int)(z-1))) && (((int)(z-1))) <= (((int)(sz-1)))
)) ? lup(nin->data, (x+1) + sx*((y-1) + sy*(z-1))) : 0.5); pid[ 7] = GETI_3(x+1, y-1, z-1)((((0) <= (((int)(x+1))) && (((int)(x+1))) <= (
((int)(sx-1)))) && ((0) <= (((int)(y-1))) &&
(((int)(y-1))) <= (((int)(sy-1)))) && ((0) <= (
((int)(z-1))) && (((int)(z-1))) <= (((int)(sz-1)))
)) ? out[(x+1) + sx*((y-1) + sy*(z-1))] : ((unsigned int)(-1)
));
218 pvl[10] = GETV_3(x+1, y, z-1)((((0) <= (((int)(x+1))) && (((int)(x+1))) <= (
((int)(sx-1)))) && ((0) <= (((int)(y))) &&
(((int)(y))) <= (((int)(sy-1)))) && ((0) <= ((
(int)(z-1))) && (((int)(z-1))) <= (((int)(sz-1))))
) ? lup(nin->data, (x+1) + sx*((y) + sy*(z-1))) : 0.5); pid[10] = GETI_3(x+1, y, z-1)((((0) <= (((int)(x+1))) && (((int)(x+1))) <= (
((int)(sx-1)))) && ((0) <= (((int)(y))) &&
(((int)(y))) <= (((int)(sy-1)))) && ((0) <= ((
(int)(z-1))) && (((int)(z-1))) <= (((int)(sz-1))))
) ? out[(x+1) + sx*((y) + sy*(z-1))] : ((unsigned int)(-1)));
219 pvl[13] = GETV_3(x+1, y+1, z-1)((((0) <= (((int)(x+1))) && (((int)(x+1))) <= (
((int)(sx-1)))) && ((0) <= (((int)(y+1))) &&
(((int)(y+1))) <= (((int)(sy-1)))) && ((0) <= (
((int)(z-1))) && (((int)(z-1))) <= (((int)(sz-1)))
)) ? lup(nin->data, (x+1) + sx*((y+1) + sy*(z-1))) : 0.5); pid[13] = GETI_3(x+1, y+1, z-1)((((0) <= (((int)(x+1))) && (((int)(x+1))) <= (
((int)(sx-1)))) && ((0) <= (((int)(y+1))) &&
(((int)(y+1))) <= (((int)(sy-1)))) && ((0) <= (
((int)(z-1))) && (((int)(z-1))) <= (((int)(sz-1)))
)) ? out[(x+1) + sx*((y+1) + sy*(z-1))] : ((unsigned int)(-1)
));
220 vl = GETV_3(x, y, z)((((0) <= (((int)(x))) && (((int)(x))) <= (((int
)(sx-1)))) && ((0) <= (((int)(y))) && (((int
)(y))) <= (((int)(sy-1)))) && ((0) <= (((int)(z
))) && (((int)(z))) <= (((int)(sz-1))))) ? lup(nin
->data, (x) + sx*((y) + sy*(z))) : 0.5);
221 p = 0;
222 if (vl == pvl[1]) {
223 id = pid[p=1];
224 }
225 TEST(3)if (vl == pvl[(3)]) { if (p) { if (id != pid[(3)]) { airEqvAdd
(eqvArr, pid[(3)], id); } } else { id = pid[p=(3)]; } };
226 TEST(9)if (vl == pvl[(9)]) { if (p) { if (id != pid[(9)]) { airEqvAdd
(eqvArr, pid[(9)], id); } } else { id = pid[p=(9)]; } };
227 if (2 <= conny) {
228 TEST(2)if (vl == pvl[(2)]) { if (p) { if (id != pid[(2)]) { airEqvAdd
(eqvArr, pid[(2)], id); } } else { id = pid[p=(2)]; } }; TEST(4)if (vl == pvl[(4)]) { if (p) { if (id != pid[(4)]) { airEqvAdd
(eqvArr, pid[(4)], id); } } else { id = pid[p=(4)]; } };
229 TEST(6)if (vl == pvl[(6)]) { if (p) { if (id != pid[(6)]) { airEqvAdd
(eqvArr, pid[(6)], id); } } else { id = pid[p=(6)]; } }; TEST(8)if (vl == pvl[(8)]) { if (p) { if (id != pid[(8)]) { airEqvAdd
(eqvArr, pid[(8)], id); } } else { id = pid[p=(8)]; } }; TEST(10)if (vl == pvl[(10)]) { if (p) { if (id != pid[(10)]) { airEqvAdd
(eqvArr, pid[(10)], id); } } else { id = pid[p=(10)]; } }; TEST(12)if (vl == pvl[(12)]) { if (p) { if (id != pid[(12)]) { airEqvAdd
(eqvArr, pid[(12)], id); } } else { id = pid[p=(12)]; } };
230 if (3 == conny) {
231 TEST(5)if (vl == pvl[(5)]) { if (p) { if (id != pid[(5)]) { airEqvAdd
(eqvArr, pid[(5)], id); } } else { id = pid[p=(5)]; } }; TEST(7)if (vl == pvl[(7)]) { if (p) { if (id != pid[(7)]) { airEqvAdd
(eqvArr, pid[(7)], id); } } else { id = pid[p=(7)]; } }; TEST(11)if (vl == pvl[(11)]) { if (p) { if (id != pid[(11)]) { airEqvAdd
(eqvArr, pid[(11)], id); } } else { id = pid[p=(11)]; } }; TEST(13)if (vl == pvl[(13)]) { if (p) { if (id != pid[(13)]) { airEqvAdd
(eqvArr, pid[(13)], id); } } else { id = pid[p=(13)]; } };
232 }
233 }
234 if (!p) {
235 /* didn't match anything previous */
236 id = *numid;
237 (*numid)++;
238 }
239 out[x + sx*(y + sy*z)] = id;
240 }
241 }
242 }
243
244 return 0;
245}
246
247int
248_nrrdCCFind_N(Nrrd *nout, unsigned int *numid, airArray *eqvArr,
249 const Nrrd *nin, unsigned int conny) {
250 static const char me[]="_nrrdCCFind_N";
251
252 AIR_UNUSED(nout)(void)(nout);
253 AIR_UNUSED(numid)(void)(numid);
254 AIR_UNUSED(eqvArr)(void)(eqvArr);
255 AIR_UNUSED(nin)(void)(nin);
256 AIR_UNUSED(conny)(void)(conny);
257 biffAddf(NRRDnrrdBiffKey, "%s: sorry, not implemented yet", me);
258 return 1;
259}
260
261/*
262******** nrrdCCFind
263**
264** finds connected components (CCs) in given integral type nrrd "nin",
265** according to connectivity "conny", putting the results in "nout".
266** The "type" argument controls what type the output will be. If
267** type == nrrdTypeDefault, the type used will be the smallest that
268** can contain the CC id values. Otherwise, the specified type "type"
269** will be used, assuming that it is large enough to hold the CC ids.
270**
271** "conny": the number of coordinates that need to varied together in
272** order to reach all the samples that are to consitute the neighborhood
273** around a sample. For 2-D, conny==1 specifies the 4 edge-connected
274** pixels, and 2 specifies the 8 edge- and corner-connected.
275**
276** The caller can get a record of the values in each CC by passing a
277** non-NULL nval, which will be allocated to an array of the same type
278** as nin, so that nval->data[I] is the value in nin inside CC #I.
279*/
280int
281nrrdCCFind(Nrrd *nout, Nrrd **nvalP, const Nrrd *nin, int type,
282 unsigned int conny) {
283 static const char me[]="nrrdCCFind", func[]="ccfind";
284 Nrrd *nfpid; /* first-pass IDs */
285 airArray *mop, *eqvArr;
286 unsigned int *fpid, numid, numsettleid, *map,
287 (*lup)(const void *, size_t), (*ins)(void *, size_t, unsigned int);
288 int ret;
289 size_t I, NN;
290 void *val;
291
292 if (!(nout && nin)) {
293 /* NULL nvalP okay */
294 biffAddf(NRRDnrrdBiffKey, "%s: got NULL pointer", me);
295 return 1;
296 }
297 if (nout == nin) {
298 biffAddf(NRRDnrrdBiffKey, "%s: nout == nin disallowed", me);
299 return 1;
300 }
301 if (!( nrrdTypeIsIntegral[nin->type]
302 && nrrdTypeIsUnsigned[nin->type]
303 && nrrdTypeSize[nin->type] <= 4 )) {
304 biffAddf(NRRDnrrdBiffKey, "%s: can only find connected components in "
305 "1, 2, or 4 byte unsigned integral values (not %s)",
306 me, airEnumStr(nrrdType, nin->type));
307 return 1;
308 }
309 if (nrrdTypeDefault != type) {
310 if (!( AIR_IN_OP(nrrdTypeUnknown, type, nrrdTypeLast)((nrrdTypeUnknown) < (type) && (type) < (nrrdTypeLast
)) )) {
311 biffAddf(NRRDnrrdBiffKey, "%s: got invalid target type %d", me, type);
312 return 1;
313 }
314 if (!( nrrdTypeIsIntegral[type]
315 && nrrdTypeIsUnsigned[nin->type]
316 && nrrdTypeSize[type] <= 4 )) {
317 biffAddf(NRRDnrrdBiffKey,
318 "%s: can only save connected components to 1, 2, or 4 byte "
319 "unsigned integral values (not %s)",
320 me, airEnumStr(nrrdType, type));
321 return 1;
322 }
323 }
324 if (!( conny <= nin->dim )) {
325 biffAddf(NRRDnrrdBiffKey, "%s: connectivity value must be in [1..%d] for %d-D "
326 "data (not %d)", me, nin->dim, nin->dim, conny);
327 return 1;
328 }
329 if (nrrdConvert(nfpid=nrrdNew(), nin, nrrdTypeUInt)) {
330 biffAddf(NRRDnrrdBiffKey, "%s: couldn't allocate fpid %s array to match input size",
331 me, airEnumStr(nrrdType, nrrdTypeUInt));
332 return 1;
333 }
334
335 mop = airMopNew();
336 airMopAdd(mop, nfpid, (airMopper)nrrdNuke, airMopAlways);
337 eqvArr = airArrayNew(NULL((void*)0), NULL((void*)0), 2*sizeof(unsigned int), _nrrdCC_EqvIncr);
338 airMopAdd(mop, eqvArr, (airMopper)airArrayNuke, airMopAlways);
339 ret = 0;
Value stored to 'ret' is never read
340 switch(nin->dim) {
341 case 1:
342 ret = _nrrdCCFind_1(nfpid, &numid, nin);
343 break;
344 case 2:
345 ret = _nrrdCCFind_2(nfpid, &numid, eqvArr, nin, conny);
346 break;
347 case 3:
348 ret = _nrrdCCFind_3(nfpid, &numid, eqvArr, nin, conny);
349 break;
350 default:
351 ret = _nrrdCCFind_N(nfpid, &numid, eqvArr, nin, conny);
352 break;
353 }
354 if (ret) {
355 biffAddf(NRRDnrrdBiffKey, "%s: initial pass failed", me);
356 airMopError(mop); return 1;
357 }
358
359 map = AIR_MALLOC(numid, unsigned int)(unsigned int*)(malloc((numid)*sizeof(unsigned int)));
360 airMopAdd(mop, map, airFree, airMopAlways);
361 numsettleid = airEqvMap(eqvArr, map, numid);
362 /* convert fpid values to final id values */
363 fpid = (unsigned int*)(nfpid->data);
364 NN = nrrdElementNumber(nfpid);
365 for (I=0; I<NN; I++) {
366 fpid[I] = map[fpid[I]];
367 }
368 if (nvalP) {
369 if (!(*nvalP)) {
370 *nvalP = nrrdNew();
371 }
372 if (nrrdMaybeAlloc_va(*nvalP, nin->type, 1,
373 AIR_CAST(size_t, numsettleid)((size_t)(numsettleid)))) {
374 biffAddf(NRRDnrrdBiffKey, "%s: couldn't allocate output value list", me);
375 airMopError(mop); return 1;
376 }
377 airMopAdd(mop, nvalP, (airMopper)airSetNull, airMopOnError);
378 airMopAdd(mop, *nvalP, (airMopper)nrrdNuke, airMopOnError);
379 val = (*nvalP)->data;
380 lup = nrrdUILookup[nin->type];
381 ins = nrrdUIInsert[nin->type];
382 /* I'm not sure if its more work to do all the redundant assignments
383 or to check whether or not to do them */
384 for (I=0; I<NN; I++) {
385 ins(val, fpid[I], lup(nin->data, I));
386 }
387 }
388
389 if (nrrdTypeDefault != type) {
390 if (numsettleid-1 > nrrdTypeMax[type]) {
391 biffAddf(NRRDnrrdBiffKey,
392 "%s: max cc id %u is too large to fit in output type %s",
393 me, numsettleid-1, airEnumStr(nrrdType, type));
394 airMopError(mop); return 1;
395 }
396 } else {
397 type = (numsettleid-1 <= nrrdTypeMax[nrrdTypeUChar]
398 ? nrrdTypeUChar
399 : (numsettleid-1 <= nrrdTypeMax[nrrdTypeUShort]
400 ? nrrdTypeUShort
401 : nrrdTypeUInt));
402 }
403 if (nrrdConvert(nout, nfpid, type)) {
404 biffAddf(NRRDnrrdBiffKey, "%s: trouble converting to final output", me);
405 airMopError(mop); return 1;
406 }
407 if (nrrdContentSet_va(nout, func, nin, "%s,%d",
408 airEnumStr(nrrdType, type), conny)) {
409 biffAddf(NRRDnrrdBiffKey, "%s:", me);
410 return 1;
411 }
412 if (nout != nin) {
413 nrrdAxisInfoCopy(nout, nin, NULL((void*)0), NRRD_AXIS_INFO_NONE0);
414 }
415 /* basic info handled by nrrdConvert */
416
417 airMopOkay(mop);
418 return 0;
419}
420
421int
422_nrrdCCAdj_1(unsigned char *out, int numid, const Nrrd *nin) {
423
424 AIR_UNUSED(out)(void)(out);
425 AIR_UNUSED(numid)(void)(numid);
426 AIR_UNUSED(nin)(void)(nin);
427 return 0;
428}
429
430int
431_nrrdCCAdj_2(unsigned char *out, unsigned int numid, const Nrrd *nin,
432 unsigned int conny) {
433 unsigned int (*lup)(const void *, size_t), x, y, sx, sy, id=0;
434 double pid[5]={0,0,0,0,0};
435
436 lup = nrrdUILookup[nin->type];
437 sx = AIR_CAST(unsigned int, nin->axis[0].size)((unsigned int)(nin->axis[0].size));
438 sy = AIR_CAST(unsigned int, nin->axis[1].size)((unsigned int)(nin->axis[1].size));
439 for (y=0; y<sy; y++) {
440 for (x=0; x<sx; x++) {
441 if (!x) {
442 pid[1] = GETV_2(-1, y)((((0) <= (((int)(-1))) && (((int)(-1))) <= (((
int)(sx-1)))) && ((0) <= (((int)(y))) && (
((int)(y))) <= (((int)(sy-1))))) ? lup(nin->data, (-1) +
sx*(y)) : 0.5);
443 pid[2] = GETV_2(-1, y-1)((((0) <= (((int)(-1))) && (((int)(-1))) <= (((
int)(sx-1)))) && ((0) <= (((int)(y-1))) &&
(((int)(y-1))) <= (((int)(sy-1))))) ? lup(nin->data, (
-1) + sx*(y-1)) : 0.5);
444 pid[3] = GETV_2(0, y-1)((((0) <= (((int)(0))) && (((int)(0))) <= (((int
)(sx-1)))) && ((0) <= (((int)(y-1))) && ((
(int)(y-1))) <= (((int)(sy-1))))) ? lup(nin->data, (0) +
sx*(y-1)) : 0.5);
445 } else {
446 pid[1] = id;
447 pid[2] = pid[3];
448 pid[3] = pid[4];
449 }
450 pid[4] = GETV_2(x+1, y-1)((((0) <= (((int)(x+1))) && (((int)(x+1))) <= (
((int)(sx-1)))) && ((0) <= (((int)(y-1))) &&
(((int)(y-1))) <= (((int)(sy-1))))) ? lup(nin->data, (
x+1) + sx*(y-1)) : 0.5);
451 id = AIR_CAST(unsigned int, GETV_2(x, y))((unsigned int)(((((0) <= (((int)(x))) && (((int)(
x))) <= (((int)(sx-1)))) && ((0) <= (((int)(y))
) && (((int)(y))) <= (((int)(sy-1))))) ? lup(nin->
data, (x) + sx*(y)) : 0.5)));
452#define TADJ(P)if (pid[(P)] != 0.5 && id != pid[(P)]) { out[id + numid
*((unsigned int)(pid[(P)]))] = out[((unsigned int)(pid[(P)]))
+ numid*id] = 1; } \
453 if (pid[(P)] != 0.5 && id != pid[(P)]) { \
454 out[id + numid*AIR_CAST(unsigned int, pid[(P)])((unsigned int)(pid[(P)]))] = \
455 out[AIR_CAST(unsigned int, pid[(P)])((unsigned int)(pid[(P)])) + numid*id] = 1; \
456 }
457 TADJ(1)if (pid[(1)] != 0.5 && id != pid[(1)]) { out[id + numid
*((unsigned int)(pid[(1)]))] = out[((unsigned int)(pid[(1)]))
+ numid*id] = 1; };
458 TADJ(3)if (pid[(3)] != 0.5 && id != pid[(3)]) { out[id + numid
*((unsigned int)(pid[(3)]))] = out[((unsigned int)(pid[(3)]))
+ numid*id] = 1; };
459 if (2 == conny) {
460 TADJ(2)if (pid[(2)] != 0.5 && id != pid[(2)]) { out[id + numid
*((unsigned int)(pid[(2)]))] = out[((unsigned int)(pid[(2)]))
+ numid*id] = 1; };
461 TADJ(4)if (pid[(4)] != 0.5 && id != pid[(4)]) { out[id + numid
*((unsigned int)(pid[(4)]))] = out[((unsigned int)(pid[(4)]))
+ numid*id] = 1; };
462 }
463 }
464 }
465
466 return 0;
467}
468
469int
470_nrrdCCAdj_3(unsigned char *out, int numid, const Nrrd *nin,
471 unsigned int conny) {
472 unsigned int (*lup)(const void *, size_t), x, y, z, sx, sy, sz, id=0;
473 double pid[14]={0,0,0,0,0,0,0,0,0,0,0,0,0,0};
474
475 lup = nrrdUILookup[nin->type];
476 sx = AIR_CAST(unsigned int, nin->axis[0].size)((unsigned int)(nin->axis[0].size));
477 sy = AIR_CAST(unsigned int, nin->axis[1].size)((unsigned int)(nin->axis[1].size));
478 sz = AIR_CAST(unsigned int, nin->axis[2].size)((unsigned int)(nin->axis[2].size));
479 for (z=0; z<sz; z++) {
480 for (y=0; y<sy; y++) {
481 for (x=0; x<sx; x++) {
482 if (!x) {
483 pid[ 1] = GETV_3(-1, y, z)((((0) <= (((int)(-1))) && (((int)(-1))) <= (((
int)(sx-1)))) && ((0) <= (((int)(y))) && (
((int)(y))) <= (((int)(sy-1)))) && ((0) <= (((int
)(z))) && (((int)(z))) <= (((int)(sz-1))))) ? lup(
nin->data, (-1) + sx*((y) + sy*(z))) : 0.5);
484 pid[ 2] = GETV_3(-1, y-1, z)((((0) <= (((int)(-1))) && (((int)(-1))) <= (((
int)(sx-1)))) && ((0) <= (((int)(y-1))) &&
(((int)(y-1))) <= (((int)(sy-1)))) && ((0) <= (
((int)(z))) && (((int)(z))) <= (((int)(sz-1))))) ?
lup(nin->data, (-1) + sx*((y-1) + sy*(z))) : 0.5);
485 pid[ 3] = GETV_3( 0, y-1, z)((((0) <= (((int)(0))) && (((int)(0))) <= (((int
)(sx-1)))) && ((0) <= (((int)(y-1))) && ((
(int)(y-1))) <= (((int)(sy-1)))) && ((0) <= (((
int)(z))) && (((int)(z))) <= (((int)(sz-1))))) ? lup
(nin->data, (0) + sx*((y-1) + sy*(z))) : 0.5);
486 pid[ 5] = GETV_3(-1, y-1, z-1)((((0) <= (((int)(-1))) && (((int)(-1))) <= (((
int)(sx-1)))) && ((0) <= (((int)(y-1))) &&
(((int)(y-1))) <= (((int)(sy-1)))) && ((0) <= (
((int)(z-1))) && (((int)(z-1))) <= (((int)(sz-1)))
)) ? lup(nin->data, (-1) + sx*((y-1) + sy*(z-1))) : 0.5);
487 pid[ 8] = GETV_3(-1, y, z-1)((((0) <= (((int)(-1))) && (((int)(-1))) <= (((
int)(sx-1)))) && ((0) <= (((int)(y))) && (
((int)(y))) <= (((int)(sy-1)))) && ((0) <= (((int
)(z-1))) && (((int)(z-1))) <= (((int)(sz-1))))) ? lup
(nin->data, (-1) + sx*((y) + sy*(z-1))) : 0.5);
488 pid[11] = GETV_3(-1, y+1, z-1)((((0) <= (((int)(-1))) && (((int)(-1))) <= (((
int)(sx-1)))) && ((0) <= (((int)(y+1))) &&
(((int)(y+1))) <= (((int)(sy-1)))) && ((0) <= (
((int)(z-1))) && (((int)(z-1))) <= (((int)(sz-1)))
)) ? lup(nin->data, (-1) + sx*((y+1) + sy*(z-1))) : 0.5);
489 pid[ 6] = GETV_3( 0, y-1, z-1)((((0) <= (((int)(0))) && (((int)(0))) <= (((int
)(sx-1)))) && ((0) <= (((int)(y-1))) && ((
(int)(y-1))) <= (((int)(sy-1)))) && ((0) <= (((
int)(z-1))) && (((int)(z-1))) <= (((int)(sz-1)))))
? lup(nin->data, (0) + sx*((y-1) + sy*(z-1))) : 0.5);
490 pid[ 9] = GETV_3( 0, y, z-1)((((0) <= (((int)(0))) && (((int)(0))) <= (((int
)(sx-1)))) && ((0) <= (((int)(y))) && (((int
)(y))) <= (((int)(sy-1)))) && ((0) <= (((int)(z
-1))) && (((int)(z-1))) <= (((int)(sz-1))))) ? lup
(nin->data, (0) + sx*((y) + sy*(z-1))) : 0.5);
491 pid[12] = GETV_3( 0, y+1, z-1)((((0) <= (((int)(0))) && (((int)(0))) <= (((int
)(sx-1)))) && ((0) <= (((int)(y+1))) && ((
(int)(y+1))) <= (((int)(sy-1)))) && ((0) <= (((
int)(z-1))) && (((int)(z-1))) <= (((int)(sz-1)))))
? lup(nin->data, (0) + sx*((y+1) + sy*(z-1))) : 0.5);
492 } else {
493 pid[ 1] = id;
494 pid[ 2] = pid[ 3];
495 pid[ 3] = pid[ 4];
496 pid[ 5] = pid[ 6];
497 pid[ 8] = pid[ 9];
498 pid[11] = pid[12];
499 pid[ 6] = pid[ 7];
500 pid[ 9] = pid[10];
501 pid[12] = pid[13];
502 }
503 pid[ 4] = GETV_3(x+1, y-1, z)((((0) <= (((int)(x+1))) && (((int)(x+1))) <= (
((int)(sx-1)))) && ((0) <= (((int)(y-1))) &&
(((int)(y-1))) <= (((int)(sy-1)))) && ((0) <= (
((int)(z))) && (((int)(z))) <= (((int)(sz-1))))) ?
lup(nin->data, (x+1) + sx*((y-1) + sy*(z))) : 0.5);
504 pid[ 7] = GETV_3(x+1, y-1, z-1)((((0) <= (((int)(x+1))) && (((int)(x+1))) <= (
((int)(sx-1)))) && ((0) <= (((int)(y-1))) &&
(((int)(y-1))) <= (((int)(sy-1)))) && ((0) <= (
((int)(z-1))) && (((int)(z-1))) <= (((int)(sz-1)))
)) ? lup(nin->data, (x+1) + sx*((y-1) + sy*(z-1))) : 0.5);
505 pid[10] = GETV_3(x+1, y, z-1)((((0) <= (((int)(x+1))) && (((int)(x+1))) <= (
((int)(sx-1)))) && ((0) <= (((int)(y))) &&
(((int)(y))) <= (((int)(sy-1)))) && ((0) <= ((
(int)(z-1))) && (((int)(z-1))) <= (((int)(sz-1))))
) ? lup(nin->data, (x+1) + sx*((y) + sy*(z-1))) : 0.5);
506 pid[13] = GETV_3(x+1, y+1, z-1)((((0) <= (((int)(x+1))) && (((int)(x+1))) <= (
((int)(sx-1)))) && ((0) <= (((int)(y+1))) &&
(((int)(y+1))) <= (((int)(sy-1)))) && ((0) <= (
((int)(z-1))) && (((int)(z-1))) <= (((int)(sz-1)))
)) ? lup(nin->data, (x+1) + sx*((y+1) + sy*(z-1))) : 0.5);
507 id = AIR_CAST(unsigned int, GETV_3(x, y, z))((unsigned int)(((((0) <= (((int)(x))) && (((int)(
x))) <= (((int)(sx-1)))) && ((0) <= (((int)(y))
) && (((int)(y))) <= (((int)(sy-1)))) && (
(0) <= (((int)(z))) && (((int)(z))) <= (((int)(
sz-1))))) ? lup(nin->data, (x) + sx*((y) + sy*(z))) : 0.5)
));
508 TADJ(1)if (pid[(1)] != 0.5 && id != pid[(1)]) { out[id + numid
*((unsigned int)(pid[(1)]))] = out[((unsigned int)(pid[(1)]))
+ numid*id] = 1; };
509 TADJ(3)if (pid[(3)] != 0.5 && id != pid[(3)]) { out[id + numid
*((unsigned int)(pid[(3)]))] = out[((unsigned int)(pid[(3)]))
+ numid*id] = 1; };
510 TADJ(9)if (pid[(9)] != 0.5 && id != pid[(9)]) { out[id + numid
*((unsigned int)(pid[(9)]))] = out[((unsigned int)(pid[(9)]))
+ numid*id] = 1; };
511 if (2 <= conny) {
512 TADJ(2)if (pid[(2)] != 0.5 && id != pid[(2)]) { out[id + numid
*((unsigned int)(pid[(2)]))] = out[((unsigned int)(pid[(2)]))
+ numid*id] = 1; }; TADJ(4)if (pid[(4)] != 0.5 && id != pid[(4)]) { out[id + numid
*((unsigned int)(pid[(4)]))] = out[((unsigned int)(pid[(4)]))
+ numid*id] = 1; };
513 TADJ(6)if (pid[(6)] != 0.5 && id != pid[(6)]) { out[id + numid
*((unsigned int)(pid[(6)]))] = out[((unsigned int)(pid[(6)]))
+ numid*id] = 1; }; TADJ(8)if (pid[(8)] != 0.5 && id != pid[(8)]) { out[id + numid
*((unsigned int)(pid[(8)]))] = out[((unsigned int)(pid[(8)]))
+ numid*id] = 1; }; TADJ(10)if (pid[(10)] != 0.5 && id != pid[(10)]) { out[id + numid
*((unsigned int)(pid[(10)]))] = out[((unsigned int)(pid[(10)]
)) + numid*id] = 1; }; TADJ(12)if (pid[(12)] != 0.5 && id != pid[(12)]) { out[id + numid
*((unsigned int)(pid[(12)]))] = out[((unsigned int)(pid[(12)]
)) + numid*id] = 1; };
514 if (3 == conny) {
515 TADJ(5)if (pid[(5)] != 0.5 && id != pid[(5)]) { out[id + numid
*((unsigned int)(pid[(5)]))] = out[((unsigned int)(pid[(5)]))
+ numid*id] = 1; }; TADJ(7)if (pid[(7)] != 0.5 && id != pid[(7)]) { out[id + numid
*((unsigned int)(pid[(7)]))] = out[((unsigned int)(pid[(7)]))
+ numid*id] = 1; }; TADJ(11)if (pid[(11)] != 0.5 && id != pid[(11)]) { out[id + numid
*((unsigned int)(pid[(11)]))] = out[((unsigned int)(pid[(11)]
)) + numid*id] = 1; }; TADJ(13)if (pid[(13)] != 0.5 && id != pid[(13)]) { out[id + numid
*((unsigned int)(pid[(13)]))] = out[((unsigned int)(pid[(13)]
)) + numid*id] = 1; };
516 }
517 }
518 }
519 }
520 }
521
522 return 0;
523}
524
525int
526_nrrdCCAdj_N(unsigned char *out, int numid, const Nrrd *nin,
527 unsigned int conny) {
528 static const char me[]="_nrrdCCAdj_N";
529
530 AIR_UNUSED(out)(void)(out);
531 AIR_UNUSED(numid)(void)(numid);
532 AIR_UNUSED(nin)(void)(nin);
533 AIR_UNUSED(conny)(void)(conny);
534 biffAddf(NRRDnrrdBiffKey, "%s: sorry, not implemented", me);
535 return 1;
536}
537
538int
539nrrdCCAdjacency(Nrrd *nout, const Nrrd *nin, unsigned int conny) {
540 static const char me[]="nrrdCCAdjacency", func[]="ccadj";
541 int ret;
542 unsigned int maxid;
543 unsigned char *out;
544
545 if (!( nout && nrrdCCValid(nin) )) {
546 biffAddf(NRRDnrrdBiffKey, "%s: invalid args", me);
547 return 1;
548 }
549 if (nout == nin) {
550 biffAddf(NRRDnrrdBiffKey, "%s: nout == nin disallowed", me);
551 return 1;
552 }
553 if (!( AIR_IN_CL(1, conny, nin->dim)((1) <= (conny) && (conny) <= (nin->dim)) )) {
554 biffAddf(NRRDnrrdBiffKey, "%s: connectivity value must be in [1..%d] for %d-D "
555 "data (not %d)", me, nin->dim, nin->dim, conny);
556 return 1;
557 }
558 maxid = nrrdCCMax(nin);
559 if (nrrdMaybeAlloc_va(nout, nrrdTypeUChar, 2,
560 AIR_CAST(size_t, maxid+1)((size_t)(maxid+1)),
561 AIR_CAST(size_t, maxid+1)((size_t)(maxid+1)))) {
562 biffAddf(NRRDnrrdBiffKey, "%s: trouble allocating output", me);
563 return 1;
564 }
565 out = (unsigned char *)(nout->data);
566
567 switch(nin->dim) {
568 case 1:
569 ret = _nrrdCCAdj_1(out, maxid+1, nin);
570 break;
571 case 2:
572 ret = _nrrdCCAdj_2(out, maxid+1, nin, conny);
573 break;
574 case 3:
575 ret = _nrrdCCAdj_3(out, maxid+1, nin, conny);
576 break;
577 default:
578 ret = _nrrdCCAdj_N(out, maxid+1, nin, conny);
579 break;
580 }
581 if (ret) {
582 biffAddf(NRRDnrrdBiffKey, "%s: trouble", me);
583 return 1;
584 }
585 /* this goofiness is just so that histo-based projections
586 return the sorts of values that we expect */
587 nout->axis[0].center = nout->axis[1].center = nrrdCenterCell;
588 nout->axis[0].min = nout->axis[1].min = -0.5;
589 nout->axis[0].max = nout->axis[1].max = maxid + 0.5;
590 if (nrrdContentSet_va(nout, func, nin, "%d", conny)) {
591 biffAddf(NRRDnrrdBiffKey, "%s:", me);
592 return 1;
593 }
594
595 return 0;
596}
597
598/*
599******** nrrdCCMerge
600**
601** Slightly-too-multi-purpose tool for merging small connected components
602** (CCs) into larger ones, according to a number of possible different
603** constraints, as explained below.
604**
605** valDir: (value direction) uses information about the original values
606** in the CC to constrain whether darker gets merged into brighter, or vice
607** versa, or neither. For non-zero valDir values, a non-NULL _nval (from
608** nrrdCCFind) must be passed.
609** valDir > 0 : merge dark CCs into bright, but not vice versa
610** valDir = 0 : merge either way, values are irrelevant
611** valDir < 0 : merge bright CCs into dark, but not vice versa
612** When merging with multiple neighbors (maxNeighbor > 1), the value
613** of the largest neighbor is considered.
614**
615** maxSize: a cap on how large "small" is- CCs any larger than maxSize are
616** not merged, as they are deemed too significant. Or, a maxSize of 0 says
617** size is no object for merging CCs.
618**
619** maxNeighbor: a maximum number of neighbors that a CC can have (either
620** bigger than the CC or not) if it is to be merged. Use 1 to merge
621** isolated islands into their surrounds, 2 to merge CC with the larger
622** of their two neighbors, etc., or 0 to allow any number of neighbors.
623**
624** conny: passed to nrrdCCAdjacency() when determining neighbors
625**
626** In order to prevent weirdness, the merging done in one call to this
627** function is not transitive: if A is merged to B, then B will not be
628** merged to anything else, even if meets all the requirements defined
629** by the given parameters. This is accomplished by working from the
630** smallest CCs to the largest. Iterated calls may be needed to acheive
631** the desired effect.
632**
633** Note: the output of this is not "settled"- the CC id values are not
634** shiftward downwards to their lowest possible values, since this would
635** needlessly invalidate the nval value store.
636*/
637int
638nrrdCCMerge(Nrrd *nout, const Nrrd *nin, Nrrd *_nval,
639 int valDir, unsigned int maxSize, unsigned int maxNeighbor,
640 unsigned int conny) {
641 static const char me[]="nrrdCCMerge", func[]="ccmerge";
642 const char *valcnt;
643 unsigned int _i, i, j, bigi=0, numid, *size, *sizeId,
644 *nn, /* number of neighbors */
645 *val=NULL((void*)0), *hit,
646 (*lup)(const void *, size_t), (*ins)(void *, size_t, unsigned int);
647 Nrrd *nadj, *nsize, *nval=NULL((void*)0), *nnn;
648 unsigned char *adj;
649 unsigned int *map, *id;
650 airArray *mop;
651 size_t I, NN;
652
653 mop = airMopNew();
654 if (!( nout && nrrdCCValid(nin) )) {
655 /* _nval can be NULL */
656 biffAddf(NRRDnrrdBiffKey, "%s: invalid args", me);
657 airMopError(mop); return 1;
658 }
659 if (valDir) {
660 airMopAdd(mop, nval = nrrdNew(), (airMopper)nrrdNuke, airMopAlways);
661 if (nrrdConvert(nval, _nval, nrrdTypeUInt)) {
662 biffAddf(NRRDnrrdBiffKey, "%s: value-directed merging needs usable nval", me);
663 airMopError(mop); return 1;
664 }
665 val = (unsigned int*)(nval->data);
666 }
667 if (nout != nin) {
668 if (nrrdCopy(nout, nin)) {
669 biffAddf(NRRDnrrdBiffKey, "%s:", me);
670 airMopError(mop); return 1;
671 }
672 }
673 airMopAdd(mop, nadj = nrrdNew(), (airMopper)nrrdNuke, airMopAlways);
674 airMopAdd(mop, nsize = nrrdNew(), (airMopper)nrrdNuke, airMopAlways);
675 airMopAdd(mop, nnn = nrrdNew(), (airMopper)nrrdNuke, airMopAlways);
676
677 if (nrrdCCSize(nsize, nin)
678 || nrrdCopy(nnn, nsize) /* just to allocate to right size and type */
679 || nrrdCCAdjacency(nadj, nin, conny)) {
680 biffAddf(NRRDnrrdBiffKey, "%s:", me);
681 airMopError(mop); return 1;
682 }
683 size = (unsigned int*)(nsize->data);
684 adj = (unsigned char*)(nadj->data);
685 nn = (unsigned int*)(nnn->data);
686 numid = AIR_CAST(unsigned int, nsize->axis[0].size)((unsigned int)(nsize->axis[0].size));
687 for (i=0; i<numid; i++) {
688 nn[i] = 0;
689 for (j=0; j<numid; j++) {
690 nn[i] += adj[j + numid*i];
691 }
692 }
693 map = AIR_MALLOC(numid, unsigned int)(unsigned int*)(malloc((numid)*sizeof(unsigned int)));
694 id = AIR_MALLOC(numid, unsigned int)(unsigned int*)(malloc((numid)*sizeof(unsigned int)));
695 hit = AIR_MALLOC(numid, unsigned int)(unsigned int*)(malloc((numid)*sizeof(unsigned int)));
696 sizeId = AIR_MALLOC(2*numid, unsigned int)(unsigned int*)(malloc((2*numid)*sizeof(unsigned int)));
697 /* we add to the mops BEFORE error checking so that anything non-NULL
698 will get airFree'd, and happily airFree is a no-op on NULL */
699 airMopAdd(mop, map, airFree, airMopAlways);
700 airMopAdd(mop, id, airFree, airMopAlways);
701 airMopAdd(mop, hit, airFree, airMopAlways);
702 airMopAdd(mop, sizeId, airFree, airMopAlways);
703 if (!(map && id && hit && sizeId)) {
704 biffAddf(NRRDnrrdBiffKey, "%s: couldn't allocate buffers", me);
705 airMopError(mop); return 1;
706 }
707
708 /* store and sort size/id pairs */
709 for (i=0; i<numid; i++) {
710 sizeId[0 + 2*i] = size[i];
711 sizeId[1 + 2*i] = i;
712 }
713 qsort(sizeId, numid, 2*sizeof(unsigned int), nrrdValCompare[nrrdTypeUInt]);
714 for (i=0; i<numid; i++) {
715 id[i] = sizeId[1 + 2*i];
716 }
717
718 /* initialize arrays */
719 for (i=0; i<numid; i++) {
720 map[i] = i;
721 hit[i] = AIR_FALSE0;
722 }
723 /* _i goes through 0 to numid-1,
724 i goes through the CC ids in ascending order of size */
725 for (_i=0; _i<numid; _i++) {
726 i = id[_i];
727 if (hit[i]) {
728 continue;
729 }
730 if (maxSize && (size[i] > maxSize)) {
731 continue;
732 }
733 if (maxNeighbor && (nn[i] > maxNeighbor)) {
734 continue;
735 }
736 /* find biggest neighbor, exploiting the fact that we already
737 sorted CC ids on size. j descends through indices of id[],
738 bigi goes through CC ids which are larger than CC i */
739 for (j=numid-1; j>_i; j--) {
740 bigi = id[j];
741 if (adj[bigi + numid*i])
742 break;
743 }
744 if (j == _i) {
745 continue; /* we had no neighbors ?!?! */
746 }
747 if (valDir && (AIR_CAST(int, val[bigi])((int)(val[bigi]))
748 - AIR_CAST(int, val[i])((int)(val[i])))*valDir < 0 ) {
749 continue;
750 }
751 /* else all criteria for merging have been met */
752 map[i] = bigi;
753 hit[bigi] = AIR_TRUE1;
754 }
755 lup = nrrdUILookup[nin->type];
756 ins = nrrdUIInsert[nout->type];
757 NN = nrrdElementNumber(nin);
758 for (I=0; I<NN; I++) {
759 ins(nout->data, I, map[lup(nin->data, I)]);
760 }
761
762 valcnt = ((_nval && _nval->content)
763 ? _nval->content
764 : nrrdStateUnknownContent);
765 if ( (valDir && nrrdContentSet_va(nout, func, nin, "%c(%s),%d,%d,%d",
766 (valDir > 0 ? '+' : '-'), valcnt,
767 maxSize, maxNeighbor, conny))
768 ||
769 (!valDir && nrrdContentSet_va(nout, func, nin, ".,%d,%d,%d",
770 maxSize, maxNeighbor, conny)) ) {
771 biffAddf(NRRDnrrdBiffKey, "%s:", me);
772 airMopError(mop); return 1;
773 }
774 /* basic info handled by nrrdCopy */
775 airMopOkay(mop);
776 return 0;
777}
778
779/*
780******** nrrdCCRevalue()
781**
782** assigns the original values back to the connected components
783** obviously, this could be subsumed by nrrdApply1DLut(), but this
784** is so special purpose that it seemed simpler to code from scratch
785*/
786int
787nrrdCCRevalue (Nrrd *nout, const Nrrd *nin, const Nrrd *nval) {
788 static const char me[]="nrrdCCRevalue";
789 size_t I, NN;
790 unsigned int (*vlup)(const void *, size_t), (*ilup)(const void *, size_t),
791 (*ins)(void *, size_t, unsigned int);
792
793 if (!( nout && nrrdCCValid(nin) && nval )) {
794 biffAddf(NRRDnrrdBiffKey, "%s: invalid args", me);
795 return 1;
796 }
797 if (nrrdConvert(nout, nin, nval->type)) {
798 biffAddf(NRRDnrrdBiffKey, "%s: couldn't initialize output", me);
799 return 1;
800 }
801 NN = nrrdElementNumber(nin);
802 vlup = nrrdUILookup[nval->type];
803 ilup = nrrdUILookup[nin->type];
804 ins = nrrdUIInsert[nout->type];
805 for (I=0; I<NN; I++) {
806 ins(nout->data, I, vlup(nval->data, ilup(nin->data, I)));
807 }
808 /* basic info handled by nrrdConvert */
809
810 return 0;
811}
812
813int
814nrrdCCSettle(Nrrd *nout, Nrrd **nvalP, const Nrrd *nin) {
815 static const char me[]="nrrdCCSettle", func[]="ccsettle";
816 unsigned int numid, maxid, jd, id, *map,
817 (*lup)(const void *, size_t), (*ins)(void *, size_t, unsigned int);
818 size_t I, NN;
819 airArray *mop;
820
821 mop = airMopNew();
822 if (!( nout && nrrdCCValid(nin) )) {
823 /* nvalP can be NULL */
824 biffAddf(NRRDnrrdBiffKey, "%s: invalid args", me);
825 airMopError(mop); return 1;
826 }
827 if (nrrdCopy(nout, nin)) {
828 biffAddf(NRRDnrrdBiffKey, "%s: initial copy failed", me);
829 airMopError(mop); return 1;
830 }
831 maxid = nrrdCCMax(nin);
832 lup = nrrdUILookup[nin->type];
833 ins = nrrdUIInsert[nin->type];
834 NN = nrrdElementNumber(nin);
835 map = AIR_CALLOC(maxid+1, unsigned int)(unsigned int*)(calloc((maxid+1), sizeof(unsigned int))); /* we do need it zeroed out */
836 if (!map) {
837 biffAddf(NRRDnrrdBiffKey, "%s: couldn't allocate internal LUT", me);
838 airMopError(mop); return 1;
839 }
840 airMopAdd(mop, map, airFree, airMopAlways);
841 for (I=0; I<NN; I++) {
842 map[lup(nin->data, I)] = 1;
843 }
844 numid = 0;
845 for (jd=0; jd<=maxid; jd++) {
846 numid += map[jd];
847 }
848
849 if (nvalP) {
850 if (!(*nvalP)) {
851 *nvalP = nrrdNew();
852 }
853 if (nrrdMaybeAlloc_va(*nvalP, nin->type, 1,
854 AIR_CAST(size_t, numid)((size_t)(numid)))) {
855 biffAddf(NRRDnrrdBiffKey, "%s: couldn't allocate output value list", me);
856 airMopError(mop); return 1;
857 }
858 airMopAdd(mop, nvalP, (airMopper)airSetNull, airMopOnError);
859 airMopAdd(mop, *nvalP, (airMopper)nrrdNuke, airMopOnError);
860 }
861
862 id = 0;
863 for (jd=0; jd<=maxid; jd++) {
864 if (map[jd]) {
865 map[jd] = id;
866 if (nvalP) {
867 ins((*nvalP)->data, id, jd);
868 }
869 id++;
870 }
871 }
872 for (I=0; I<NN; I++) {
873 ins(nout->data, I, map[lup(nin->data, I)]);
874 }
875
876 if (nrrdContentSet_va(nout, func, nin, "")) {
877 biffAddf(NRRDnrrdBiffKey, "%s:", me);
878 airMopError(mop); return 1;
879 }
880 /* basic info handled by nrrdCopy */
881 airMopOkay(mop);
882 return 0;
883}