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 send a time-series (fixed index, variable ival)
12 into a dataset - inverse function to THD_extract_series.
13 ind = spatial index of desired voxel
14 = ix + jy * n1 + kz * n1*n2
15 npt = # of points in array
16 typ = type of data in array [MRI_byte, MRI_short, MRI_float, MRI_complex]
17 dar = data array
18 raw = 0 if input is to be scaled, 1 if it is not
19
20 -- RWCox - 29 Oct 1999
21 -----------------------------------------------------------------*/
22
THD_insert_series(int ind,THD_3dim_dataset * dset,int npt,int typ,void * dar,int raw)23 void THD_insert_series( int ind , THD_3dim_dataset * dset ,
24 int npt , int typ , void * dar , int raw )
25 {
26 int nv , ival , dtyp , ii ;
27 float *far=NULL , *fac=NULL ;
28
29 if( ind < 0 ||
30 npt <= 0 ||
31 dar == NULL ||
32 !ISVALID_DSET(dset) ||
33 ind >= DSET_NVOX(dset) ||
34 !DSET_IS_MALLOC(dset) ) return ;
35
36
37 nv = DSET_NVALS(dset) ; if( npt > nv ) npt = nv; /* truncate input? */
38
39 if( !DSET_LOADED(dset) ) DSET_load(dset) ; /* read from disk? */
40 if( !DSET_LOADED(dset) ) return ;
41
42 dtyp = DSET_BRICK_TYPE(dset,0) ; /* dataset array type */
43
44 /* convert input to a floating point type */
45
46 if( !raw && THD_need_brick_factor(dset) )
47 fac = dset->dblk->brick_fac ;
48 else
49 raw = 1 ;
50
51 #define FAC(q) ( (fac != NULL && fac[q] != 0.0) ? 1.0/fac[q] : 1.0 )
52
53 if( dtyp == MRI_complex ){ /* complex output ?! */
54
55 complex *car = (complex *) malloc( sizeof(complex) * npt ) ;
56 complex *bar ;
57
58 switch( typ ){
59 default:
60 free(car) ; return ; /* bad input */
61
62 case MRI_complex:
63 memcpy( car , dar , sizeof(complex)*npt ) ;
64 break ;
65
66 case MRI_float:{
67 float *a = (float *) dar ;
68 for( ii=0 ; ii < npt ; ii++ ){ car[ii].r = a[ii] ; car[ii].i = 0.0 ; }
69 }
70 break ;
71
72 case MRI_short:{
73 short *a = (short *) dar ;
74 for( ii=0 ; ii < npt ; ii++ ){ car[ii].r = a[ii] ; car[ii].i = 0.0 ; }
75 }
76 break ;
77
78 case MRI_byte:{
79 byte *a = (byte *) dar ;
80 for( ii=0 ; ii < npt ; ii++ ){ car[ii].r = a[ii] ; car[ii].i = 0.0 ; }
81 }
82 break ;
83 }
84
85 /* can now copy car into dataset, and exit */
86
87 if( !raw )
88 for( ii=0 ; ii < npt ; ii++ ){ car[ii].r *= FAC(ii) ; car[ii].i *= FAC(ii) ; }
89
90 for( ii=0 ; ii < npt ; ii++ ){
91 bar = (complex *) DSET_ARRAY(dset,ii) ;
92 bar[ind] = car[ii] ;
93 }
94
95 free(car) ; return ;
96 }
97
98 /* otherwise, compute a temporary float array */
99
100 far = (float *) malloc( sizeof(float) * npt ) ;
101 switch( typ ){
102 default:
103 free(far) ; return ; /* bad input */
104
105 case MRI_complex:{
106 complex * a = (complex *) dar ;
107 for( ii=0 ; ii < npt ; ii++ ) far[ii] = CABS(a[ii]) ;
108 }
109 break ;
110
111 case MRI_float:
112 memcpy( far , dar , sizeof(float)*npt ) ;
113 break ;
114
115 case MRI_short:{
116 short * a = (short *) dar ;
117 for( ii=0 ; ii < npt ; ii++ ) far[ii] = a[ii] ;
118 }
119 break ;
120
121 case MRI_byte:{
122 byte * a = (byte *) dar ;
123 for( ii=0 ; ii < npt ; ii++ ) far[ii] = a[ii] ;
124 }
125 break ;
126 }
127
128 if( !raw ) for( ii=0 ; ii < npt ; ii++ ) far[ii] *= FAC(ii) ;
129
130 /* now convert this to the output */
131
132 switch( dtyp ){
133
134 default:
135 free(far) ; return ; /* bad dataset? */
136
137 case MRI_float:{
138 float * bar ;
139 for( ii=0 ; ii < npt ; ii++ ){
140 bar = (float *) DSET_ARRAY(dset,ii) ;
141 bar[ind] = far[ii] ;
142 }
143 }
144 break ;
145
146 case MRI_short:{
147 short * bar ;
148 for( ii=0 ; ii < npt ; ii++ ){
149 bar = (short *) DSET_ARRAY(dset,ii) ;
150 bar[ind] = SHORTIZE(far[ii]) ;
151 }
152 }
153 break ;
154
155 case MRI_byte:{
156 byte * bar ;
157 for( ii=0 ; ii < npt ; ii++ ){
158 bar = (byte *) DSET_ARRAY(dset,ii) ;
159 bar[ind] = BYTEIZE(far[ii]) ;
160 }
161 }
162 break ;
163 }
164
165 free(far) ; return ;
166 }
167