1 /*****************************************************************************
2    Major portions of this software are copyrighted by the Medical College
3    of Wisconsin, 1994-2000, and are released under the Gnu General Public
4    License, Version 2.  See the file README.Copyright for details.
5 ******************************************************************************/
6 
7 #include "mrilib.h"
8 #include "thd.h"
9 
10 /*---------------------------------------------------------------
11   Routine to extract a time-series (fixed index, variable ival)
12     from a previously set up FD_brick structure.
13   ixyz = spatial index of desired voxel (in brick coordinates)
14        = ix + jy * n1 + kz * n1*n2
15   Return value is an image of the type of the dataset (assumed
16   uniform).
17 -----------------------------------------------------------------*/
18 
FD_brick_to_series(int ixyz,FD_brick * br)19 MRI_IMAGE * FD_brick_to_series( int ixyz , FD_brick *br )
20 {
21    MRI_IMAGE *im ;  /* output */
22    int nv , ival ;
23    char *iar ;      /* brick in the input */
24    MRI_TYPE typ ;
25    int ix,jy,kz , ind ;
26    THD_ivec3 ind_fd , ind_ds ;
27 
28 ENTRY("FD_brick_to_series") ;
29 
30    if( ixyz < 0 || ixyz >= br->n1 * br->n2 * br->n3 ) RETURN(NULL) ;
31 
32    /** otherwise, get ready for a real image **/
33 
34    ix  = ixyz % br->n1 ;
35    jy  = ( ixyz % (br->n1 * br->n2) ) / br->n1 ;
36    kz  = ixyz / (br->n1 * br->n2) ;
37    LOAD_IVEC3( ind_fd , ix,jy,kz ) ; ind_ds = THD_fdind_to_3dind( br , ind_fd ) ;
38    ix  = ind_ds.ijk[0] ;
39    jy  = ind_ds.ijk[1] ;
40    kz  = ind_ds.ijk[2] ;
41    ind = (kz * br->dset->daxes->nyy + jy) * br->dset->daxes->nxx + ix ;
42 
43    nv = br->dset->dblk->nvals ;
44 
45    iar = DSET_ARRAY(br->dset,0) ;
46    if( iar == NULL ){  /* if data needs to be loaded from disk */
47       (void) THD_load_datablock( br->dset->dblk ) ;
48       iar = DSET_ARRAY(br->dset,0) ;
49       if( iar == NULL ) RETURN(NULL) ;
50    }
51 
52    /* 15 Sep 2004: allow for nonconstant datum */
53 
54    if( !DSET_datum_constant(br->dset) ){  /* only for stupid users */
55      float *ar ;
56      im = mri_new( nv , 1 , MRI_float ) ; ar = MRI_FLOAT_PTR(im) ;
57      for( ival = 0 ; ival < nv ; ival++ )
58        ar[ival] = THD_get_voxel( br->dset , ind , ival ) ;
59      goto image_done ;
60    }
61 
62    /* the older (more efficient) way */
63 
64    typ = DSET_BRICK_TYPE(br->dset,0) ;
65    im  = mri_new( nv , 1 , typ ) ;
66 #if 0
67    mri_zero_image(im) ;             /* 18 Oct 2001 */
68 #endif
69 
70    switch( typ ){
71 
72       default:             /* don't know what to do --> return nada */
73          mri_free( im ) ;
74          RETURN(NULL) ;
75 
76       case MRI_byte:{
77          byte *ar  = MRI_BYTE_PTR(im) , *bar ;
78          for( ival=0 ; ival < nv ; ival++ ){
79             bar = (byte *) DSET_ARRAY(br->dset,ival) ;
80             if( bar != NULL ) ar[ival] = bar[ind] ;
81          }
82       }
83       break ;
84 
85       case MRI_short:{
86          short *ar  = MRI_SHORT_PTR(im) , *bar ;
87          for( ival=0 ; ival < nv ; ival++ ){
88             bar = (short *) DSET_ARRAY(br->dset,ival) ;
89             if( bar != NULL ) ar[ival] = bar[ind] ;
90          }
91       }
92       break ;
93 
94       case MRI_float:{
95          float *ar  = MRI_FLOAT_PTR(im) , *bar ;
96          for( ival=0 ; ival < nv ; ival++ ){
97             bar = (float *) DSET_ARRAY(br->dset,ival) ;
98             if( bar != NULL ) ar[ival] = bar[ind] ;
99          }
100       }
101       break ;
102 
103       case MRI_int:{
104          int *ar  = MRI_INT_PTR(im) , *bar ;
105          for( ival=0 ; ival < nv ; ival++ ){
106             bar = (int *) DSET_ARRAY(br->dset,ival) ;
107             if( bar != NULL ) ar[ival] = bar[ind] ;
108          }
109       }
110       break ;
111 
112       case MRI_double:{
113          double *ar  = MRI_DOUBLE_PTR(im) , *bar ;
114          for( ival=0 ; ival < nv ; ival++ ){
115             bar = (double *) DSET_ARRAY(br->dset,ival) ;
116             if( bar != NULL ) ar[ival] = bar[ind] ;
117          }
118       }
119       break ;
120 
121       case MRI_complex:{
122          complex *ar  = MRI_COMPLEX_PTR(im) , *bar ;
123          for( ival=0 ; ival < nv ; ival++ ){
124             bar = (complex *) DSET_ARRAY(br->dset,ival) ;
125             if( bar != NULL ) ar[ival] = bar[ind] ;
126          }
127       }
128       break ;
129 
130       /* 15 Apr 2002: RGB types */
131 
132       case MRI_rgb:{
133          rgbyte *ar  = (rgbyte *) MRI_RGB_PTR(im) , *bar ;
134          for( ival=0 ; ival < nv ; ival++ ){
135             bar = (rgbyte *) DSET_ARRAY(br->dset,ival) ;
136             if( bar != NULL ) ar[ival] = bar[ind] ;
137          }
138       }
139       break ;
140 
141       case MRI_rgba:{
142          rgba *ar  = (rgba *) MRI_RGBA_PTR(im) , *bar ;
143          for( ival=0 ; ival < nv ; ival++ ){
144             bar = (rgba *) DSET_ARRAY(br->dset,ival) ;
145             if( bar != NULL ) ar[ival] = bar[ind] ;
146          }
147       }
148       break ;
149 
150    }
151 
152    if( THD_need_brick_factor(br->dset) ){
153       MRI_IMAGE *qim ;
154       qim = mri_mult_to_float( br->dset->dblk->brick_fac , im ) ;
155       mri_free(im) ; im = qim ;
156    }
157 
158    /* at this point, the image is ready to ship out;
159       but first, maybe attach a time origin and spacing */
160 
161 image_done:
162    if( br->dset->taxis != NULL ){  /* 21 Oct 1996 */
163       float zz , tt ;
164 
165       zz = br->dset->daxes->zzorg + kz * br->dset->daxes->zzdel ;
166       tt = THD_timeof( 0 , zz , br->dset->taxis ) ;
167 
168       im->xo = tt ; im->dx = br->dset->taxis->ttdel ;   /* origin and delta */
169 
170       if( br->dset->taxis->units_type == UNITS_MSEC_TYPE ){ /* convert to sec */
171          im->xo *= 0.001 ; im->dx *= 0.001 ;
172       }
173    } else {
174       im->xo = 0.0 ; im->dx = 1.0 ;  /* 08 Nov 1996 */
175    }
176 
177    RETURN(im) ;
178 }
179