1 #include "mrilib.h"
2
main(int argc,char * argv[])3 int main( int argc , char * argv[] )
4 {
5 int iarg=1 , dcode=0 , maxgap=9 , nftot=0, i=0;
6 char * prefix="rowfillin" , * dstr=NULL;
7 THD_3dim_dataset * inset , * outset ;
8 MRI_IMAGE * brim ;
9 int verb=0, bin=0;
10
11 if( argc < 2 || strcmp(argv[1],"-help") == 0 ){
12 printf("Usage: 3dRowFillin [options] dataset\n"
13 "Extracts 1D rows in the given direction from a 3D dataset,\n"
14 "searches for blank (zero) regions, and fills them in if\n"
15 "the blank region isn't too large and it is flanked by\n"
16 "the same value on either edge. For example:\n"
17 " input row = 0 1 2 0 0 2 3 0 3 0 0 4 0\n"
18 " output row = 0 1 2 2 2 2 3 3 3 0 0 4 0\n"
19 "\n"
20 "OPTIONS:\n"
21 " -maxgap N = set the maximum length of a blank region that\n"
22 " will be filled in to 'N' [default=9].\n"
23 " -dir D = set the direction of fill to 'D', which can\n"
24 " be one of the following:\n"
25 " A-P, P-A, I-S, S-I, L-R, R-L, x, y, z, \n"
26 " XYZ.OR, XYZ.AND\n"
27 " The first 6 are anatomical directions;\n"
28 " x,y, and z, are reference to the dataset\n"
29 " internal axes. \n"
30 " XYZ.OR means do the fillin in x, followed by y,\n"
31 " followed by z directions.\n"
32 " XYZ.AND is like XYZ.OR but only accept voxels that\n"
33 " would have been filled in each of the three fill\n"
34 " calls. \n"
35 " Note that with XYZ* options, the fill value depends on\n"
36 " on the axis orientation. So you're better off sticking\n"
37 " to single valued dsets when using them. \n"
38 " See also -binary option below\n"
39 " -binary: Turn input dataset to 0 and 1 before filling in.\n"
40 " Output will also be a binary valued dataset.\n"
41 " -prefix P = set the prefix to 'P' for the output dataset.\n"
42 "\n"
43 "N.B.: If the input dataset has more than one sub-brick,\n"
44 " only the first one will be processed.\n"
45 "\n"
46 "* The intention of this program is to let you fill in slice gaps\n"
47 " made when drawing ROIs with the 'Draw Dataset' plugin. If you\n"
48 " draw every 5th coronal slice, say, then you could fill in using\n"
49 " 3dRowFillin -maxgap 4 -dir A-P -prefix fredfill fred+orig\n"
50 "* This program is moderately obsolescent, since I later added\n"
51 " the 'Linear Fillin' controls to the 'Draw Dataset' plugin.\n"
52 "\n"
53 ) ;
54 PRINT_COMPILE_DATE ; exit(0) ;
55 }
56
57 mainENTRY("3dRowFillin main"); machdep();
58 AFNI_logger("3dRowFillin",argc,argv);
59 PRINT_VERSION("3dRowFillin") ;
60
61 /*-- scan args --*/
62
63 while( iarg < argc && argv[iarg][0] == '-' ){
64
65 if( strncmp(argv[iarg],"-verb",5) == 0 ){
66 verb++ ; iarg++ ; continue ;
67 }
68
69 if( strncmp(argv[iarg],"-bin",4) == 0 ){
70 bin=1; iarg++ ; continue ;
71 }
72
73 if( strcmp(argv[iarg],"-prefix") == 0 ){
74 prefix = argv[++iarg] ;
75 if( !THD_filename_ok(prefix) ){
76 fprintf(stderr,"*** Illegal string after -prefix!\n"); exit(1) ;
77 }
78 iarg++ ; continue ;
79 }
80
81 if( strcmp(argv[iarg],"-maxgap") == 0 ){
82 maxgap = strtol( argv[++iarg] , NULL , 10 ) ;
83 if( maxgap < 1 ){
84 fprintf(stderr,"*** Illegal value after -maxgap!\n"); exit(1);
85 }
86 iarg++ ; continue ;
87 }
88
89 if( strcmp(argv[iarg],"-dir") == 0 ){
90 dstr = argv[++iarg] ;
91 iarg++ ; continue ;
92 }
93
94 fprintf(stderr,"*** Illegal option: %s\n",argv[iarg]) ; exit(1) ;
95 }
96
97 if( dstr == NULL ){
98 fprintf(stderr,"*** No -dir option on command line!\n"); exit(1);
99 }
100 if( iarg >= argc ){
101 fprintf(stderr,"*** No input dataset on command line!\n"); exit(1);
102 }
103
104 inset = THD_open_dataset( argv[iarg] ) ;
105 CHECK_OPEN_ERROR(inset,argv[iarg]) ;
106
107 outset = EDIT_empty_copy( inset ) ;
108 EDIT_dset_items( outset , ADN_prefix , prefix , ADN_none ) ;
109 if( THD_deathcon() && THD_is_file( DSET_HEADNAME(outset) ) ){
110 fprintf(stderr,"** Output file %s exists -- cannot overwrite!\n",
111 DSET_HEADNAME(outset) ) ;
112 exit(1) ;
113 }
114
115 tross_Copy_History( inset , outset ) ;
116 tross_Make_History( "3dRowFillin" , argc,argv , outset ) ;
117
118 if( DSET_NVALS(inset) > 1 ){
119 fprintf(stderr,"++ WARNING: input dataset has more than one sub-brick!\n");
120 EDIT_dset_items( outset ,
121 ADN_ntt , 0 ,
122 ADN_nvals , 1 ,
123 ADN_none ) ;
124 }
125
126 if( DSET_BRICK_TYPE(outset,0) != MRI_short &&
127 DSET_BRICK_TYPE(outset,0) != MRI_byte ){
128
129 fprintf(stderr,"*** This program only works on byte and short dataset!\n");
130 exit(1) ;
131 }
132
133 if (strcasecmp(dstr,"XYZ.OR") ==0) dcode = 4;
134 else if (strcasecmp(dstr,"XYZ.AND")==0) dcode = 5;
135 else {
136 switch( *dstr ){
137 case 'X': case 'x': dcode = 1 ; break ;
138 case 'Y': case 'y': dcode = 2 ; break ;
139 case 'Z': case 'z': dcode = 3 ; break ;
140 default:
141 if( *dstr == ORIENT_tinystr[outset->daxes->xxorient][0] ||
142 *dstr == ORIENT_tinystr[outset->daxes->xxorient][1] ) dcode = 1 ;
143
144 if( *dstr == ORIENT_tinystr[outset->daxes->yyorient][0] ||
145 *dstr == ORIENT_tinystr[outset->daxes->yyorient][1] ) dcode = 2 ;
146
147 if( *dstr == ORIENT_tinystr[outset->daxes->zzorient][0] ||
148 *dstr == ORIENT_tinystr[outset->daxes->zzorient][1] ) dcode = 3 ;
149 break ;
150 }
151 }
152 if( dcode == 0 ){
153 ERROR_exit("Illegal -dir direction '%s'",dstr) ;
154 }
155 if( verb )
156 INFO_message("Direction = axis %d in dataset",dcode) ;
157
158 if (bin) DSET_mallocize(inset);
159 DSET_load(inset) ; CHECK_LOAD_ERROR(inset) ;
160
161
162 /* binarize input */
163 if (bin) {
164 switch(DSET_BRICK_TYPE(inset,0)) {
165 case MRI_short: {
166 short *dp;
167 dp = DSET_ARRAY(inset, 0);
168 for (i=0; i<DSET_NVOX(inset); ++i) {
169 if (dp[i]) dp[i]=1;
170 }
171 break;
172 }
173 case MRI_byte: {
174 byte *dp;
175 dp = DSET_ARRAY(inset, 0);
176 for (i=0; i<DSET_NVOX(inset); ++i) {
177 if (dp[i]) dp[i]=1;
178 }
179 break;
180 }
181 default:
182 fprintf(stderr,"*** Illegal brick type!\n") ; exit(1) ;
183 }
184 }
185
186 if (dcode < 4) {
187 brim = mri_copy( DSET_BRICK(inset,0) ) ;
188 EDIT_substitute_brick( outset , 0 , brim->kind , mri_data_pointer(brim) ) ;
189 DSET_unload(inset) ;
190 nftot = THD_dataset_rowfillin( outset , 0 , dcode , maxgap ) ;
191 fprintf(stderr,"++ Number of voxels filled = %d\n",nftot) ;
192 } else if (dcode == 4) {
193 brim = mri_copy( DSET_BRICK(inset,0) ) ;
194 EDIT_substitute_brick( outset , 0 , brim->kind , mri_data_pointer(brim) ) ;
195 DSET_unload(inset) ;
196 nftot = 0; dcode=0;
197 for (dcode=1; dcode<4; ++dcode) {
198 nftot = nftot +
199 THD_dataset_rowfillin( outset , 0 , dcode , maxgap ) ;
200 }
201 fprintf(stderr,"++ Number of voxels OR filled = %d\n",nftot) ;
202 } else if (dcode == 5) {
203 int nf=0;
204 THD_3dim_dataset *otmp[3]={NULL, NULL, NULL};
205 MRI_IMAGE *brimtmp=NULL;
206
207 brim = mri_copy( DSET_BRICK(inset,0) ) ;
208 EDIT_substitute_brick( outset , 0 , brim->kind , mri_data_pointer(brim) ) ;
209 for (i=0; i<3; ++i) {
210 otmp[i] = EDIT_empty_copy( outset ) ;
211 brimtmp = mri_copy( DSET_BRICK(inset,0) ) ;
212 EDIT_substitute_brick( otmp[i] , 0 ,
213 brimtmp->kind, mri_data_pointer(brimtmp) ) ;
214 THD_dataset_rowfillin( otmp[i] , 0 , i+1 , maxgap);
215 }
216 DSET_unload(inset) ;
217
218 switch(DSET_BRICK_TYPE(outset,0)) {
219 case MRI_short: {
220 short * dp[3], *doo;
221 nftot = 0;
222 doo = DSET_ARRAY(outset, 0);
223 for (i=0; i<3; ++i) dp[i] = DSET_ARRAY(otmp[i], 0);
224 for (i=0; i<DSET_NVOX(outset); ++i) {
225 if ( !doo[i] && (dp[0][i] == dp[1][i]) &&
226 (dp[0][i] == dp[2][i]) && dp[0][i] ) {
227 doo[i] = dp[0][i]; ++nftot;
228 }
229 }
230 break;
231 }
232 case MRI_byte: {
233 byte * dp[3], *doo;
234 nftot = 0;
235 doo = DSET_ARRAY(outset, 0);
236 for (i=0; i<3; ++i) dp[i] = DSET_ARRAY(otmp[i], 0);
237 for (i=0; i<DSET_NVOX(outset); ++i) {
238 if ( !doo[i] && (dp[0][i] == dp[1][i]) &&
239 (dp[0][i] == dp[2][i]) && dp[0][i] ) {
240 doo[i] = dp[0][i]; ++nftot;
241 }
242 }
243 break;
244 }
245 default:
246 fprintf(stderr,"*** Illegal brick type!\n") ; exit(1) ;
247 }
248
249 for (i=0; i<3; ++i) DSET_delete(otmp[i]); otmp[i] = NULL;
250 fprintf(stderr,"++ Number of voxels AND filled = %d\n",nftot) ;
251 }
252
253 DSET_write(outset) ;
254 exit(0) ;
255 }
256