Home My Page Projects Code Snippets Project Openings diderot
Summary Activity Tracker Tasks SCM

SCM Repository

[diderot] Annotation of /branches/vis12/src/dnorm/dnorm.c
ViewVC logotype

Annotation of /branches/vis12/src/dnorm/dnorm.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1739 - (view) (download) (as text)

1 : jhr 114 /*
2 :     Teem: Tools to process and visualize scientific data and images
3 :     Copyright (C) 2008, 2007, 2006, 2005 Gordon Kindlmann
4 :     Copyright (C) 2004, 2003, 2002, 2001, 2000, 1999, 1998 University of Utah
5 :    
6 :     This library is free software; you can redistribute it and/or
7 :     modify it under the terms of the GNU Lesser General Public License
8 :     (LGPL) as published by the Free Software Foundation; either
9 :     version 2.1 of the License, or (at your option) any later version.
10 :     The terms of redistributing and/or modifying this software also
11 :     include exceptions to the LGPL that facilitate static linking.
12 :    
13 :     This library is distributed in the hope that it will be useful,
14 :     but WITHOUT ANY WARRANTY; without even the implied warranty of
15 :     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
16 :     Lesser General Public License for more details.
17 :    
18 :     You should have received a copy of the GNU Lesser General Public License
19 :     along with this library; if not, write to Free Software Foundation, Inc.,
20 :     51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
21 :     */
22 :    
23 :     #ifdef DIDEROT
24 :     # include <teem/air.h>
25 :     # include <teem/biff.h>
26 :     # include <teem/nrrd.h>
27 :     #else
28 :     # include "../nrrd.h"
29 :     #endif
30 :    
31 :     char *dnormInfo = ("Normalizes nrrd representation for Diderot. "
32 :     "Forces information about kind and orientation into "
33 :     "a consistent form, and nixes various other fields. ");
34 :    
35 :     int
36 : glk 1234 main(int argc, const char **argv) {
37 : glk 1235 const char *me;
38 :     char *outS;
39 : jhr 114 hestOpt *hopt;
40 :     hestParm *hparm;
41 :     airArray *mop;
42 :    
43 :     char *err;
44 :     Nrrd *nin, *nout;
45 :     NrrdIoState *nio;
46 : glk 1739 int kindIn, kindOut, headerOnly, haveMM, trivialOrient;
47 : jhr 114 unsigned int kindAxis, axi, si, sj;
48 :    
49 :     me = argv[0];
50 :     mop = airMopNew();
51 :     hparm = hestParmNew();
52 :     hopt = NULL;
53 :     airMopAdd(mop, hparm, (airMopper)hestParmFree, airMopAlways);
54 : glk 125 hestOptAdd(&hopt, "h,header", NULL, airTypeInt, 0, 0, &headerOnly, NULL,
55 :     "output header of nrrd file only, not the data itself");
56 : glk 1739 hestOptAdd(&hopt, "to", NULL, airTypeInt, 0, 0, &trivialOrient, NULL,
57 :     "(*t*rivial *o*rientation) "
58 :     "even if the input nrrd comes with full orientation or "
59 :     "per-axis min-max info, ignore it and instead assert the "
60 :     "most trivial mapping between index and world space");
61 : jhr 114 hestOptAdd(&hopt, "i", "nin", airTypeOther, 1, 1, &nin, NULL,
62 :     "input image", NULL, NULL, nrrdHestNrrd);
63 :     hestOptAdd(&hopt, "o", "nout", airTypeString, 1, 1, &outS, "-",
64 :     "output filename", NULL);
65 :    
66 :     hestParseOrDie(hopt, argc-1, argv+1, hparm,
67 :     me, dnormInfo, AIR_TRUE, AIR_TRUE, AIR_TRUE);
68 :     airMopAdd(mop, hopt, (airMopper)hestOptFree, airMopAlways);
69 :     airMopAdd(mop, hopt, (airMopper)hestParseFree, airMopAlways);
70 :    
71 :     /* can't deal with block type */
72 :     if (nrrdTypeBlock == nin->type) {
73 :     fprintf(stderr, "%s: can only have scalar kinds (not %s)\n", me,
74 :     airEnumStr(nrrdType, nrrdTypeBlock));
75 :     airMopError(mop); exit(1);
76 :     }
77 :    
78 :     /* make sure all kinds are set to something */
79 :     /* see if there's a range kind, verify that there's only one */
80 : jhr 1116 /* set haveMM */
81 :     haveMM = AIR_TRUE;
82 : jhr 114 kindIn = nrrdKindUnknown;
83 :     kindAxis = 0;
84 :     for (axi=0; axi<nin->dim; axi++) {
85 :     if (nrrdKindUnknown == nin->axis[axi].kind
86 :     || nrrdKindIsDomain(nin->axis[axi].kind)) {
87 : jhr 1116 haveMM &= AIR_EXISTS(nin->axis[axi].min);
88 :     haveMM &= AIR_EXISTS(nin->axis[axi].max);
89 :     } else {
90 :     if (nrrdKindUnknown != kindIn) {
91 :     fprintf(stderr, "%s: got non-domain kind %s on axis %u, but already "
92 :     "have %s from axis %u\n", me,
93 :     airEnumStr(nrrdKind, nin->axis[axi].kind), axi,
94 :     airEnumStr(nrrdKind, kindIn), kindAxis);
95 :     airMopError(mop); exit(1);
96 :     }
97 :     kindIn = nin->axis[axi].kind;
98 :     kindAxis = axi;
99 : jhr 114 }
100 :     }
101 :     /* see if the non-domain kind is something we can interpret as a tensor */
102 :     if (nrrdKindUnknown != kindIn) {
103 :     switch (kindIn) {
104 :     /* ======= THESE are the kinds that we can possibly output ======= */
105 :     case nrrdKind2Vector:
106 :     case nrrdKind3Vector:
107 :     case nrrdKind4Vector:
108 :     case nrrdKind2DSymMatrix:
109 :     case nrrdKind2DMatrix:
110 :     case nrrdKind3DSymMatrix:
111 :     case nrrdKind3DMatrix:
112 :     /* =============================================================== */
113 :     kindOut = kindIn;
114 :     break;
115 :     case nrrdKind3Color:
116 :     case nrrdKindRGBColor:
117 :     kindOut = nrrdKind3Vector;
118 :     break;
119 :     case nrrdKind4Color:
120 :     case nrrdKindRGBAColor:
121 :     kindOut = nrrdKind4Vector;
122 :     break;
123 :     default:
124 :     fprintf(stderr, "%s: got non-conforming kind %s on axis %u\n", me,
125 :     airEnumStr(nrrdKind, kindIn), kindAxis);
126 :     airMopError(mop); exit(1);
127 :     break;
128 :     }
129 :     } else {
130 :     kindOut = nrrdKindUnknown;
131 :     }
132 :    
133 : glk 147 /* initialize output by copying */
134 : jhr 114 nout = nrrdNew();
135 :     airMopAdd(mop, nout, (airMopper)nrrdNuke, airMopAlways);
136 :     if (nrrdCopy(nout, nin)) {
137 :     airMopAdd(mop, err = biffGet(NRRD), airFree, airMopAlways);
138 :     fprintf(stderr, "%s: trouble copying:\n%s", me, err);
139 :     airMopError(mop); exit(1);
140 :     }
141 :    
142 :     /* no comments, either advertising the format URL or anything else */
143 :     nio = nrrdIoStateNew();
144 :     airMopAdd(mop, nio, (airMopper)nrrdIoStateNix, airMopAlways);
145 :     nio->skipFormatURL = AIR_TRUE;
146 : glk 125 if (headerOnly) {
147 :     nio->skipData = AIR_TRUE;
148 :     }
149 : jhr 114 nrrdCommentClear(nout);
150 :    
151 :     /* no measurement frame */
152 :     for (si=0; si<NRRD_SPACE_DIM_MAX; si++) {
153 :     for (sj=0; sj<NRRD_SPACE_DIM_MAX; sj++) {
154 :     nout->measurementFrame[si][sj] = AIR_NAN;
155 :     }
156 :     }
157 :    
158 :     /* no key/value pairs */
159 :     nrrdKeyValueClear(nout);
160 :    
161 :     /* no content field */
162 :     nout->content = airFree(nout->content);
163 :    
164 :     /* normalize domain kinds to "space" */
165 :     /* turn off centers (perhaps Diderot should assume cell-centered) */
166 :     /* turn off thickness */
167 :     /* turn off labels and units */
168 :     for (axi=0; axi<nout->dim; axi++) {
169 :     if (nrrdKindUnknown == kindOut) {
170 :     nout->axis[axi].kind = nrrdKindSpace;
171 :     } else {
172 :     nout->axis[axi].kind = (kindAxis == axi
173 :     ? kindOut
174 :     : nrrdKindSpace);
175 :     }
176 :     nout->axis[axi].center = nrrdCenterUnknown;
177 :     nout->axis[axi].thickness = AIR_NAN;
178 :     nout->axis[axi].label = airFree(nout->axis[axi].label);
179 :     nout->axis[axi].units = airFree(nout->axis[axi].units);
180 :     nout->axis[axi].min = AIR_NAN;
181 :     nout->axis[axi].max = AIR_NAN;
182 :     nout->axis[axi].spacing = AIR_NAN;
183 :     }
184 :    
185 :     /* logic of orientation definition:
186 :     if space dimension is known:
187 :     set origin to zero if not already set
188 :     set space direction to unit vector if not already set
189 :     else:
190 :     set origin to zero and all space directions to units
191 :     might be nice to use gage's logic for mapping from world to index,
192 :     but we have to accept a greater variety of kinds and dimensions
193 :     than gage ever has to process.
194 :     */
195 : glk 1739 if (nout->spaceDim && !trivialOrient) {
196 : jhr 114 int saxi = 0;
197 : glk 147 /* we use only the space dimension, not any named space */
198 :     nout->space = nrrdSpaceUnknown;
199 : jhr 114 if (!nrrdSpaceVecExists(nout->spaceDim, nout->spaceOrigin)) {
200 :     nrrdSpaceVecSetZero(nout->spaceOrigin);
201 :     }
202 :     for (axi=0; axi<nout->dim; axi++) {
203 :     if (nrrdKindUnknown == kindOut || kindAxis != axi) {
204 :     if (!nrrdSpaceVecExists(nout->spaceDim,
205 :     nout->axis[axi].spaceDirection)) {
206 :     nrrdSpaceVecSetZero(nout->axis[axi].spaceDirection);
207 :     nout->axis[axi].spaceDirection[saxi] = 1.0;
208 :     }
209 :     saxi++;
210 :     } else {
211 :     nrrdSpaceVecSetNaN(nout->axis[axi].spaceDirection);
212 :     }
213 :     }
214 : glk 1739 } else if (haveMM && !trivialOrient) {
215 : jhr 1116 int saxi = 0;
216 :     for (axi=0; axi<nout->dim; axi++) {
217 :     if (nrrdKindUnknown == kindOut || kindAxis != axi) {
218 :     nrrdSpaceVecSetZero(nout->axis[axi].spaceDirection);
219 :     nout->axis[axi].spaceDirection[saxi]
220 :     = (nin->axis[axi].max - nin->axis[axi].min)/(nin->axis[axi].size-1);
221 :     nout->spaceOrigin[saxi] = nin->axis[axi].min;
222 :     saxi++;
223 :     } else {
224 :     nrrdSpaceVecSetNaN(nout->axis[axi].spaceDirection);
225 :     }
226 :     }
227 :     nout->spaceDim = saxi;
228 : jhr 114 } else {
229 :     int saxi = 0;
230 :     nrrdSpaceVecSetZero(nout->spaceOrigin);
231 :     for (axi=0; axi<nout->dim; axi++) {
232 :     if (nrrdKindUnknown == kindOut || kindAxis != axi) {
233 :     nrrdSpaceVecSetZero(nout->axis[axi].spaceDirection);
234 : glk 125 nout->axis[axi].spaceDirection[saxi]
235 :     = (AIR_EXISTS(nin->axis[axi].spacing)
236 :     ? nin->axis[axi].spacing
237 :     : 1.0);
238 : jhr 114 saxi++;
239 :     } else {
240 :     nrrdSpaceVecSetNaN(nout->axis[axi].spaceDirection);
241 :     }
242 :     }
243 :     nout->spaceDim = saxi;
244 :     }
245 : jhr 1116
246 :     /* space dimension has to match the number of domain axes */
247 : jhr 114 if (nout->dim != nout->spaceDim + !!kindOut) {
248 :     fprintf(stderr, "%s: output dim %d != spaceDim %d + %d %s%s%s\n",
249 :     me, nout->dim, nout->spaceDim, !!kindOut,
250 :     kindOut ? "for non-scalar (" : "(scalar data)",
251 :     kindOut ? airEnumStr(nrrdKind, kindOut) : "",
252 :     kindOut ? ") data" : "");
253 :     airMopError(mop); exit(1);
254 :     }
255 :    
256 :     if (nrrdSave(outS, nout, nio)) {
257 :     airMopAdd(mop, err = biffGet(NRRD), airFree, airMopAlways);
258 :     fprintf(stderr, "%s: trouble saving \"%s\":\n%s",
259 :     me, outS, err);
260 :     airMopError(mop); exit(1);
261 :     }
262 :    
263 :     airMopOkay(mop);
264 :     exit(0);
265 :     }

root@smlnj-gforge.cs.uchicago.edu
ViewVC Help
Powered by ViewVC 1.0.0