KUL_dcm2bids.sh 47 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840841842843844845846847848849850851852853854855856857858859860861862863864865866867868869870871872873874875876877878879880881882883884885886887888889890891892893894895896897898899900901902903904905906907908909910911912913914915916917918919920921922923924925926927928929930931932933934935936937938939940941942943944945946947948949950951952953954955956957958959960961962963964965966967968969970971972973974975976977978979980981982983984985986987988989990991992993994995996997998999100010011002100310041005100610071008100910101011101210131014101510161017101810191020102110221023102410251026102710281029103010311032103310341035103610371038103910401041104210431044104510461047104810491050105110521053105410551056105710581059106010611062106310641065106610671068106910701071107210731074107510761077107810791080108110821083108410851086108710881089109010911092109310941095109610971098109911001101110211031104110511061107110811091110111111121113111411151116111711181119112011211122112311241125112611271128112911301131113211331134113511361137113811391140114111421143114411451146114711481149115011511152115311541155115611571158115911601161116211631164116511661167116811691170117111721173117411751176117711781179118011811182118311841185118611871188118911901191119211931194119511961197119811991200120112021203120412051206120712081209121012111212121312141215121612171218121912201221122212231224122512261227122812291230123112321233123412351236123712381239124012411242124312441245124612471248124912501251125212531254125512561257125812591260126112621263126412651266126712681269127012711272127312741275127612771278127912801281128212831284128512861287128812891290129112921293129412951296129712981299130013011302130313041305130613071308130913101311131213131314131513161317131813191320132113221323132413251326132713281329133013311332133313341335133613371338133913401341134213431344134513461347134813491350135113521353135413551356135713581359136013611362136313641365136613671368136913701371137213731374137513761377137813791380138113821383138413851386138713881389139013911392139313941395139613971398139914001401140214031404140514061407140814091410141114121413141414151416141714181419142014211422142314241425142614271428142914301431143214331434143514361437143814391440144114421443144414451446144714481449145014511452145314541455145614571458145914601461146214631464146514661467146814691470147114721473147414751476147714781479148014811482148314841485148614871488148914901491149214931494149514961497149814991500150115021503150415051506150715081509151015111512151315141515
  1. #!/bin/bash -e
  2. # set -x
  3. # Bash shell script to convert dicoms to bids format
  4. #
  5. # Requires dcm2bids, dcm2niix, Mrtrix3
  6. #
  7. # @ Stefan Sunaert - UZ/KUL - stefan.sunaert@uzleuven.be
  8. # @ Ahmed Radwan - KUL - ahmed.radwan@kuleuven.be
  9. #
  10. v="v0.8 - dd 08/09/2020"
  11. # Notes
  12. # - NOW USES https://github.com/UNFmontreal/Dcm2Bids
  13. # - works for GE/Siemens/Philips
  14. # - wrap around for multiple subjects: use KUL_multisubjects_dcm2bids
  15. # ----------------------------------- MAIN ---------------------------------------------
  16. # this script defines a few functions:
  17. # - Usage (for information to the novice user)
  18. # - kul_e2cl from KUL_main_functions (for logging)
  19. # - kul_dcmtags (for reading specific parameters from dicom header)
  20. # source general functions
  21. kul_main_dir=`dirname "$0"`
  22. script=`basename "$0"`
  23. source "${kul_main_dir}/KUL_main_functions.sh"
  24. cwd=$(pwd)
  25. # BEGIN LOCAL FUNCTIONS --------------
  26. # --- function Usage ---
  27. function Usage {
  28. cat <<USAGE
  29. `basename $0` converts dicoms to bids format
  30. Usage:
  31. `basename $0` -d dicom_dir -p subject (participant) -c config_file -o bids_dir <OPT_ARGS>
  32. Depends on a config file that defines parameters with sequence information, e.g.
  33. For Philips dicom we need to manually specify the mb and pe_dir:
  34. # Identifier,search-string,fmritask/inteded_for,mb,pe_dir,acq_label
  35. # Structural scans
  36. T1w,T1_PRE
  37. cT1w,T1_Post
  38. FLAIR,3D_FLAIR
  39. T2w,3D_T2
  40. SWI,SWI
  41. MTI,mtc_2dyn
  42. # functional scans
  43. func,rsfMRI_MB6,rest,6,j,singleTE
  44. sbref,rsfMRI_SBREF,rest,1,j,singleTE
  45. func,MB_mTE,rest,4,j,multiTE
  46. sbref,mTE_SBREF,rest,4,j,multiTE
  47. func,MB2_hand,HAND,2,j
  48. func,MB2_lip,LIP,2,j
  49. func,MB2_nback,nback,2,j
  50. sbref,MB2_SBREF,HAND nback,1,j
  51. # fmap: 'task' is now 'IntendedFor' (the func images/tasks that the B0_map is used for SDC)
  52. fmap,B0_map,[HAND LIP nback]
  53. # dMRI
  54. dwi,p1_b1200,-,3,j-,b1200
  55. dwi,p2_b0,-,3,j-,b0
  56. dwi,p3_b2500,-,3,j-,b2500
  57. dwi,p4_b2500,-,3,j,rev
  58. # ASL support is very limited (not BIDS compliant for now)
  59. # Indentifier, search-string
  60. ASL,pCASL
  61. explains that the T1w scan should be found by the search string "T1_PRE"
  62. func by rsfMRI, and has multiband_factor 6, and pe_dir = j
  63. the sbref will be used for the tb_fMRI for both tasks (hands and nback)
  64. For Siemens and GE dicom it can be as simple as:
  65. # Identifier,search-string,fmritask/inteded_for,mb,pe_dir,acq_label
  66. # Structural scans
  67. T1w,T1_PRE
  68. cT1w,T1_Post
  69. FLAIR,3D_FLAIR
  70. T2w,3D_T2
  71. # functional scans
  72. func,rsfMRI_MB6,rest,-,-,singleTE
  73. sbref,rsfMRI_SBREF,rest,-,-,singleTE
  74. func,MB_mTE,rest,-,-,multiTE
  75. sbref,mTE_SBREF,rest,-,-,multiTE
  76. func,MB2_hand,HAND,-,-
  77. func,MB2_lip,LIP,-,-
  78. func,MB2_nback,nback,-,-
  79. sbref,MB2_SBREF,HAND nback,-,-
  80. # dMRI
  81. dwi,p1_b1200,-,-,-,b1200
  82. dwi,p2_b0,-,-,-,b0
  83. dwi,p3_b2500,-,-,-,b2500
  84. dwi,p4_b2500,-,-,-,rev
  85. # fmap, SWI, MTC and ASL have not been tested on Siemens or GE
  86. Example:
  87. `basename $0` -p pat001 -d pat001.zip -c definitions_of_sequences.txt -o BIDS
  88. Required arguments:
  89. -d: dicom_zip_file (the zip or tar.gz containing all your dicoms, or directory containing dicoms)
  90. -p: participant (anonymised name of the subject in bids convention)
  91. -c: definitions of sequences (T1w=MPRAGE,dwi=seq, etc..., see above)
  92. Optional arguments:
  93. -o: bids directory
  94. -s: session (for longitudinal study with multiple timepoints)
  95. -t: temporary directory (default = /tmp)
  96. -e: copy task-*_events.tsv from config to BIDS dir
  97. -a: further anonymise the subject by using pydeface (takes much longer)
  98. -v: verbose
  99. USAGE
  100. exit 1
  101. }
  102. # check if jsontool is installed and install it if not
  103. if [[ $(which jsontool) ]]; then
  104. echo " jsontool already installed, good" $log
  105. else
  106. echo " jsontool not installed, installing it with pip using pip install jsontool" $log
  107. pip install jsontool
  108. fi
  109. # check if pydeface is installed and install it if not
  110. if [[ $(which pydeface) ]]; then
  111. echo " pydeface already installed, good" $log
  112. else
  113. echo " pydeface not installed, installing it with pip using pip install jsontool" $log
  114. pip install pydeface
  115. fi
  116. # check if the correct dcm2bids is installed and install it if not
  117. if [[ $(which dcm2bids_scaffold) ]]; then
  118. echo " dcm2bids already installed, good" $log
  119. else
  120. echo " dcm2bids not installed, installing it with pip using pip install dcm2bids" $log
  121. pip uninstall Dcm2Bids
  122. pip install dcm2bids
  123. fi
  124. # --- function kul_dcmtags (for reading specific parameters from dicom header & calculating missing BIDS parameters) ---
  125. function kul_dcmtags {
  126. local dcm_file=$1
  127. local out=$final_dcm_tags_file
  128. # 0/ Read out standard tags for logging
  129. local seriesdescr=$(dcminfo "$dcm_file" -tag 0008 103E | cut -c 13-)
  130. local manufacturer=$(dcminfo "$dcm_file" -tag 0008 0070 | cut -c 13-)
  131. local software=$(dcminfo "$dcm_file" -tag 0018 1020 | cut -c 13-)
  132. local imagetype=$(dcminfo "$dcm_file" -tag 0008 0008 2>/dev/null | cut -c 13- | head -n 1)
  133. local patid=$(dcminfo "$dcm_file" -tag 0010 0020 | cut -c 13-)
  134. local pixelspacing=$(dcminfo "$dcm_file" -tag 0028 0030 2>/dev/null | cut -c 13- | head -n 1)
  135. local slicethickness=$(dcminfo "$dcm_file" -tag 0018 0088 2>/dev/null | cut -c 13- | head -n 1)
  136. local acquisitionMatrix=$(dcminfo "$dcm_file" -tag 0018 1310 2>/dev/null | cut -c 13- | head -n 1)
  137. local FovAP=$(dcminfo "$dcm_file" -tag 2005 1074 | cut -c 13-)
  138. local FovFH=$(dcminfo "$dcm_file" -tag 2005 1075 | cut -c 13-)
  139. local FovRL=$(dcminfo "$dcm_file" -tag 2005 1076 | cut -c 13-)
  140. # local echonumber=$(dcminfo "$dcm_file" -tag 0018 0086 2>/dev/null | cut -c 13- | head -n 1)
  141. # need to add local echonumber or something similar for mTE (0018,0086)
  142. # Now we need to determine what vendor it is.
  143. # Philips needs all the following calculations
  144. # Siemens works out of the box
  145. # GE most recent version also seem to work fine
  146. #echo $manufacturer
  147. if [ "$manufacturer" = "SIEMENS" ]; then
  148. slicetime_provided_by_vendor=1
  149. ees_trt_provided_by_vendor=1
  150. #elif [ "$manufacturer" = 'GE ]
  151. # need to be tested
  152. else
  153. slicetime_provided_by_vendor=0
  154. ees_trt_provided_by_vendor=0
  155. fi
  156. # 1/ Calculate ess/trt; needed are : FieldStrength, WaterFatShift, EPIFactor
  157. # Note: only calculate it when it is provided (sometimes this has been thrown away by anonymising the dicom-data)
  158. # We check whether needed tags exist
  159. tags_are_present=1
  160. test_waterfatshift=$(dcminfo "$dcm_file" -tag 2001 1022)
  161. if [ -z "$test_waterfatshift" ]; then
  162. tags_are_present=0
  163. local waterfatshift="empty"
  164. else
  165. local waterfatshift=$(dcminfo "$dcm_file" -tag 2001 1022 | awk '{print $(NF)}')
  166. fi
  167. test_fieldstrength=$(dcminfo "$dcm_file" -tag 0018 0087)
  168. if [ -z "$test_fieldstrength" ]; then
  169. tags_are_present=0
  170. local fieldstrength="empty"
  171. else
  172. local fieldstrength=$(dcminfo "$dcm_file" -tag 0018 0087 | awk '{print $(NF)}')
  173. fi
  174. test_epifactor=$(dcminfo "$dcm_file" -tag 2001 1013)
  175. if [ -z "$test_epifactor" ]; then
  176. tags_are_present=0
  177. local epifactor="empty"
  178. else
  179. local epifactor=$(dcminfo "$dcm_file" -tag 2001 1013 | awk '{print $(NF)}')
  180. fi
  181. if [ $tags_are_present -eq 0 ]; then
  182. ees_sec="empty"
  183. trt_sec="empty"
  184. else
  185. local water_fat_diff_ppm=3.3995
  186. local resonance_freq_mhz_tesla=42.576
  187. local water_fat_shift_hz=$(echo $fieldstrength $water_fat_diff_ppm $resonance_freq_mhz_tesla echo $fieldstrength $water_fat_diff_ppm | awk '{print $1 * $2 * $3}')
  188. #effective_echo_spacing_msec = 1000 * WFS_PIXEL/(water_fat_shift_hz * (EPI_FACTOR + 1))
  189. ees_sec=$(echo $waterfatshift $water_fat_shift_hz $epifactor | awk '{print $1 / ($2 * ($3 + 1))}')
  190. #total_readout_time_fsl_msec = EPI_FACTOR * effective_echo_spacing_msec;
  191. trt_sec=$(echo $epifactor $ees_sec | awk '{print $1 * $2 }')
  192. fi
  193. # 2/ Calculate slice SliceTiming
  194. #function SliceTime=KUL_slicetiming(MB, NS, TR)
  195. #%MB = 1; % Multiband factor
  196. #%NS = 30; % Number of Slices
  197. #%TR = 1.7; % TR in seconds
  198. #st = repmat(0:TR/(NS/MB):TR-.0000001,1,MB);
  199. #SliceTime = sprintf('%.8f,' , st);
  200. #SliceTime = ['"SliceTiming": [' SliceTime(1:end-1) ']'];
  201. #end
  202. multiband_factor=$mb
  203. # if mb is not found in the config file or is set to 0
  204. # simply set it to 1 and run
  205. if [[ -z ${multiband_factor} ]] || [[ ${multiband_factor} == 0 ]]; then
  206. multiband_factor=1
  207. fi
  208. tags_are_present=1
  209. test_number_of_slices=$(dcminfo "$dcm_file" -tag 2001 1018)
  210. if [ -z "$test_number_of_slices" ]; then
  211. tags_are_present=0
  212. local number_of_slices="empty"
  213. else
  214. local number_of_slices=$(dcminfo "$dcm_file" -tag 2001 1018 | awk '{print $(NF)}')
  215. fi
  216. test_repetion_time_msec=$(dcminfo "$dcm_file" -tag 0018 0080 )
  217. if [ -z "$test_repetion_time_msec" ]; then
  218. tags_are_present=0
  219. local repetion_time_msec="empty"
  220. else
  221. local repetion_time_msec=$(dcminfo "$dcm_file" -tag 0018 0080 | awk '{print $(NF)}')
  222. fi
  223. test_slice_scan_order=$(dcminfo "$dcm_file" -tag 2005 1081 )
  224. if [ -z "$test_slice_scan_order" ]; then
  225. #tags_are_present=0
  226. local slice_scan_order="empty"
  227. else
  228. local slice_scan_order=$(dcminfo "$dcm_file" -tag 2005 1081 | head -n 1 | awk '{print $(NF)}')
  229. fi
  230. if [ $tags_are_present -eq 1 ]; then
  231. if [ ! $multiband_factor = "" ];then
  232. #single_slice_time (in seconds)
  233. local single_slice_time=$(echo $repetion_time_msec $number_of_slices $multiband_factor | awk '{print $1 / ($2 / $3) / 1000}')
  234. # number of excitations given multiband
  235. # e = n. of excitations/slices per band
  236. local e=$(echo $number_of_slices $multiband_factor | awk '{print ($1 / $2) -1 }')
  237. local spb=$((${e}+1));
  238. echo $e
  239. echo $spb
  240. echo "${slice_scan_order}"
  241. # here we need to adapt to account for different slice orders
  242. # e.g.
  243. if [[ "${slice_scan_order}" == "rev. central" ]]; then
  244. # this is a bit different from interleaved... namely we split it into 2 gps
  245. # lower group is regular ascending and second group is regular descending
  246. half_e=$(echo "scale=2;(${spb}/2)" | bc | awk '{print int($1+0.5)}')
  247. for (( zc=0; zc<${half_e}; zc++ )); do
  248. sl1=$(echo $zc $single_slice_time | awk '{print $1 * $2}')
  249. if [[ -z ${slit1} ]]; then
  250. slit1="${sl1}"
  251. else
  252. slit1="${slit1}, ${sl1}"
  253. fi
  254. done
  255. for (( zx=${e}; zx>=${half_e}; zx-- )); do
  256. sl2=$(echo $zx $single_slice_time | awk '{print $1 * $2}')
  257. if [[ -z ${slit2} ]]; then
  258. slit2="${sl2}"
  259. else
  260. slit2="${slit2}, ${sl2}"
  261. fi
  262. done
  263. slit="${slit1}, ${slit2}"
  264. echo "${slit}"
  265. for iz in ${!tmp_order[@]}; do
  266. sl=$(echo $((${tmp_order[$iz]})) ${single_slice_time} | awk '{print $1 * $2}');
  267. echo ${sl}
  268. if [[ -z ${slit} ]]; then
  269. slit="${sl}"
  270. echo ${slit}
  271. else
  272. slit="${slit}, ${sl}"
  273. echo ${slit}
  274. fi
  275. done
  276. # interleaved is ready!
  277. elif [[ "${slice_scan_order}" == "interleaved" ]]; then
  278. declare -a tmp_order
  279. step=$(echo "sqrt(${spb})" | bc)
  280. unset tmp_order curr bh ik slgp;
  281. slgp=0;
  282. declare -a tmp_order;
  283. tmp_order[0]=0;
  284. for ik in $(seq 1 ${e}); do
  285. bh=$((${ik}-1));
  286. curr=$((${tmp_order[$bh]}+${step}));
  287. if [[ ${curr} -gt ${e} ]]; then
  288. ((slgp++));
  289. curr=${slgp};
  290. fi;
  291. tmp_order[$ik]=${curr};
  292. done;
  293. for iz in ${!tmp_order[@]}; do
  294. sl=$(echo $((${tmp_order[$iz]})) ${single_slice_time} | awk '{print $1 * $2}');
  295. echo ${sl}
  296. if [[ -z ${slit} ]]; then
  297. slit="${sl}"
  298. echo ${slit}
  299. else
  300. slit="${slit}, ${sl}"
  301. echo ${slit}
  302. fi
  303. done
  304. elif [[ "${slice_scan_order}" == "FH" ]]; then
  305. for (( c=0; c<=$e; c++ )); do
  306. sl=$(echo $c $single_slice_time | awk '{print $1 * $2}')
  307. if [[ -z ${slit} ]]; then
  308. slit="${sl}"
  309. else
  310. slit="${slit}, ${sl}"
  311. fi
  312. done
  313. elif [[ "${slice_scan_order}" == "HF" ]]; then
  314. for (( c=${e}; c>=0; c-- )); do
  315. sl=$(echo $c $single_slice_time | awk '{print $1 * $2}')
  316. if [[ -z ${slit} ]]; then
  317. slit="${sl}"
  318. else
  319. slit="${slit}, ${sl}"
  320. fi
  321. done
  322. elif [[ "${slice_scan_order}" == "default" ]]; then
  323. echo ${slice_scan_order}
  324. a=$((${spb} %2));
  325. echo ${a}
  326. # this is still untested
  327. if [[ "${spb}" -le 6 ]]; then
  328. step=2;
  329. hlpp=$(($((${e}+1))/${step}));
  330. lpsl=0;
  331. lpal=${lpsl};
  332. hpsl=${e};
  333. hpal=${hpsl};
  334. order=0;
  335. for ii in $(seq 0 1 ${spb}); do
  336. if [[ ${lpal} -lt $((${hlpp}-1)) ]]; then
  337. tmp_order[${order}]=${lpal};
  338. lpal=$((${lpal}+${step}));
  339. ((order++))
  340. elif [[ ${hpal} -ge $((${hlpp}-1)) ]]; then
  341. tmp_order[${order}]=${hpal};
  342. hpal=$((${hpal}-${step}));
  343. ((order++))
  344. else
  345. lpal=$((${lpsl}+1))
  346. hpal=$((${hpsl}-1))
  347. fi
  348. done
  349. # We will not add a 1 as done in the matlab version but we iterate over spb not e also
  350. for iz in ${!tmp_order[@]}; do
  351. sl=$(echo $((${tmp_order[$iz]})) ${single_slice_time} | awk '{print $1 * $2}');
  352. if [[ -z ${slit} ]]; then
  353. slit="${sl}"
  354. else
  355. slit="${slit}, ${sl}"
  356. fi
  357. done
  358. # this is still untested
  359. elif [[ "${spb}" == 8 ]]; then
  360. declare -a tmp_order
  361. step=$(echo "sqrt(${spb})" | bc)
  362. echo "step is ${step}"
  363. unset tmp_order curr bh ik slgp;
  364. slgp=0;
  365. declare -a tmp_order;
  366. tmp_order[0]=0;
  367. echo ${tmp_order[@]};
  368. for ik in $(seq 1 ${e}); do
  369. echo " ik is ${ik}";
  370. bh=$((${ik}-1));
  371. echo "bh is ${bh}";
  372. curr=$((${tmp_order[$bh]}+${step}));
  373. echo "tmp_order of bh is ${tmp_order[$bh]};
  374. echo "initially tmp_order of ik is ${tmp_order[$ik]};
  375. echo "step is ${step}";
  376. if [[ ${curr} -gt ${e} ]]; then
  377. echo "slgp is ${slgp}";
  378. ((slgp++));
  379. echo "inceremented slgp is ${slgp}";
  380. curr=${slgp};
  381. echo "now curr = ${slgp}";
  382. fi;
  383. tmp_order[$ik]=${curr};
  384. echo "tmp_order of ik is ${tmp_order[$ik]}";
  385. done;
  386. echo ${tmp_order[@]}
  387. for iz in ${!tmp_order[@]}; do
  388. sl=$(echo $((${tmp_order[$iz]})) ${single_slice_time} | awk '{print $1 * $2}');
  389. if [[ -z ${slit} ]]; then
  390. slit="${sl}"
  391. else
  392. slit="${slit}, ${sl}"
  393. fi
  394. done
  395. else
  396. # if we have an odd no. of slices per band
  397. if [[ ${a} == 1 ]]; then
  398. for zc in $(seq 0 2 ${e} ); do
  399. sl1=$(echo $zc $single_slice_time | awk '{print $1 * $2}')
  400. if [[ -z ${slit1} ]]; then
  401. slit1="${sl1}"
  402. else
  403. slit1="${slit1}, ${sl1}"
  404. fi
  405. done
  406. for zx in $(seq 1 2 ${e} ); do
  407. sl2=$(echo $zx $single_slice_time | awk '{print $1 * $2}')
  408. if [[ -z ${slit2} ]]; then
  409. slit2="${sl2}"
  410. else
  411. slit2="${slit2}, ${sl2}"
  412. fi
  413. done
  414. slit="${slit1}, ${slit2}"
  415. echo "${slit}"
  416. # if we have an odd no. of slices per band
  417. elif [[ ${a} == 0 ]]; then
  418. step=2;
  419. hlpp=$(($((${e}+1))/${step}));
  420. lpsl=0;
  421. lpal=${lpsl};
  422. hpsl=${e};
  423. hpal=${hpsl};
  424. order=0;
  425. for ii in $(seq 0 1 ${spb}); do
  426. if [[ ${lpal} -lt $((${hlpp}-1)) ]]; then
  427. tmp_order[${order}]=${lpal};
  428. lpal=$((${lpal}+${step}));
  429. ((order++))
  430. elif [[ ${hpal} -ge $((${hlpp}-1)) ]]; then
  431. tmp_order[${order}]=${hpal};
  432. hpal=$((${hpal}-${step}));
  433. ((order++))
  434. else
  435. lpal=$((${lpsl}+1))
  436. hpal=$((${hpsl}-1))
  437. fi
  438. done
  439. # We will not add a 1 as done in the matlab version but we iterate over spb not e also
  440. for iz in ${!tmp_order[@]}; do
  441. sl=$(echo $((${tmp_order[$iz]})) ${single_slice_time} | awk '{print $1 * $2}');
  442. if [[ -z ${slit} ]]; then
  443. slit="${sl}"
  444. else
  445. slit="${slit}, ${sl}"
  446. fi
  447. done
  448. fi
  449. fi
  450. fi
  451. slit2=$slit
  452. rep=$(echo $multiband_factor | awk '{print $1 - 1}')
  453. for (( c=1; c<=$rep; c++ )); do
  454. slit2="$slit2, $slit"
  455. done
  456. slice_time=[$slit2]
  457. fi
  458. else
  459. slice_time="empty"
  460. fi
  461. if [ $silent -eq 0 ]; then
  462. echo " patid = $patid"
  463. echo " the dicom file we are reading = $dcm_file"
  464. echo " manufacturer = $manufacturer"
  465. echo " software version = $software"
  466. echo " imagetype = $imagetype"
  467. echo " acquisitionMatrix = $acquisitionMatrix"
  468. echo " FovAP = $FovAP"
  469. echo " FovFH = $FovFH"
  470. echo " FovRL = $FovRL"
  471. echo " pixelspacing = $pixelspacing"
  472. echo " slicethickness = $slicethickness"
  473. echo " series = $seriesdescr"
  474. echo " fieldstrength = $fieldstrength"
  475. echo " waterfatshift = $waterfatshift"
  476. echo " epifactor = $epifactor"
  477. echo " calulated ees = $ees_sec"
  478. echo " calculated trt = $trt_sec"
  479. echo " number of slices = $number_of_slices"
  480. echo " repetion_time_msec = $repetion_time_msec"
  481. echo " slice_scan_order = $slice_scan_order"
  482. if [ ! $multiband_factor = "" ];then
  483. echo " multiband_factor = $mb"
  484. echo " calculated single_slice_time = $single_slice_time"
  485. echo " number of excitations - 1 = $e"
  486. echo " for 1 multiband = $slit"
  487. echo " complete slice_time = $slice_time"
  488. fi
  489. fi
  490. if [ ! -f $out ]; then
  491. echo -e "participant,session,dcm_file,manufacturer,software_version,series_descr,imagetype,fieldstrength,acquisitionMatrix,FovAP,FovFH,FovRL,pixelspacing,slicethickness,epifactor,wfs,ees_sec,trt_sec,nr_slices,slice_scan_order,repetion_time_msec,multiband_factor" > $out
  492. fi
  493. echo -e "$subj,${sess},$dcm_file,$manufacturer,$software,$seriesdescr,$imagetype,$fieldstrength,$acquisitionMatrix,$FovAP,$FovFH,$FovRL,$pixelspacing,$slicethickness,$epifactor,$waterfatshift,$ees_sec,$trt_sec,$number_of_slices,$slice_scan_order,$repetion_time_msec,$multiband_factor" >> $out
  494. }
  495. function kul_find_relevant_dicom_file {
  496. kul_e2cl " Searching for ${identifier} using search_string $search_string" $log
  497. # find the search_string in the dicom dump_file
  498. # search for search_string in dump_file, find ORIGINAL, remove dicom tags, sort, take first line, remove trailing space
  499. seq_file=$(grep $search_string $dump_file | grep ORIGINAL - | cut -f1 -d"[" | sort | head -n 1 | sed -e 's/[[:space:]]*$//')
  500. if [ "$seq_file" = "" ]; then
  501. kul_e2cl " ${identifier} dicoms are NOT FOUND" $log
  502. seq_found=0
  503. else
  504. seq_found=1
  505. kul_e2cl " a relevant ${identifier} dicom is $(basename "${seq_file}") " $log
  506. fi
  507. }
  508. # END LOCAL FUNCTIONS --------------
  509. # CHECK COMMAND LINE OPTIONS -------------
  510. #
  511. # Set defaults
  512. sess=""
  513. bids_output=BIDS
  514. # Set flags
  515. subj_flag=0
  516. sess_flag=0
  517. dcm_flag=0
  518. conf_flag=0
  519. bids_flag=0
  520. tmp_flag=0
  521. events_flag=0
  522. n_sbref_tasks=0
  523. silent=1
  524. anon=0
  525. if [ "$#" -lt 4 ]; then
  526. Usage >&2
  527. exit 1
  528. else
  529. while getopts "c:d:p:o:s:t:aveh" OPT; do
  530. case $OPT in
  531. d) #dicom_zip_file
  532. dcm_flag=1
  533. dcm=$OPTARG
  534. ;;
  535. p) #participant
  536. subj_flag=1
  537. subj=$OPTARG
  538. ;;
  539. c) #config_file
  540. conf_flag=1
  541. conf=$OPTARG
  542. ;;
  543. o) #bids output directory
  544. bids_flag=1
  545. bids_output=$OPTARG
  546. ;;
  547. s) #session
  548. sess_flag=1
  549. sess=$OPTARG
  550. ;;
  551. a) #pydeface
  552. anon=1
  553. ;;
  554. v) #verbose
  555. silent=0
  556. ;;
  557. e) #events for task based fmri
  558. events_flag=1
  559. ;;
  560. t) #temporary directory
  561. tmp_flag=1
  562. tempo=$OPTARG
  563. ;;
  564. h) #help
  565. Usage >&2
  566. exit 0
  567. ;;
  568. \?)
  569. echo "Invalid option: -$OPTARG" >&2
  570. echo
  571. Usage >&2
  572. exit 1
  573. ;;
  574. :)
  575. echo "Option -$OPTARG requires an argument." >&2
  576. echo
  577. Usage >&2
  578. exit 1
  579. ;;
  580. esac
  581. done
  582. fi
  583. # check for required options
  584. if [ $dcm_flag -eq 0 ] ; then
  585. echo
  586. echo "Option -d is required: give the file with your raw dicoms either .zip of tar.gz" >&2
  587. echo
  588. exit 2
  589. fi
  590. if [ $subj_flag -eq 0 ] ; then
  591. echo
  592. echo "Option -p is required: give the anonymised name of a subject this will create a directory subject_preproc with results." >&2
  593. echo
  594. exit 2
  595. fi
  596. if [ $conf_flag -eq 0 ] ; then
  597. echo
  598. echo "Option -c is required: give the path to the file that describes the sequences" >&2
  599. echo
  600. exit 2
  601. elif [ ! -f $conf ] ; then
  602. echo
  603. echo "The config file $conf does not exist"
  604. echo
  605. exit 2
  606. fi
  607. # INITIATE ---
  608. # main log file naming
  609. d=$(date "+%Y-%m-%d_%H-%M-%S")
  610. log=$log_dir/${subj}_main_log_${d}.txt
  611. # file for initial dicom tags
  612. dump_file=${log_dir}/${subj}_${sess}_initial_dicom_info.txt
  613. # file with final dicom tags
  614. final_dcm_tags_file=${log_dir}/${subj}_${sess}_final_dicom_info.csv
  615. # location of bids_config_json_file
  616. bids_config_json_file=${log_dir}/${subj}_${sess}_bids_config.json
  617. # location of dcm2niix_log_file
  618. dcm2niix_log_file=$log_dir/${subj}_${sess}_dcm2niix_log_file.txt
  619. if [[ $tmp_flag -eq 1 ]] ; then
  620. tmp=${cwd}/${tempo}
  621. else
  622. tmp="/tmp/${subj}"
  623. fi
  624. rm -fr ${tmp}
  625. # exit
  626. # remove previous existances to start fresh
  627. rm -f $dump_file
  628. rm -f $final_dcm_tags_file
  629. rm -f $bids_config_json_file
  630. # rm -fr ${cwd}/${tmp}/$subj
  631. # ----------- SAY HELLO ----------------------------------------------------------------------------------
  632. if [[ $silent -eq 0 ]]; then
  633. echo " The script you are running has basename `basename "$0"`, located in dirname $kul_main_dir"
  634. echo " The present working directory is `pwd`"
  635. fi
  636. # uncompress the zip file with dicoms or link the directory to tmp
  637. # clear the /tmp directory
  638. mkdir -p ${tmp}
  639. if [[ -d "$dcm" ]]; then
  640. # it is a directory
  641. kul_e2cl " you gave the directory $dcm as input; linking to to $tmp/$subj" $log
  642. ln -s "${cwd}/${dcm}" $tmp/$subj
  643. else
  644. kul_e2cl " uncompressing the zip file $dcm to $tmp/$subj" $log
  645. # Check the extention of the archive
  646. arch_ext="${dcm##*.}"
  647. #echo $arch_ext
  648. if [[ $arch_ext = "zip" ]]; then
  649. unzip -q -o ${dcm} -d ${tmp}
  650. elif [[] $arch_ext = "tar" ]]; then
  651. tar --strip-components=5 -C ${tmp} -xzf ${dcm}
  652. fi
  653. fi
  654. # dump the dicom tags of all dicoms in a file
  655. kul_e2cl " brute force extraction of some relevant dicom tags of all dicom files of subject $subj into file $dump_file" $log
  656. # check if a DICOMIR file exists and archive it (we do a brute force extraction and don't need the DICOMDIR file)
  657. test_DICOMDIR=($(find -L ${tmp} -type f -name "DICOMDIR" ))
  658. #echo "DICOMDIR = $test_DICOMDIR"
  659. if [[ $test_DICOMDIR == "" ]]; then
  660. echo " OK there is no DICOMDIR"
  661. else
  662. gzip $test_DICOMDIR
  663. fi
  664. # Do the bruce force extract
  665. echo hello > $dump_file
  666. task(){
  667. dcm1=$(dcminfo "$dcm_file" -tag 0008 103E -tag 0008 0008 -tag 0008 0070 -nthreads 4 2>/dev/null | tr -s '\n' ' ')
  668. echo "$dcm_file" $dcm1 >> $dump_file
  669. }
  670. (
  671. find -L ${tmp} -type f |
  672. while IFS= read -r dcm_file; do
  673. task
  674. done
  675. )
  676. kul_e2cl " done reading dicom tags of $dcm" $log
  677. # create empty bids description
  678. # bids=""
  679. # we read the config file
  680. declare -a sub_bids
  681. while IFS=, read identifier search_string task mb pe_dir acq_label; do
  682. bs=$(( $bs + 1))
  683. if [[ ! ${identifier} == \#* ]]; then
  684. if [[ ${identifier} == "T1w" ]]; then
  685. kul_find_relevant_dicom_file
  686. if [ $seq_found -eq 1 ]; then
  687. # read the relevant dicom tags
  688. kul_dcmtags "${seq_file}"
  689. sub_bids_T1='{"dataType": "anat", "modalityLabel": "T1w", "criteria": {
  690. "SeriesDescription": "*'${search_string}'*"}}'
  691. sub_bids_[$bs]=$(echo ${sub_bids_T1} | python -m json.tool )
  692. fi
  693. fi
  694. if [[ ${identifier} == "cT1w" ]]; then
  695. kul_find_relevant_dicom_file
  696. if [ $seq_found -eq 1 ]; then
  697. # read the relevant dicom tags
  698. kul_dcmtags "${seq_file}"
  699. sub_bids_T1='{"dataType": "anat", "modalityLabel": "T1w", "criteria": {
  700. "SeriesDescription": "*'${search_string}'*"},
  701. "customLabels": "ce-gadolinium",
  702. "sidecarChanges": {"KUL_dcm2bids": "yes","ContrastBolusIngredient": "gadolinium"}
  703. }'
  704. sub_bids_[$bs]=$(echo ${sub_bids_T1} | python -m json.tool )
  705. fi
  706. fi
  707. if [[ ${identifier} == "T2w" ]]; then
  708. kul_find_relevant_dicom_file
  709. if [ $seq_found -eq 1 ]; then
  710. # read the relevant dicom tags
  711. kul_dcmtags "${seq_file}"
  712. sub_bids_T2='{"dataType": "anat", "modalityLabel": "T2w", "criteria": {
  713. "SeriesDescription": "*'${search_string}'*"}}'
  714. sub_bids_[$bs]=$(echo ${sub_bids_T2} | python -m json.tool)
  715. fi
  716. fi
  717. if [[ ${identifier} == "PD" ]]; then
  718. kul_find_relevant_dicom_file
  719. The complete ASL
  720. if [ $seq_found -eq 1 ]; then
  721. # read the relevant dicom tags
  722. kul_dcmtags "${seq_file}"
  723. sub_bids_PD='{"dataType": "anat", "modalityLabel": "PD", "criteria": {
  724. "SeriesDescription": "*'${search_string}'*"}}'
  725. sub_bids_[$bs]=$(echo ${sub_bids_PD} | python -m json.tool)
  726. fi
  727. fi
  728. if [[ ${identifier} == "SWI" ]]; then
  729. kul_find_relevant_dicom_file
  730. if [ $seq_found -eq 1 ]; then
  731. # read the relevant dicom tags
  732. kul_dcmtags "${seq_file}"
  733. sub_bids_SWI='{"dataType": "anat", "modalityLabel": "SWI", "criteria": {
  734. "SeriesDescription": "*'${search_string}'*"}}'
  735. sub_bids_[$bs]=$(echo ${sub_bids_SWI} | python -m json.tool)
  736. fi
  737. fi
  738. if [[ ${identifier} == "MTI" ]]; then
  739. kul_find_relevant_dicom_file
  740. if [ $seq_found -eq 1 ]; then
  741. # read the relevant dicom tags
  742. kul_dcmtags "${seq_file}"
  743. sub_bids_SWI='{"dataType": "anat", "modalityLabel": "MTI", "criteria": {
  744. "SeriesDescription": "*'${search_string}'*","ImageType": [
  745. "ORIGINAL","PRIMARY","M","FFE","M","FFE"]}}'
  746. sub_bids_[$bs]=$(echo ${sub_bids_SWI} | python -m json.tool)
  747. fi
  748. fi
  749. if [[ ${identifier} == "ASL" ]]; then
  750. kul_find_relevant_dicom_file
  751. if [ $seq_found -eq 1 ]; then
  752. # read the relevant dicom tags
  753. kul_dcmtags "${seq_file}"
  754. sub_bids_SWI='{"dataType": "perf", "modalityLabel": "asl",
  755. "criteria": {
  756. "SeriesDescription": "*'${search_string}'*",
  757. "ImageType": [
  758. "ORIGINAL","PRIMARY","PERFUSION","NONE","REAL"
  759. ]}}'
  760. sub_bids_[$bs]=$(echo ${sub_bids_SWI} | python -m json.tool)
  761. fi
  762. fi
  763. if [[ ${identifier} == "FLAIR" ]]; then
  764. kul_find_relevant_dicom_file
  765. if [ $seq_found -eq 1 ]; then
  766. # read the relevant dicom tags
  767. kul_dcmtags "${seq_file}"
  768. sub_bids_FL='{"dataType": "anat", "modalityLabel": "FLAIR", "criteria": {
  769. "SeriesDescription": "*'${search_string}'*"}}'
  770. sub_bids_[$bs]=$(echo ${sub_bids_FL} | python -m json.tool)
  771. fi
  772. fi
  773. if [[ ${identifier} == "fmap" ]]; then
  774. kul_find_relevant_dicom_file
  775. if [ $seq_found -eq 1 ]; then
  776. # read the relevant dicom tags
  777. kul_dcmtags "${seq_file}"
  778. sub_bids_fm1='
  779. {
  780. "dataType": "fmap",
  781. "modalityLabel": "magnitude",
  782. "criteria":
  783. {
  784. "SeriesDescription": "*'${search_string}'*",
  785. "EchoNumber": 1
  786. }
  787. }'
  788. #echo $sub_bids_fm1
  789. sub_bids_fm1a=$(echo ${sub_bids_fm1} | python -m json.tool)
  790. sub_bids_fm2='
  791. {
  792. "dataType": "fmap",
  793. "modalityLabel": "fieldmap",
  794. "criteria":
  795. {
  796. "SeriesDescription": "'*${search_string}'*",
  797. "EchoNumber": 2
  798. },
  799. "sidecarChanges":
  800. {"Units": "Hz","IntendedFor": "##REPLACE_ME_INTENDED_FOR##"}
  801. }'
  802. #echo $sub_bids_fm2
  803. sub_bids_fm2a=$(echo ${sub_bids_fm2} | python -m json.tool)
  804. sub_bids_[$bs]="${sub_bids_fm1a},${sub_bids_fm2a}"
  805. echo $sub_bids_[$bs]
  806. fmap_task=$task
  807. intended_tasks_array=($task)
  808. fi
  809. fi
  810. if [[ ${identifier} == "func" ]]; then
  811. kul_find_relevant_dicom_file
  812. if [ $seq_found -eq 1 ]; then
  813. # read the relevant dicom tags
  814. kul_dcmtags "${seq_file}"
  815. # remove any whitespaces
  816. task_nospace="$(echo -e "${task}" | tr -d '[:space:]')"
  817. sub_bids_fu1='{"dataType": "func","modalityLabel":
  818. "bold","criteria": {"SeriesDescription": "*'${search_string}'*"},
  819. "customLabels": "task-'${task_nospace}''
  820. # add an acq_label if any
  821. if [ "$acq_label" = "" ];then
  822. sub_bids_fu1b='"',
  823. else
  824. sub_bids_fu1b='_acq-'${acq_label}'",'
  825. fi
  826. sub_bids_fu1c='"sidecarChanges": {"KUL_dcm2bids": "yes","TaskName": "'${task}'"'
  827. # for siemens (& ge?) ess/trt is not necessary as sidecarChanges
  828. # also not for philips, when it cannot be calculated
  829. #echo "ess_trt_provided_by_vendor: $ees_trt_provided_by_vendor"
  830. #echo "ees_sec: $ees_sec"
  831. if [ $ees_trt_provided_by_vendor -eq 1 ] || [ $ees_sec = "empty" ]; then
  832. if [ $ees_trt_provided_by_vendor -eq 1 ]; then
  833. kul_e2cl " It's a SIEMENS, ees/trt are in the dicom-header" $log
  834. else
  835. kul_e2cl " It's NOT original dicom data (anonymised?): ees/trt could not be calculated" $log
  836. fi
  837. sub_bids_fu2=""
  838. else
  839. kul_e2cl " It's a PHILIPS, ees/trt are were calculated" $log
  840. sub_bids_fu2=',"EffectiveEchoSpacing": '${ees_sec}',"TotalReadoutTime":
  841. '${trt_sec}',"MultibandAccelerationFactor": '${mb}',"PhaseEncodingDirection": "'${pe_dir}'"'
  842. fi
  843. # for siemens (& ge?) slicetiming is not necessary as sidecarChanges
  844. # also not for philips, when it cannot be calculated
  845. #echo "slicetime_provided_by_vendor: $slicetime_provided_by_vendor"
  846. #echo "slice_time: $slice_time"
  847. if [ $slicetime_provided_by_vendor -eq 1 ] || [ "$slice_time" = "empty" ]; then
  848. if [ $slicetime_provided_by_vendor -eq 1 ]; then
  849. kul_e2cl " It's a SIEMENS, slicetiming is in the dicom-header" $log
  850. else
  851. kul_e2cl " It's NOT original dicom data (anonymised?): slicetiming could not be calculated" $log
  852. fi
  853. sub_bids_fu3='}}'
  854. else
  855. kul_e2cl " It's a PHILIPS, slicetiming was calculated" $log
  856. sub_bids_fu3=',"SliceTiming": '${slice_time}'}}'
  857. fi
  858. sub_bids_[$bs]=$(echo ${sub_bids_fu1}${sub_bids_fu1b}${sub_bids_fu1c}${sub_bids_fu2}${sub_bids_fu3} | python -m json.tool)
  859. fi
  860. fi
  861. if [[ ${identifier} == "sbref" ]]; then
  862. kul_find_relevant_dicom_file
  863. if [ $seq_found -eq 1 ]; then
  864. # read the relevant dicom tags
  865. kul_dcmtags "${seq_file}"
  866. # Take care of fact that 1 sbref could be used for multiple funcs
  867. # we convert the first here, and copy later (see below)
  868. sbref_tasks=($task)
  869. n_sbref_tasks=${#sbref_tasks[@]}
  870. sbref_task1=${sbref_tasks[0]}
  871. #echo $sbref_task1
  872. sub_bids_sb1='{"dataType": "func","modalityLabel":
  873. "sbref","criteria": {"SeriesDescription": "*'${search_string}'*"},
  874. "customLabels": "task-'${sbref_task1}''
  875. if [ "$acq_label" = "" ];then
  876. sub_bids_sb1b='"',
  877. else
  878. sub_bids_sb1b='_acq-'${acq_label}'",'
  879. fi
  880. sub_bids_sb1c='"sidecarChanges": {"KUL_dcm2bids": "yes","TaskName": "'${sbref_task1}'"'
  881. # for siemens (& ge?) ess/trt is not necessary as sidecarChanges
  882. # also not for philips, when it cannot be calculated
  883. #echo "ess_trt_provided_by_vendor: $ees_trt_provided_by_vendor"
  884. #echo "ees_sec: $ees_sec"
  885. if [ $ees_trt_provided_by_vendor -eq 1 ] || [ $ees_sec = "empty" ]; then
  886. if [ $ees_trt_provided_by_vendor -eq 1 ]; then
  887. kul_e2cl " It's a SIEMENS, ees/trt are in the dicom-header" $log
  888. else
  889. kul_e2cl " It's NOT original dicom data (anonymised?): ees/trt could not be calculated" $log
  890. fi
  891. sub_bids_sb2=""
  892. else
  893. kul_e2cl " It's a PHILIPS, ees/trt are were calculated" $log
  894. sub_bids_sb2=',"EffectiveEchoSpacing": '${ees_sec}',"TotalReadoutTime":
  895. '${trt_sec}',"MultibandAccelerationFactor": '${mb}',"PhaseEncodingDirection": "'${pe_dir}'"'
  896. fi
  897. # for siemens (& ge?) slicetiming is not necessary as sidecarChanges
  898. # also not for philips, when it cannot be calculated
  899. #echo "slicetime_provided_by_vendor: $slicetime_provided_by_vendor"
  900. #echo "slice_time: $slice_time"
  901. if [ $slicetime_provided_by_vendor -eq 1 ] || [ "$slice_time" = "empty" ]; then
  902. if [ $slicetime_provided_by_vendor -eq 1 ]; then
  903. kul_e2cl " It's a SIEMENS, slicetiming is in the dicom-header" $log
  904. else
  905. kul_e2cl " It's NOT original dicom data (anonymised?): slicetiming could not be calculated" $log
  906. fi
  907. sub_bids_sb3='}}'
  908. else
  909. kul_e2cl " It's a PHILIPS, slicetiming was calculated" $log
  910. sub_bids_sb3=',"SliceTiming": '${slice_time}'}}'
  911. fi
  912. sub_bids_[$bs]=$(echo ${sub_bids_sb1}${sub_bids_sb1b}${sub_bids_sb1c}${sub_bids_sb2}${sub_bids_sb3} | python -m json.tool)
  913. fi
  914. fi
  915. if [[ ${identifier} == "dwi" ]]; then
  916. kul_find_relevant_dicom_file
  917. if [ $seq_found -eq 1 ]; then
  918. # read the relevant dicom tags
  919. kul_dcmtags "${seq_file}"
  920. sub_bids_dw1='{"dataType": "dwi","modalityLabel": "dwi",
  921. "criteria": {"SeriesDescription": "*'${search_string}'*"},'
  922. if [ "$acq_label" = "" ];then
  923. sub_bids_dw1b=""
  924. else
  925. sub_bids_dw1b='"customLabels": "acq-'${acq_label}'",'
  926. fi
  927. sub_bids_dw1c='"sidecarChanges": {"KUL_dcm2bids": "yes"'
  928. # for siemens (& ge?) ess/trt is not necessary as sidecarChanges
  929. # also not for philips, when it cannot be calculated
  930. #echo "ees_trt_provided_by_vendor: $ees_trt_provided_by_vendor"
  931. #echo "ees_sec: $ees_sec"
  932. if [ $ees_trt_provided_by_vendor -eq 1 ] || [ $ees_sec = "empty" ]; then
  933. if [ $ees_trt_provided_by_vendor -eq 1 ]; then
  934. kul_e2cl " It's a SIEMENS, ees/trt are in the dicom-header" $log
  935. else
  936. kul_e2cl " It's NOT original dicom data (anonymised?): ees/trt could not be calculated" $log
  937. fi
  938. sub_bids_dw2=""
  939. else
  940. kul_e2cl " It's a PHILIPS, ees/trt were calculated" $log
  941. sub_bids_dw2=',"EffectiveEchoSpacing": '${ees_sec}', "TotalReadoutTime": '${trt_sec}',
  942. "MultibandAccelerationFactor": '${mb}', "PhaseEncodingDirection": "'${pe_dir}'"'
  943. fi
  944. # for siemens (& ge?) slicetiming is not necessary as sidecarChanges
  945. # also not for philips, when it cannot be calculated
  946. #echo "slicetime_provided_by_vendor: $slicetime_provided_by_vendor"
  947. #echo "slice_time: $slice_time"
  948. if [ $slicetime_provided_by_vendor -eq 1 ] || [ "$slice_time" = "empty" ]; then
  949. if [ $slicetime_provided_by_vendor -eq 1 ]; then
  950. kul_e2cl " It's a SIEMENS, slicetiming is in the dicom-header" $log
  951. else
  952. kul_e2cl " It's NOT original dicom data (anonymised?): slicetiming could not be calculated" $log
  953. fi
  954. sub_bids_dw3='}}'
  955. else
  956. kul_e2cl " It's a PHILIPS, slicetiming was calculated" $log
  957. sub_bids_dw3=', "SliceTiming": '${slice_time}'}}'
  958. fi
  959. sub_bids_[$bs]=$(echo ${sub_bids_dw1}${sub_bids_dw1b}${sub_bids_dw1c}${sub_bids_dw2}${sub_bids_dw3} | python -m json.tool)
  960. fi
  961. fi
  962. fi
  963. done < $conf
  964. # make the full bids_conf string and write it to file
  965. bids_conf=""
  966. for bf in ${!sub_bids_[@]}; do
  967. #echo "now generating dcm2bids config files"
  968. if [[ ! ${bids_conf} ]]; then
  969. #echo "first field"
  970. bids_conf="${sub_bids_[$bf]}"
  971. else
  972. #echo "Series ${bf}"
  973. bids_conf="${bids_conf}, ${sub_bids_[$bf]}"
  974. fi
  975. done
  976. # WRITE THE .JSON FILE USED BY dcm2bids
  977. # Note: we also set dcm2niix options here and set pydeface
  978. json_anon=""
  979. if [ $anon -eq 1 ]; then
  980. $json_anon='"defaceTpl": "pydeface --outfile {dstFile} {srcFile}",'
  981. fi
  982. bids_conf_str="{${json_anon}
  983. \"dcm2niixOptions\": \"-b y -ba y -z y -i y -f '%3s_%f_%p_%t'\",
  984. \"descriptions\":[ ${bids_conf} ] }"
  985. #echo ${bids_conf_str}
  986. echo ${bids_conf_str} | python -m json.tool > ${bids_config_json_file}
  987. # MAIN HERE - WE RUN dcm2bids - HERE
  988. # invoke dcm2bids
  989. kul_e2cl " Calling dcm2bids... (for the output see $dcm2niix_log_file)" $log
  990. if [ ! -d BIDS ];then
  991. mkdir -p BIDS
  992. cd BIDS
  993. dcm2bids_scaffold
  994. echo "tmp_dcm2bids/*" > .bidsignore
  995. echo "*/anat/*SWI*" >> .bidsignore
  996. echo "*/anat/*MTI*" >> .bidsignore
  997. echo "*/perf/*asl*" >> .bidsignore
  998. cd ..
  999. fi
  1000. if [ $sess_flag -eq 1 ]; then
  1001. dcm2bids_session=" -s ${sess} "
  1002. else
  1003. dcm2bids_session=""
  1004. fi
  1005. dcm2bids -d "${tmp}" -p $subj $dcm2bids_session -c $bids_config_json_file \
  1006. -o $bids_output -l DEBUG --clobber > $dcm2niix_log_file
  1007. # Multi Echo func needs extra work. dcm2bids does not convert these correctly. "run" needs to be "echo"
  1008. if [[ ${sess} = "" ]] ; then
  1009. ses_long=""
  1010. else
  1011. ses_long="/ses-${sess}"
  1012. fi
  1013. me_file=($(grep EchoNumber ${bids_output}/sub-${subj}${ses_long}/func/*.json 2> /dev/null | awk -F ':' '{print $1}'))
  1014. me_echo=($(grep EchoNumber ${bids_output}/sub-${subj}${ses_long}/func/*.json 2> /dev/null | awk -F ':' '{print $3}' | cut -c 2 ))
  1015. if [[ ${me_file} = "" ]] ; then
  1016. echo " No Multiecho fMRI data found "
  1017. else
  1018. echo " Multiecho fMRI data found "
  1019. n_multi_echo=${#me_echo[@]}
  1020. for echo_number in $(seq 0 $(($n_multi_echo-1)) ) ; do
  1021. me_file_before_run=$(echo ${me_file[$echo_number]} | awk -F '_run-' '{print $1}')
  1022. me_file_after_run=$(echo ${me_file[$echo_number]} | awk -F '_run-' '{print $2}')
  1023. me_file_after_run=${me_file_after_run:2}
  1024. cmd_json="mv ${me_file[$echo_number]} ${me_file_before_run}_echo-${me_echo[$echo_number]}${me_file_after_run}"
  1025. cmd_nii=$(echo $cmd_json | perl -p -e 's/json/nii.gz/g')
  1026. eval $cmd_json
  1027. eval $cmd_nii
  1028. done
  1029. fi
  1030. # Update the Intended For of the fmaps
  1031. # Here we define the Intended For according to the BIDS specs
  1032. # It tells fmriprep how to use the fmap, notably for which func(s)
  1033. # fmap_task variable (1 strings or space-separated string) is used
  1034. if [[ ${fmap_task} ]] ; then
  1035. for intended_task in "${intended_tasks_array[@]}"; do
  1036. #echo $intended_task
  1037. #echo $cwd
  1038. search_runs_of_task=($(find ${cwd}/${bids_output}/sub-${subj}/func -type f | grep task-${intended_task} | grep nii.gz))
  1039. #echo ${search_runs_of_task[@]}
  1040. n_runs=${#search_runs_of_task[@]}
  1041. echo " we found $n_runs of task $intended_task"
  1042. for run_func in ${search_runs_of_task[@]}; do
  1043. if [[ $sess = "" ]]; then
  1044. intended_for_string="func\/$(basename $run_func)"
  1045. else
  1046. intended_for_string="ses-${sess}\/func\/$(basename $run_func)"
  1047. fi
  1048. full_intended_for_string="$full_intended_for_string, \"${intended_for_string}\""
  1049. done
  1050. done
  1051. full_intended_for_string="[ ${full_intended_for_string:1} ]"
  1052. # Now we replace it in the json file
  1053. echo " NOTE: we set the following string as Intended_For in the fieldmap: ${full_intended_for_string}"
  1054. if [[ $sess = "" ]]; then
  1055. perl -pi -e "s/\"##REPLACE_ME_INTENDED_FOR##\"/${full_intended_for_string}/g" ${bids_output}/sub-${subj}/fmap/sub-${subj}_fieldmap.json
  1056. else
  1057. perl -pi -e "s/\"##REPLACE_ME_INTENDED_FOR##\"/${full_intended_for_string}/g" ${bids_output}/sub-${subj}/ses-${sess}/fmap/sub-${subj}_fieldmap.json
  1058. fi
  1059. else
  1060. echo " No fmap tasks given "
  1061. fi
  1062. # Take care of fact that 1 sbref could be used for multiple funcs
  1063. # we convert the first above, now copy for each task
  1064. i=$((n_bref_tasks-1))
  1065. if [ $n_sbref_tasks -gt 1 ];then
  1066. for t in ${sbref_tasks[@]:0-$i}; do
  1067. cp ${bids_output}/sub-${subj}/func/sub-${subj}_task-${sbref_tasks[0]}_sbref.json \
  1068. ${bids_output}/sub-${subj}/func/sub-${subj}_task-${t}_sbref.json
  1069. cp ${bids_output}/sub-${subj}/func/sub-${subj}_task-${sbref_tasks[0]}_sbref.nii.gz \
  1070. ${bids_output}/sub-${subj}/func/sub-${subj}_task-${t}_sbref.nii.gz
  1071. done
  1072. fi
  1073. # copying task based events.tsv to BIDS directory
  1074. if [ $events_flag -eq 1 ]; then
  1075. test_events_exist=$(ls -l *conf*/task-*_events.tsv | grep "No such file")
  1076. #echo $test_events_exist
  1077. if [ "$test_events_exist" = "" ]; then
  1078. kul_e2cl "Copying task based events.tsv to BIDS directory" $log
  1079. cp *conf*/task-*_events.tsv $bids_output
  1080. fi
  1081. fi
  1082. # clean up
  1083. cleanup="rm -fr ${tmp}"
  1084. echo ${cleanup}
  1085. eval ${cleanup}
  1086. # Fix BIDS validation
  1087. echo "This BIDS was made using KUL_NeuroImagingTools" >> ${bids_output}/README
  1088. sed -i.bck 's/"Funding": ""/"Funding": [""]/' ${bids_output}/dataset_description.json
  1089. rm ${bids_output}/dataset_description.json.bck
  1090. # Run BIDS validation
  1091. docker run -ti --rm -v ${cwd}/${bids_output}:/data:ro bids/validator /data
  1092. kul_e2cl "Finished $script" $log