1 #include "mrilib.h"
2
3 /*-------------------------------------------------------------------
4 Fill in the blanks (zeros) in a row of shorts that are at most
5 maxgap in length, and are bounded by the same value.
6 Return value is number of filled-in blanks.
7 ---------------------------------------------------------------------*/
8
THD_rowfillin_short(int nrow,short * row,int maxgap)9 static int THD_rowfillin_short( int nrow , short * row , int maxgap )
10 {
11 int ii , nfill=0 , jj ;
12 short vbot ;
13
14 /*-- skip zeros at start --*/
15
16 for( ii=0 ; ii < nrow && row[ii] == 0 ; ii++ ) ; /* nada */
17 if( ii == nrow ) return nfill ; /*** was all zeros ***/
18
19 /*-- row[ii] is not zero here; now find a gap at or above this --*/
20
21 #if 0
22 { int kk ; fprintf(stderr,"maxgap=%d: row in:",maxgap);
23 for( kk=0 ; kk < nrow ; kk++ ) fprintf(stderr," %d",row[kk]) ;
24 fprintf(stderr,"\n") ; }
25 #endif
26
27 while( 1 ){
28 /*-- find next zero value (start of blank region) --*/
29
30 for( ; ii < nrow && row[ii] != 0 ; ii++ ) ; /* nada */
31 if( ii == nrow ) return nfill ; /*** didn't find any zero ***/
32
33 #if 0
34 fprintf(stderr," gap start at %d\n",ii) ;
35 #endif
36
37 vbot = row[ii-1] ; /* value just before gap */
38
39 /*-- find next nonzero value above this zero --*/
40
41 for( jj=ii+1 ; jj < nrow && row[jj] == 0 ; jj++ ) ; /* nada */
42 if( jj == nrow ) return nfill ; /*** was all zeros to the end ***/
43
44 #if 0
45 fprintf(stderr," gap end at %d: vbot=%d row[jj]=%d\n",jj,vbot,row[jj]) ;
46 #endif
47
48 if( row[jj] == vbot && jj-ii <= maxgap ){ /*** fill in this gap!!! ***/
49 nfill += (jj-ii) ;
50 for( ; ii < jj ; ii++ ) row[ii] = vbot ;
51 #if 0
52 { int kk ; fprintf(stderr,"filled:");
53 for( kk=0 ; kk < nrow ; kk++ ) fprintf(stderr," %d",row[kk]) ;
54 fprintf(stderr,"\n") ; }
55 #endif
56 }
57
58 ii = jj ; /* loop back and look for another gap */
59
60 } /* endless loop */
61 }
62
63 /*------------------------------------------------------------------------*/
64
THD_rowfillin_byte(int nrow,byte * row,int maxgap)65 static int THD_rowfillin_byte( int nrow , byte * row , int maxgap )
66 {
67 int ii , nfill=0 , jj ;
68 byte vbot ;
69
70 /*-- skip zeros --*/
71
72 for( ii=0 ; ii < nrow && row[ii] == 0 ; ii++ ) ; /* nada */
73 if( ii == nrow ) return nfill ; /*** was all zeros ***/
74
75 /*-- row[ii] is not zero here; now find a gap at or above this --*/
76
77 while( 1 ){
78 /*-- find next zero value --*/
79
80 for( ; ii < nrow && row[ii] != 0 ; ii++ ) ; /* nada */
81 if( ii == nrow ) return nfill ; /*** didn't find any zero ***/
82
83 vbot = row[ii-1] ; /* value at start of gap */
84
85 /*-- find next nonzero value above this zero --*/
86
87 for( jj=ii+1 ; jj < nrow && row[jj] == 0 ; jj++ ) ; /* nada */
88 if( jj == nrow ) return nfill ; /*** was all zeros to the end ***/
89
90 if( row[jj] == vbot && jj-ii <= maxgap ){ /*** fill in this gap!!! ***/
91 nfill += (jj-ii) ;
92 for( ; ii < jj ; ii++ ) row[ii] = vbot ;
93 }
94
95 ii = jj ; /* loop back and look for another gap */
96
97 } /* endless loop */
98 }
99
100 /*------------------------------------------------------------------------*/
101
THD_rowfillin_float(int nrow,float * row,int maxgap)102 static int THD_rowfillin_float( int nrow , float * row , int maxgap )
103 {
104 int ii , nfill=0 , jj ;
105 float vbot ;
106
107 /*-- skip zeros --*/
108
109 for( ii=0 ; ii < nrow && row[ii] == 0 ; ii++ ) ; /* nada */
110 if( ii == nrow ) return nfill ; /*** was all zeros ***/
111
112 /*-- row[ii] is not zero here; now find a gap at or above this --*/
113
114 while( 1 ){
115 /*-- find next zero value --*/
116
117 for( ; ii < nrow && row[ii] != 0 ; ii++ ) ; /* nada */
118 if( ii == nrow ) return nfill ; /*** didn't find any zero ***/
119
120 vbot = row[ii-1] ; /* value at start of gap */
121
122 /*-- find next nonzero value above this zero --*/
123
124 for( jj=ii+1 ; jj < nrow && row[jj] == 0 ; jj++ ) ; /* nada */
125 if( jj == nrow ) return nfill ; /*** was all zeros to the end ***/
126
127 if( row[jj] == vbot && jj-ii <= maxgap ){ /*** fill in this gap!!! ***/
128 nfill += (jj-ii) ;
129 for( ; ii < jj ; ii++ ) row[ii] = vbot ;
130 }
131
132 ii = jj ; /* loop back and look for another gap */
133
134 } /* endless loop */
135 }
136
137 /*---------------------------------------------------------------------------
138 Do the fillin thing on each 1D row from a dataset sub-brick.
139 Return value is number of voxels filled in (-1 if an error transpired).
140 -----------------------------------------------------------------------------*/
141
THD_dataset_rowfillin(THD_3dim_dataset * dset,int ival,int dcode,int maxgap)142 int THD_dataset_rowfillin( THD_3dim_dataset * dset, int ival, int dcode, int maxgap )
143 {
144 int kind , xx,yy,zz , nrow , nx,ny,nz ;
145 int xtop,ytop,ztop , nff , nfftot=0 ;
146
147 ENTRY("THD_dataset_rowfillin") ;
148
149 if( !ISVALID_DSET(dset) ||
150 ival < 0 ||
151 ival >= DSET_NVALS(dset) ||
152 maxgap < 1 ) RETURN(-1) ; /* bad things */
153
154 kind = DSET_BRICK_TYPE(dset,ival) ;
155 if( kind != MRI_short && kind != MRI_byte && kind != MRI_float ) RETURN(-1) ; /* bad */
156
157 nrow = THD_get_dset_rowcount( dset , dcode ) ;
158 if( nrow < 1 ) RETURN(-1) ; /* bad */
159
160 nx = DSET_NX(dset) ;
161 ny = DSET_NY(dset) ;
162 nz = DSET_NZ(dset) ;
163
164 xtop = ytop = ztop = 1 ;
165 switch( dcode ){
166 case 1: case -1: ytop=ny ; ztop=nz ; break ;
167 case 2: case -2: xtop=nx ; ztop=nz ; break ;
168 case 3: case -3: xtop=nx ; ytop=ny ; break ;
169 }
170
171 switch( kind ){
172
173 case MRI_short:{
174 short * row ;
175 for( zz=0 ; zz < ztop ; zz++ )
176 for( yy=0 ; yy < ytop ; yy++ )
177 for( xx=0 ; xx < xtop ; xx++ ){
178 row = THD_get_dset_row( dset,ival , dcode,xx,yy,zz ) ;
179 nff = THD_rowfillin_short( nrow , row , maxgap ) ;
180 if( nff > 0 ){
181 THD_put_dset_row( dset,ival , dcode,xx,yy,zz , row ) ;
182 nfftot += nff ;
183 }
184 free(row) ;
185 }
186 }
187 break ;
188
189 case MRI_byte:{
190 byte * row ;
191 for( zz=0 ; zz < ztop ; zz++ )
192 for( yy=0 ; yy < ytop ; yy++ )
193 for( xx=0 ; xx < xtop ; xx++ ){
194 row = THD_get_dset_row( dset,ival , dcode,xx,yy,zz ) ;
195 nff = THD_rowfillin_byte( nrow , row , maxgap ) ;
196 if( nff > 0 ){
197 THD_put_dset_row( dset,ival , dcode,xx,yy,zz , row ) ;
198 nfftot += nff ;
199 }
200 free(row) ;
201 }
202 }
203 break ;
204
205 case MRI_float:{
206 float * row ;
207 for( zz=0 ; zz < ztop ; zz++ )
208 for( yy=0 ; yy < ytop ; yy++ )
209 for( xx=0 ; xx < xtop ; xx++ ){
210 row = THD_get_dset_row( dset,ival , dcode,xx,yy,zz ) ;
211 nff = THD_rowfillin_float( nrow , row , maxgap ) ;
212 if( nff > 0 ){
213 THD_put_dset_row( dset,ival , dcode,xx,yy,zz , row ) ;
214 nfftot += nff ;
215 }
216 free(row) ;
217 }
218 }
219 break ;
220
221 }
222
223 RETURN(nfftot) ;
224 }
225