User Tools

Site Tools


biac:analysis:topup_correction

Differences

This shows you the differences between two versions of the page.

Link to this comparison view

Both sides previous revision Previous revision
Next revision
Previous revision
Next revision Both sides next revision
biac:analysis:topup_correction [2019/02/21 18:41]
cmp12
biac:analysis:topup_correction [2019/09/04 17:45]
cmp12
Line 7: Line 7:
 There will be two series of single timepoint EPI images.  You can grab relevant info to create the acq_params.txt files from your XML header. There will be two series of single timepoint EPI images.  You can grab relevant info to create the acq_params.txt files from your XML header.
  
-You need to calculate readout time in seconds ( the physical time it takes to get the acquisition matrix of a single slice ) and get the polarity direction ( phase encode direction ).+You need to calculate readout time in seconds of the PEPOLAR images ( the physical time it takes to get the acquisition matrix of a single slice ) and get the polarity direction ( phase encode direction ).
  
 the readout time in seconds for the parameter file will be: the readout time in seconds for the parameter file will be:
Line 15: Line 15:
 </code> </code>
  
 +unfortunately at this point the polarity of the images will have to be determined from visual inspection. we aren't provided enough information in the metadata to give an entry into the BXH file ( yet ).
 +
 +Here is a rough guide to help with inspection.  Take note of the eyeballs being crushed in for AP and stretched out for PA.  Below shows the PEPOLAR images in red/blue on top of a mean functional.
 +{{:biac:analysis:connneuro_func.png?800|}}
  
-the polarity for the entry will have to be determined from the **seriesdescription**, which is typically "field map reverse polarity" or "field map regular" 
  
 reverse will be "-1" in the acq_params.txt and regular will be "1" reverse will be "-1" in the acq_params.txt and regular will be "1"
Line 62: Line 65:
 applytopup --imain=../bia6_00186_009_01.nii.gz --inindex=1 --method=jac --datain=acq_params.txt --topup=rs_topup --out=run009 --verbose applytopup --imain=../bia6_00186_009_01.nii.gz --inindex=1 --method=jac --datain=acq_params.txt --topup=rs_topup --out=run009 --verbose
 applytopup --imain=../bia6_00186_010_01.nii.gz --inindex=1 --method=jac --datain=acq_params.txt --topup=rs_topup --out=run010 --verbose applytopup --imain=../bia6_00186_010_01.nii.gz --inindex=1 --method=jac --datain=acq_params.txt --topup=rs_topup --out=run010 --verbose
 +</code>
 +
 +If you have the scenario where the functional data was acquired with acceleration, but the PEPOLAR images were not.  You can create an additional acq_params file with a modified readout time to help prevent over correction.  You will need to divide the echospacing by the sensefactor from your BXH header of the functional run.  If your functional data and PEPOLAR images have the same sensefactor, this compensation is not needed.
 +
 +<code>readout = ((echospacing/sensefactor) * (acquisitionmatrix[0] * (percentsampling/100))) / 1e6 
 +echospacing in the BXH header is in microseconds
 </code> </code>
  
Line 180: Line 189:
  
 [[https://mrtrix.readthedocs.io/en/latest/|mrtrix3]] has a nice wrapper for FSL's eddy/topup/apply_topup that can do the above corrections and along with eddy corrections [[https://mrtrix.readthedocs.io/en/latest/|mrtrix3]] has a nice wrapper for FSL's eddy/topup/apply_topup that can do the above corrections and along with eddy corrections
 +
 +If you are running dwidenoise, do it BEFORE dwipreproc. 
  
 ==== Calculate your readout time ==== ==== Calculate your readout time ====
Line 190: Line 201:
  
 there are 2 scenarios that typically apply here there are 2 scenarios that typically apply here
-1) a short acquisition with RPE B0s +1) DWI and a short acquisition with RPE B0s 
-2) an entire acquisition with the same gradient table and reverse phase encoding +2) entire acquisitions with the same gradient table and reverse phase encoding  
 + 
 +In both scenarios it is important to create your data with normal phase encoding direction first, followed by reversed.
  
 ===== Scenario 1 ===== ===== Scenario 1 =====
 +
 +==== Prepare your input datasets ====
 +
 +Specify that a set of images (typically b=0 volumes) will be provided for use in inhomogeneity field estimation only
  
 For scenario 1, create your blip up / blip down B0 data the same way as above. For scenario 1, create your blip up / blip down B0 data the same way as above.
 <code> <code>
-bxhselect --timeselect 0:1 bia6_00197_012.bxh bu  +bxhselect --timeselect 0 bi 
-bxhselect --timeselect 0:1 bia6_00197_013.bxh bd+a6_00197_012.bxh bu  
 +bxhselect --timeselect 0 bia6_00197_013.bxh bd
 fslmerge -t bud bu.nii.gz bd.nii.gz fslmerge -t bud bu.nii.gz bd.nii.gz
-<code>+rm bu.* bd.* 
 +</code>
  
 Prepare your dwi series for input to mrtrix3.  If you use their mif data format, then the gradient tables will be automatically passed throughout their software with additional flags. Prepare your dwi series for input to mrtrix3.  If you use their mif data format, then the gradient tables will be automatically passed throughout their software with additional flags.
Line 211: Line 230:
 mrconvert -fslgrad bvecs bvals bia6_00414_013.nii.gz dwi.mif  mrconvert -fslgrad bvecs bvals bia6_00414_013.nii.gz dwi.mif 
  
-#you can verify the table  +#export the table back out in their format  
-tesla:00414 cmp12$ mrinfo dwi.mif -dwgrad +mrinfo dwi.mif -export_grad_mrtrix grads.b 
-          0                               0 + 
-          0                               0 +#put THOSE back into the mif header 
-          1                            1500 +mrconvert -grad grads.b dwi.mif dwi.mif -force
-   0.456784    0.889578               1501.42 +
-  -0.264866    0.484757    0.833581     1501.51 +
-   0.937268   -0.284081    0.202058     1499.14 +
-  0.0630209   -0.806262    0.588192     1499.02 +
-  -0.679047  -0.0940064     0.72805     1499.79 +
-   0.686982    0.278993    0.670983     1500.08 +
-  -0.232946    0.949777    0.208951      1500.7 +
-   0.922569    0.178916    0.341841      1501.4 +
-  -0.111021   -0.980188   -0.164032     1499.43 +
-  -0.132016   -0.854107   -0.503064     1499.62+
  
 </code> </code>
  
 +The last step seems redundant, but there is an issue with how mrtrix3 tools are handling the multi-shell GE data produced on our scanners.  If you data is NOT multi-shell, then you don't need it.  If your data IS multi-shell, then putting it back in can workaround other issues until their RC4 code is released.  [[https://github.com/MRtrix3/mrtrix3/issues/1487|discussion here]]
  
- +==== Run dwipreproc -rpe_pair ====
  
 +Run [[https://mrtrix.readthedocs.io/en/latest/reference/scripts/dwipreproc.html|dwipreproc]] with the "-rpe_pair" option
 +
 +<code>
 +dwipre
 +proc dwi.mif dwi_corr.mif -rpe_pair -se_epi bud.nii.gz -pe_dir AP -readout_time 0.10656  -debug
 +
 +-rpe_pair specifies you're providing a pair of B0s ( regular, reversed )
 +-se_epi is you bud image
 +-pe_dir is your dwi's phase encode direction
 +-readout_time is from calculation above
 +</code>
    
 +dwipreproc will run topup, eddy and apply topup.  eddy will output update diffusion directions and insert them into your dwi_corr.mif output header.
 +
 +You can generate a BXH header for a mif file with:
 +<code>
 +mif2bxh dwi_corr.mif dwi_corr.bxh
 +</code>
 +
 +===== Scenario 2 =====
 +
 +You have an entire acquisition with the same gradient table as your DWI sequence, but with reverse phase encoding.  This information will be used to perform a recombination of image volumes (each pair of volumes with the same b-vector but different phase encoding directions will be combined together into a single volume).
 +
 +Calculate the readout time in the same way as above
 +
 +==== Prepare your input datasets ====
 +
 +This time you'll concatenate the 2 DWI acquisitions.  regular, then reversed 
 +
 +<code>
 +#concat the 2 series ( regular, reversed ) 
 +bxh_concat bia6_00260_00
 +7_LAS.bxh bia6_00260_008_LAS.bxh both
 +
 +#extract the gradients
 +extractdiffdirs --fsl both.bxh bvecs bvals
 +
 +#convert to mif
 +mrconvert -fslgrad bvecs bvals both.nii.gz dwi.mif
 +
 +#re-extract table
 +mrinfo dwi.mif -export_grad_mrtrix grads.b
 +
 +#put corrected table back in
 +mrconvert -grad grads.b dwi.mif dwi.mif -force
 +
 +#remove the temporary files ( bvecs, bvals, both.* )
 +</code>
 +
 +
 +==== Run dwipreproc -rpe_all ====
 +
 +Run dwipreproc with the "-rpe_all" option
 +
 +<code>
 +dwipreproc dwi.mif dwi_corr.mif -rpe_all -pe_dir AP -readout_time 0.10656 -debug
 +
 +-rpe_all signals that you've replicated 
 +A
 +LL the directions with a rpe acquisition
 +-pe_dir is the phase encode direction of your regular acquisition
 +-readout_time from above  
 +</code>
 +
 +
 +{{ :biac:analysis:trace_raw_corr.png?800 |}}
 +
biac/analysis/topup_correction.txt · Last modified: 2023/02/23 18:43 (external edit)