1#!/bin/tcsh
2
3@global_parse `basename $0` "$*" ; if ($status) exit 0
4
5# [PT: July 2, 2018]
6# + Edited to allow for inclusion of path when naming output JPGs
7#
8# [PT: Mar 6, 2020]
9# + note about updating PBAR_ALL, which doesn't appear necessary at the
10#   moment, just based on how outputs look currently
11# + added '-no1D' to the afni command, to remove annoying warnings if there
12#   happens to be a *.1D file lurking
13# + add in disable IMSAVE warnings with env var
14# + make the GUI run quietly (-q)
15# + set many other env vars to be quieter.
16
17# help?
18
19if ( $#argv < 2 ) then
20  echo
21  echo '-----------------------------------------------------------------'
22  echo
23  echo 'This script will make a JPEG image showing the edges of an'
24  echo 'EPI dataset overlay-ed on an anatomical dataset.  The purpose is'
25  echo 'to let the user (you) judge the quality of the 3D registration.'
26  echo
27  echo 'Three images from each of the coronal, axial, and sagittal'
28  echo 'AFNI image viewers are used, laid out in a 3x3 grid.'
29  echo
30  echo '@snapshot_volreg works by running the AFNI GUI inside a "virtual"'
31  echo 'X11 display server program named "Xvfb", and saving images from'
32  echo 'that copy of AFNI.  The script also uses programs from the netpbm11'
33  echo 'software library to put the saved images together into a pleasing'
34  echo 'layout.  If the script cannot find the netpbm11 software, it will'
35  echo 'not run :('
36  echo
37  echo '-----------------------------------------------------------------'
38  echo 'Usage: @snapshot_volreg ANATdataset EPIdataset [jname] [xdisplay]'
39  echo
40  echo 'Sample (from an afni_proc.py results directory):'
41  echo
42  echo '  @snapshot_volreg anat_final.sub-10506+tlrc      \'
43  echo '                   pb02.sub-10506.r01.volreg+tlrc sub-10506'
44  echo
45  echo 'The output file from this example is "sub-10506.jpg".'
46  echo '-----------------------------------------------------------------'
47  echo
48  echo 'Do NOT put a sub-brick index (such as "[0]") on the EPIdataset'
49  echo 'name -- the script will automatically only use the "[0]" volume.'
50  echo
51  echo '(( Although the original use was for visualizing how well EPI ))'
52  echo '(( and anatomical datasets were aligned by align_epi_anat.py, ))'
53  echo '(( it is also useful to see how well 3dQwarp aligned an       ))'
54  echo '(( anatomical dataset to a template dataset.                  ))'
55  echo
56  echo 'The optional third argument is the name of the output JPEG'
57  echo 'file -- if it does not end in ".jpg", that suffix will be added.'
58  echo 'If you do NOT supply a 3rd argument, the script will invent a name:'
59  echo 'it is probably better for you to supply a 3rd argument.'
60  echo 'It is now permitted to include an output path as part of the third'
61  echo 'argument.'
62  echo
63  echo 'The fourth (and very optional) argument is the display number'
64  echo 'of an ALREADY RUNNING copy of Xvfb, as in'
65  echo '  Xvfb :88 -screen 0 1024x768x24 &'
66  echo 'If you do NOT supply this number (88 in the example), then'
67  echo 'the script will start its own Xvfb (on a display of its choosing),'
68  echo 'use it once, and then stop it. If you are going to run this script'
69  echo 'many times in a row, starting and stopping your own Xvfb'
70  echo 'instance will speed things up a little. Normally, you do not'
71  echo 'need to use this 4th argument.'
72  echo
73  echo '-----------------------------------------------------------------'
74  echo
75  echo 'The edges from a typical EPI dataset are usually broken up and'
76  echo 'do not completely outline sulci, ventricles, etc.  In judging'
77  echo 'the quality of alignment, I usually start by looking at the'
78  echo 'outlines of the large lateral ventricles -- if those are very'
79  echo 'wrong, the alignment is not good.  After that, I look at the'
80  echo 'sulci in the superior part of the brain -- if the EPI edges'
81  echo 'there seem to be mostly aligned with the sulci, then I am'
82  echo 'usually happy.  The base of the brain, where lots of EPI'
83  echo 'dropout happens, often does not not show good edge alignment'
84  echo 'even when the rest of the brain alignment looks good.'
85  echo
86  echo '-----------------------------------------------------------------'
87  echo
88  echo 'If this script crashes, then it might leave behind files with'
89  echo 'names that start with "zzerm".  Delete these files.'
90  echo 'It is also possible that the Xvfb program will still be running'
91  echo 'if this script crashes.  A command such as that below can'
92  echo 'be used to see if you have any stray Xvfb programs running:'
93  echo
94  echo '  ps X | grep Xvfb | grep -v grep'
95  echo
96  echo 'If there are any such programs, the command below can be used'
97  echo 'to kill all of them:'
98  echo
99  echo '  killall Xvfb'
100  echo
101  echo '-------------- Author: The Madd Allineator ----------------------'
102  echo
103  exit 0
104endif
105
106# control variables
107
108setenv AFNI_DONT_LOGFILE YES
109unset noclobber
110
111# check if all the needed auxiliary programs exist
112
113set nerr = 0
114set errm = "** ERROR:"
115
116set plist = ( Xvfb djpeg cjpeg pnmcat pbmtext pnmscale pbmtopgm )
117foreach pppp ( $plist )
118  set wwww = `which $pppp`
119  if ( $status != 0 ) then
120    @ nerr++
121    set errm = "$errm $pppp"
122  endif
123end
124
125# we can use either pamcomp or pnmcomp, so only need to
126# find one of the twain
127
128set pcprog = pamcomp
129set wwww = `which $pcprog`
130if ( $status != 0 ) then
131  set pcprog = pnmcomp
132  set wwww = `which $pcprog`
133  if ( $status != 0 ) then
134    @ nerr++
135    set errm = "$errm (pnmcomp OR pamcomp)"
136  endif
137endif
138
139if ( $nerr > 0 ) then
140  echo "$errm -- not found in path -- @snapshot_volreg fails"
141  echo "** WARNING: this script cannot run without installing package netpbm11"
142  exit 1
143endif
144
145# set the prefix for the anat and epi datasets
146
147set adset = $argv[1]
148set abase = $adset:t
149set anat  = `basename $abase .gz`
150set anat  = `basename $anat  .nii`
151set anat  = `basename $anat  .HEAD`
152set anat  = `basename $anat  .BRIK`
153set anat  = `basename $anat  +tlrc`
154set anat  = `basename $anat  +orig`
155set anat  = `basename $anat  +acpc`
156set anat  = `basename $anat  +tlrc.`
157set anat  = `basename $anat  +orig.`
158set anat  = `basename $anat  +acpc.`
159if ( $abase == $anat.nii.gz ) then
160  set asuff = .nii.gz
161else if ( $abase == $anat.nii ) then
162  set asuff = .nii
163else
164  set asuff = ""
165endif
166
167set edset = $argv[2]
168set epi   = `basename $edset .gz`
169set epi   = `basename $epi   .nii`
170set epi   = `basename $epi   .HEAD`
171set epi   = `basename $epi   .BRIK`
172set epi   = `basename $epi   +orig`
173set epi   = `basename $epi   +acpc`
174set epi   = `basename $epi   +tlrc`
175set epi   = `basename $epi   +tlrc.`
176set anat  = `basename $anat  +tlrc.`
177set anat  = `basename $anat  +tlrc`
178set anat  = `basename $anat  +orig.`
179set anat  = `basename $anat  +acpc.`
180
181# set output image prefix
182
183if ( $#argv > 2 ) then
184  # [PT: July 2, 2018] Include output path with filename
185  set jdir = `dirname $argv[3]`
186
187  set jnam = $argv[3]
188  set jnam = `basename "$jnam" .jpg`
189  set jnam = `basename "$jnam" .JPG`
190  set jnam = `basename "$jnam" .jpeg`
191  set jnam = `basename "$jnam" .JPEG`
192else
193  set jdir = "."
194  set jnam = $anat.$epi
195endif
196
197# Are we re-using Xvfb?
198
199unset xdisplay
200unset killX
201if ( $#argv > 3 ) then
202  set xdisplay = $argv[4]
203  if ( ! -e /tmp/.X${xdisplay}-lock ) then
204    echo "** ERROR: it doesn't look like Xvfb is running with display $xdisplay"
205    exit 1
206  endif
207endif
208
209set exad = `3dinfo -exists $adset`
210set exep = `3dinfo -exists $edset`
211
212if ( $exad == 0 || $exep == 0 ) then
213  if ( $exad == 0 ) echo "** ERROR: @snapshot_volreg can't find $adset"
214  if ( $exep == 0 ) echo "** ERROR: @snapshot_volreg can't find $edset"
215  exit 1
216endif
217
218# set some AFNI GUI things
219
220setenv AFNI_NOSPLASH             YES
221setenv AFNI_SPLASH_MELT          NO
222setenv AFNI_LEFT_IS_LEFT         YES
223setenv AFNI_IMAGE_LABEL_MODE     5
224setenv AFNI_IMAGE_LABEL_SIZE     2
225setenv AFNI_ENVIRON_WARNINGS     NO
226setenv AFNI_COMPRESSOR           NONE
227setenv AFNI_IMSAVE_WARNINGS      NO
228setenv AFNI_NEVER_SAY_GOODBYE    YES
229setenv AFNI_STARTUP_WARNINGS     NO
230setenv AFNI_MOTD_CHECK           NO
231
232# start the X virtual frame buffer on display #xdisplay
233
234set ranval = `count -dig 1 1 999999 R1`
235
236if ( $?xdisplay == 0 ) then
237  set killX     = 1
238  set ntry      = 1
239  set Xnotfound = 1
240  while( $Xnotfound )
241    set xdisplay = `count -dig 1 3 999 R1`
242    if ( -e /tmp/.X${xdisplay}-lock ) continue
243    echo " -- trying to start Xvfb :${xdisplay}"
244    Xvfb :${xdisplay} -screen 0 1024x768x24 >& /dev/null &
245    sleep 2
246    jobs > zzerm.$ranval.txt
247    grep -q Xvfb zzerm.$ranval.txt
248    set Xnotfound = $status
249    \rm -f zzerm.$ranval.txt
250    if ( $Xnotfound == 0 )break ;
251    killall Xvfb ; sleep 1
252    @ ntry++
253    if ( $ntry > 9 ) then
254      echo "** ERROR: can't start Xvfb -- exiting"
255      exit 1
256    endif
257  end
258endif
259
260setenv DISPLAY :${xdisplay}
261
262# quasi-random temp filename prefix
263
264set zpref = zzerm.X${xdisplay}-${ranval}
265
266# crop the input anat to a standard size
267
268set cdset = ${zpref}.acrop.nii
269
2703dAutobox -npad 17 -prefix $cdset $adset
271
272# resample the EPI to the anat grid
273
2743dAllineate -input  ${edset}'[0]'     \
275            -master ${cdset}          \
276            -prefix ${zpref}.epiR.nii \
277            -1Dparam_apply IDENTITY   \
278            -final cubic
279
280set Rnam = ${zpref}.epiR.nii
281
282# aniso- or median-smooth the resampled EPI
283
284if ( 0 ) then
285  3danisosmooth -iters 2 -prefix ${zpref}.epiRS.nii ${zpref}.epiR.nii
286  set Rnam = ${zpref}.epiRS.nii
287else
288  3dMedianFilter -irad 1.01 -iter 1 -prefix ${zpref}.epiRS.nii ${zpref}.epiR.nii
289  set Rnam = ${zpref}.epiRS.nii
290endif
291
292# find edges in the EPI
293
2943dedge3 -input $Rnam -prefix ${zpref}.epiE.nii
295
296# get the EPI automask and apply it to the edgized EPI
297
2983dAutomask -q -prefix ${zpref}.epiM.nii -clfrac 0.333 -dilate 5 $Rnam
2993dcalc -a ${zpref}.epiE.nii -b ${zpref}.epiM.nii -expr 'a*b' -prefix ${zpref}.epiEM.nii
300
301# get the lower and upper thresholds for the edgized EPI
302
303set epp = `3dBrickStat -non-zero -percentile 20 60 80 ${zpref}.epiEM.nii`
304set eth = $epp[2]
305set emx = $epp[4]
306
307# run AFNI to make the 3 overlay images
308
309set anatCM = ( `3dCM $cdset` )
310
311set anatNN = `3dinfo -nijk $cdset`
312
313set astep = `ccalc -int "cbrt($anatNN)/6.111"`
314
315if ( $astep < 2 ) then
316    set astep = 2
317endif
318
319# [PT] NB: at some point, it may become necessary to adjust this line:
320##     -com "SET_PBAR_ALL +99 1 Plasma"                            \
321# ... but at present, things tend to look good as they are.  For that
322# future generation, you could perhaps use 3dBrickStat to get an
323# appropriate top_val for the PBAR_ALL
324
325echo "++ Start GUI to snapshot"
326
327afni -q -noplugins -no_detach -no1D                              \
328     -DAFNI_DEFAULT_OPACITY=9                                    \
329     -DAFNI_FUNC_ALPHA=Linear                                    \
330     -com "SWITCH_UNDERLAY ${cdset}"                             \
331     -com "SWITCH_OVERLAY ${zpref}.epiEM.nii"                    \
332     -com "SET_DICOM_XYZ $anatCM"                                \
333     -com "SET_PBAR_ALL +99 1 Plasma"                            \
334     -com "SET_FUNC_RANGE $emx"                                  \
335     -com "SET_THRESHNEW $eth *"                                 \
336     -com "SEE_OVERLAY +"                                        \
337     -com "SET_XHAIRS OFF"                                       \
338     -com "OPEN_WINDOW sagittalimage opacity=6 mont=3x1:$astep"  \
339     -com "OPEN_WINDOW axialimage opacity=6 mont=3x1:$astep"     \
340     -com "OPEN_WINDOW coronalimage opacity=6 mont=3x1:$astep"   \
341     -com "SAVE_JPEG sagittalimage ${zpref}.sag.jpg blowup=2"    \
342     -com "SAVE_JPEG coronalimage  ${zpref}.cor.jpg blowup=2"    \
343     -com "SAVE_JPEG axialimage    ${zpref}.axi.jpg blowup=2"    \
344     -com "QUITT"                                                \
345     ${cdset} ${zpref}.epiEM.nii
346
347# convert the output JPEGs to PNMs for manipulation
348
349djpeg ${zpref}.sag.jpg > ${zpref}.sag.pnm
350djpeg ${zpref}.cor.jpg > ${zpref}.cor.pnm
351djpeg ${zpref}.axi.jpg > ${zpref}.axi.pnm
352
353# cat them together, make a text label, overlay it, make output JPEG
354
355pnmcat -tb -jcenter -black ${zpref}.sag.pnm ${zpref}.axi.pnm ${zpref}.cor.pnm           > ${zpref}.pnm
356pbmtext -builtin fixed "$jnam"   | pnmcrop | pbmtopgm 1 1 | pnmscale 2                  > ${zpref}.t1.pgm
357pbmtext -builtin fixed "<- Left"           | pbmtopgm 1 1                               > ${zpref}.t2.pgm
358$pcprog -align=right -valign=bottom ${zpref}.t1.pgm ${zpref}.pnm                        > ${zpref}.t3.pnm
359$pcprog -align=right -valign=top    ${zpref}.t2.pgm ${zpref}.t3.pnm | cjpeg -quality 95 > "$jdir/$jnam.jpg"
360
361# delete the trash data
362
363\rm -f ${zpref}.*
364
365echo "@snapshot_volreg output image = $jdir/$jnam.jpg"
366
367# stop Xvfb if we started it ourselves
368
369if ( $?killX ) kill %1
370
371exit 0
372