update to grasp script (courtesy of tobias block)
Martin Uecker
8 years ago
18 | 18 | export ITER=30 |
19 | 19 | export REG=0.05 |
20 | 20 | SCALE=0.6 |
21 | LOGFILE=/dev/null | |
21 | LOGFILE=/dev/stdout | |
22 | 22 | MAXPROC=4 |
23 | 23 | MAXTHREADS=4 |
24 | 24 | |
100 | 100 | exit 1 |
101 | 101 | fi |
102 | 102 | |
103 | if [ ! -e $TOOLBOX_PATH/fft ] ; then | |
103 | if [ ! -e $TOOLBOX_PATH/bart ] ; then | |
104 | 104 | echo "\$TOOLBOX_PATH is not set correctly!" >&2 |
105 | 105 | exit 1 |
106 | 106 | fi |
116 | 116 | { |
117 | 117 | |
118 | 118 | # read TWIX file |
119 | twixread -A $input grasp | |
119 | bart twixread -A $input grasp | |
120 | 120 | |
121 | 121 | export READ=$(tail -n1 grasp.hdr | cut -f1 -d" ") |
122 | 122 | export COILS=$(tail -n1 grasp.hdr | cut -f4 -d" ") |
132 | 132 | #rm grasp.* grasp2.* |
133 | 133 | |
134 | 134 | # inverse FFT along 3rd dimension |
135 | fft -i -u $(bitmask 2) grasp grasp_hybrid | |
135 | bart fft -i -u $(bart bitmask 2) grasp grasp_hybrid | |
136 | 136 | rm grasp.cfl grasp.hdr |
137 | 137 | |
138 | 138 | SLICES=$(tail -n1 grasp_hybrid.hdr | cut -f3 -d" ") |
139 | 139 | |
140 | 140 | |
141 | 141 | # create trajectory with 400 spokes and 2x oversampling |
142 | traj -G -x$READ -y$CALIB r | |
143 | scale $SCALE r rcalib | |
142 | bart traj -G -x$READ -y$CALIB r | |
143 | bart scale $SCALE r rcalib | |
144 | 144 | |
145 | 145 | # create trajectory with 2064 spokes and 2x oversampling |
146 | traj -G -x$READ -y$(($SPOKES * $PHASES)) r | |
147 | scale $SCALE r r2 | |
146 | bart traj -G -x$READ -y$(($SPOKES * $PHASES)) r | |
147 | bart scale $SCALE r r2 | |
148 | 148 | |
149 | 149 | # split off time dimension into index 10 |
150 | reshape $(bitmask 2 10) $SPOKES $PHASES r2 rfull | |
150 | bart reshape $(bart bitmask 2 10) $SPOKES $PHASES r2 rfull | |
151 | 151 | |
152 | 152 | # number of threads per slice |
153 | 153 | export OMP_NUM_THREADS=$MAXTHREADS |
155 | 155 | calib_slice() |
156 | 156 | { |
157 | 157 | # extract slice |
158 | slice 2 $1 grasp_hybrid grasp1-$1 | |
158 | bart slice 2 $1 grasp_hybrid grasp1-$1 | |
159 | 159 | |
160 | 160 | # extract first $CALIB spokes |
161 | extract 1 $(($SKIP + 0)) $(($SKIP + $CALIB - 1)) grasp1-$1 grasp2-$1 | |
161 | bart extract 1 $(($SKIP + 0)) $(($SKIP + $CALIB - 1)) grasp1-$1 grasp2-$1 | |
162 | 162 | |
163 | 163 | # reshape dimensions |
164 | reshape $(bitmask 0 1 2 3) 1 $READ $CALIB $COILS grasp2-$1 grasp3-$1 | |
164 | bart reshape $(bart bitmask 0 1 2 3) 1 $READ $CALIB $COILS grasp2-$1 grasp3-$1 | |
165 | 165 | |
166 | 166 | # apply inverse nufft to first $CALIB spokes |
167 | nufft -i -t rcalib grasp3-$1 img-$1.coo | |
167 | bart nufft -i -t rcalib grasp3-$1 img-$1.coo | |
168 | 168 | } |
169 | 169 | |
170 | 170 | recon_slice() |
171 | 171 | { |
172 | 172 | # extract sensitivities for slice |
173 | slice 2 $1 sens sens-$1 | |
173 | bart slice 2 $1 sens sens-$1 | |
174 | 174 | |
175 | 175 | # extract spokes and split-off time dim |
176 | extract 1 $(($SKIP + 0)) $(($SKIP + $SPOKES * $PHASES - 1)) grasp1-$1 grasp2-$1 | |
177 | reshape $(bitmask 1 2) $SPOKES $PHASES grasp2-$1 grasp1-$1 | |
176 | bart extract 1 $(($SKIP + 0)) $(($SKIP + $SPOKES * $PHASES - 1)) grasp1-$1 grasp2-$1 | |
177 | bart reshape $(bart bitmask 1 2) $SPOKES $PHASES grasp2-$1 grasp1-$1 | |
178 | 178 | |
179 | 179 | # move time dimensions to dim 10 and reshape |
180 | transpose 2 10 grasp1-$1 grasp2-$1 | |
181 | reshape $(bitmask 0 1 2) 1 $READ $SPOKES grasp2-$1 grasp1-$1 | |
180 | bart transpose 2 10 grasp1-$1 grasp2-$1 | |
181 | bart reshape $(bart bitmask 0 1 2) 1 $READ $SPOKES grasp2-$1 grasp1-$1 | |
182 | 182 | rm grasp2-$1.cfl grasp2-$1.hdr |
183 | 183 | |
184 | 184 | # reconstruction with tv penality along dimension 10 |
185 | 185 | # old (v0.2.08): |
186 | 186 | # pics -S -d5 -lv -u10. -r$REG -R$(bitmask 10) -i$ITER -t rfull grasp1-$1 sens-$1 i-$1.coo |
187 | 187 | # new (v0.2.09): |
188 | pics -S -d5 -u10. -RT:$(bitmask 10):0:$REG -i$ITER -t rfull grasp1-$1 sens-$1 i-$1.coo | |
188 | bart pics -S -d5 -u10. -RT:$(bart bitmask 10):0:$REG -i$ITER -t rfull grasp1-$1 sens-$1 i-$1.coo | |
189 | 189 | |
190 | 190 | # clean up temp files |
191 | 191 | rm *-$1.cfl *-$1.hdr |
198 | 198 | seq -w 0 $(($SLICES - 1)) | xargs -I {} -P $MAXPROC bash -c "calib_slice {}" |
199 | 199 | |
200 | 200 | # transform back to k-space and compute sensitivities |
201 | join 2 img-*.coo img | |
202 | fft -u $(bitmask 0 1 2) img ksp | |
201 | bart join 2 img-*.coo img | |
202 | bart fft -u $(bart bitmask 0 1 2) img ksp | |
203 | 203 | |
204 | 204 | #ecalib -S -c0.8 -m1 -r20 ksp sens |
205 | 205 | |
206 | 206 | # transpose because we already support off-center calibration region |
207 | 207 | # in dim 0 but here we might have it in 2 |
208 | transpose 0 2 ksp ksp2 | |
209 | ecalib -S -c0.8 -m1 -r20 ksp2 sens2 | |
210 | transpose 0 2 sens2 sens | |
208 | bart transpose 0 2 ksp ksp2 | |
209 | bart ecalib -S -c0.8 -m1 -r20 ksp2 sens2 | |
210 | bart transpose 0 2 sens2 sens | |
211 | 211 | |
212 | 212 | # loop over slices |
213 | 213 | seq -w 0 $(($SLICES - 1)) | xargs -I {} -P $MAXPROC bash -c "recon_slice {}" |
214 | 214 | #echo 20 | xargs -i --max-procs=$MAXPROC bash -c "recon_slice {}" |
215 | 215 | |
216 | 216 | # join slices back together |
217 | join 2 i-*.coo $output | |
217 | bart join 2 i-*.coo $output | |
218 | 218 | |
219 | 219 | # generate dicoms |
220 | 220 | #for s in $(seq -w 0 $(($SLICES - 1))) ; do |
221 | 221 | # for p in $(seq -w 0 $(($PHASES - 1))) ; do |
222 | # | |
223 | # slice 10 $p i-$s.coo i-$p-$s.coo | |
224 | # toimg i-$p-$s.coo $output-$p$-$s.dcm | |
222 | # bart slice 10 $p i-$s.coo i-$p-$s.coo | |
223 | # bart toimg i-$p-$s.coo $output.series$p.slice$s.dcm | |
225 | 224 | # done |
226 | 225 | #done |
227 | 226 | |
228 | 227 | } > $LOGFILE |
229 | 228 | |
230 | 229 | exit 0 |
231 | ||
232 |