search

brainlife / app-mrtrix3-preproc

Run the recommended preprocessing procedure provided by mrtrix3. The options available mostly reflect the optimal DESIGNER pipeline that was recently proposed. This App runs for >15 on topup if both PA and AP dwi files are provided. It detects bvecs flipping (dwigradcheck) and update the gradient table accordingly.
3 stars 6 forks source link

a possible bug with converting from odd sliced b0 to even sliced b0 #14

Open soichih opened 3 years ago

soichih commented 3 years ago

I started seeing this error message on the latest version of app-mrtrix3-preproc (mrtrix 3.0.2)

+ mrconvert raw2.mif -grad cor2.b rpe_dwi.mif -nthreads 8 -quiet -force
mrconvert: [ERROR] number of studies in base image (6) does not match number of rows in diffusion gradient table (7)
mrconvert: [ERROR] error importing diffusion gradient table for image "raw2.mif"

I see 7 entries in cor2.b

# command_history: mrconvert -fslgrad ../6021a7e98b397412ae22be8e/5e849e7f952fef5e197a1771/dwi.bvecs ../6021a7e98b397412ae22be8e/5e849e7f952fef5e197a1771/dwi.bvals ../6021a7e98b397412ae22be8e/5e849e7f952fef5e197a1771/dwi.nii.gz raw2.mif --export_grad_mrtrix raw2.b -nthreads 8 -quiet -force  (version=3.0.2)
-0.5773502692 0.5773502692 0.5773502692 5
-0.5773502692 0.5773502692 0.5773502692 5
-0.5773502692 0.5773502692 0.5773502692 5
-0.5773502692 0.5773502692 0.5773502692 5
-0.5773502692 0.5773502692 0.5773502692 5
-0.5773502692 0.5773502692 0.5773502692 5
-4.697983962e-20 -3.666057883e-19 -1 30

However, the raw2.mif contains only 6 directions.

╰─$ mrinfo raw2.mif 
************************************************
Image:               "raw2.mif"
************************************************
  Dimensions:        140 x 128 x 92 x 6
  Voxel size:        1.5 x 1.5 x 1.5 x 3.495

Although the raw2.mif had initially had 7 volumes ...

++ mrinfo -size raw2.mif
++ grep -oE '[^[:space:]]+$'
+ nb0R=7
+ echo 'Reverse phase encoded dwi volume has 7 volumes.'

The script then run dwiextract -bzero to extract only the b0 volumes

+ echo 'The RPE file has an odd number of volumes. Only the b0 volumes were extracted.'
+ dwiextract -bzero raw2.mif raw2.mif -nthreads 8 -quiet -force
++ mrinfo -size raw2.mif
++ grep -oE '[^[:space:]]+$'
+ ob0=6
+ echo 'This should be an even number: 6'
+ echo 'RPE assigned as: pairs'

I believe we also need to update raw1.b (just re-export it from updated raw2.mif?)

soichih commented 3 years ago

Actually.. I am not sure why we are doing this

mrconvert raw2.mif -grad cor2.b rpe_dwi.mif -nthreads 8 -quiet -force

Can't we just remove -grad core2.b and let mrtrix use the correct bvector stored inside the raw2.mif? And if we do, what is the point of converting to rpe_dwi.mif?

bcmcpher commented 3 years ago

The convert appending cor2.b is replacing the gradient in raw2.mif with the flip checked gradient table.

When is this error occurring? All this text output indicates to me is you're passing a rpe input with the wrong bval/bvec.

soichih commented 3 years ago

Here is the stdout 

OMP_NUM_THREADS=
Converting input files to mrtrix format...
Forward phase encoded dwi volume has 80 volumes.
Reverse phase encoded dwi volume has 7 volumes.
The RPE file has an odd number of volumes. Only the b0 volumes were extracted.
This should be an even number: 6
RPE assigned as: pairs
creating dwimask (dwi2mask) from raw1.mif ...
Identifying correct gradient orientation...

and stderr

+ set -e
+ export SINGULARITYENV_OMP_NUM_THREADS=
+ SINGULARITYENV_OMP_NUM_THREADS=
+ export SINGULARITYENV_CUDA91PATH=/ocean/projects/dbs190008p/shayashi/cuda-9.1
+ SINGULARITYENV_CUDA91PATH=/ocean/projects/dbs190008p/shayashi/cuda-9.1
+ singularity exec -e --nv docker://brainlife/mrtrix3:3.0.2 ./mrtrix3_preproc.sh
INFO:    Using cached SIF image
+ set -e
+ echo OMP_NUM_THREADS=
+ '[' -z '' ']'
+ export OMP_NUM_THREADS=8
+ OMP_NUM_THREADS=8
+ export LD_LIBRARY_PATH=/usr/local/fsl/lib:/.singularity.d/libs:/ocean/projects/dbs190008p/shayashi/cuda-9.1/lib64
+ LD_LIBRARY_PATH=/usr/local/fsl/lib:/.singularity.d/libs:/ocean/projects/dbs190008p/shayashi/cuda-9.1/lib64
+ out=proc
+ difm=dwi
+ mask=b0_dwi_brain_mask
+ rm -rf ./tmp ./eddyqc cor1.b cor2.b corr.b
+ mkdir -p ./tmp
+ common='-nthreads 8 -quiet -force'
++ jq -r .diff config.json
+ DIFF=../6021a7e98b397412ae22be8e/5e849e7e952fef74e37a1721/dwi.nii.gz
++ jq -r .bval config.json
+ BVAL=../6021a7e98b397412ae22be8e/5e849e7e952fef74e37a1721/dwi.bvals
++ jq -r .bvec config.json
+ BVEC=../6021a7e98b397412ae22be8e/5e849e7e952fef74e37a1721/dwi.bvecs
++ jq -r .anat config.json
+ ANAT=../6021a7e98b397412ae22be8e/5e849e7e952fefcef47a170f/t1.nii.gz
++ jq -r .rdif config.json
+ RDIF=../6021a7e98b397412ae22be8e/5e849e7f952fef5e197a1771/dwi.nii.gz
++ jq -r .rbvl config.json
+ RBVL=../6021a7e98b397412ae22be8e/5e849e7f952fef5e197a1771/dwi.bvals
++ jq -r .rbvc config.json
+ RBVC=../6021a7e98b397412ae22be8e/5e849e7f952fef5e197a1771/dwi.bvecs
++ jq -r .round_bvals config.json
+ ROUND_BVALS=false
++ jq -r .acqd config.json
+ ACQD=AP
++ jq -r .denoise config.json
+ DO_DENOISE=true
++ jq -r .degibbs config.json
+ DO_DEGIBBS=true
++ jq -r .eddy config.json
+ DO_EDDY=true
++ jq -r .bias config.json
+ DO_BIAS=true
++ jq -r .ricn config.json
+ DO_RICN=true
++ jq -r .norm config.json
+ DO_NORM=false
++ jq -r .acpc config.json
+ DO_ACPC=true
++ jq -r .reslice config.json
+ NEW_RES=
++ jq -r .nval config.json
+ NORM=1000
++ jq -r .prct config.json
+ PRCT=null
++ jq -r .bias_method config.json
+ BIAS_METHOD=ants
++ jq -r .antsb config.json
+ ANTSB='[150,3]'
++ jq -r .antsc config.json
+ ANTSC='[200x200,1e-6]'
++ jq -r .antss config.json
+ ANTSS=2
+ common_fslpreproc='-eddy_mask b0_dwi_brain_mask.mif -eddyqc_all ./eddyqc -scratch ./tmp'
++ jq -r .eddy_data_is_shelled config.json
+ eddy_data_is_shelled=true
++ jq -r .eddy_slm config.json
+ eddy_slm=linear
++ jq -r .eddy_niter config.json
+ eddy_niter=5
++ jq -r .eddy_repol config.json
+ eddy_repol=true
++ jq -r .eddy_mporder config.json
+ eddy_mporder=0
+ eddy_options=' '
+ '[' true == true ']'
+ eddy_options='  --repol'
+ '[' true == true ']'
+ eddy_options='  --repol --data_is_shelled'
+ eddy_options='  --repol --data_is_shelled --slm=linear'
+ eddy_options='  --repol --data_is_shelled --slm=linear --niter=5'
+ '[' 0 '!=' 0 ']'
+ jq -rj .eddy_slspec config.json
+ '[' -s slspec.txt ']'
++ jq -r .topup_lambda config.json
+ topup_lambda=0.005,0.001,0.0001,0.000015,0.000005,0.0000005,0.00000005,0.0000000005,0.00000000001
+ topup_options=' '
+ '[' 0.005,0.001,0.0001,0.000015,0.000005,0.0000005,0.00000005,0.0000000005,0.00000000001 '!=' 0.005,0.001,0.0001,0.000015,0.000005,0.0000005,0.00000005,0.0000000005,0.00000000001 ']'
+ '[' -z ']'
+ DO_RESLICE=false
+ cp ../6021a7e98b397412ae22be8e/5e849e7e952fefcef47a170f/t1.nii.gz t1_acpc.nii.gz
+ chmod +w t1_acpc.nii.gz
+ ANAT=t1_acpc
+ echo 'Converting input files to mrtrix format...'
+ '[' false == true ']'
+ mrconvert -fslgrad ../6021a7e98b397412ae22be8e/5e849e7e952fef74e37a1721/dwi.bvecs ../6021a7e98b397412ae22be8e/5e849e7e952fef74e37a1721/dwi.bvals ../6021a7e98b397412ae22be8e/5e849e7e952fef74e37a1721/dwi.nii.gz raw1.mif --export_grad_mrtrix raw1.b -nthreads 8 -quiet -force
+ '[' -e ../6021a7e98b397412ae22be8e/5e849e7f952fef5e197a1771/dwi.nii.gz ']'
+ '[' false == true ']'
+ mrconvert -fslgrad ../6021a7e98b397412ae22be8e/5e849e7f952fef5e197a1771/dwi.bvecs ../6021a7e98b397412ae22be8e/5e849e7f952fef5e197a1771/dwi.bvals ../6021a7e98b397412ae22be8e/5e849e7f952fef5e197a1771/dwi.nii.gz raw2.mif --export_grad_mrtrix raw2.b -nthreads 8 -quiet -force
+ '[' ../6021a7e98b397412ae22be8e/5e849e7f952fef5e197a1771/dwi.nii.gz == null ']'
++ mrinfo -size raw1.mif
++ grep -oE '[^[:space:]]+$'
+ nb0F=80
++ mrinfo -size raw2.mif
++ grep -oE '[^[:space:]]+$'
+ nb0R=7
+ echo 'Forward phase encoded dwi volume has 80 volumes.'
+ echo 'Reverse phase encoded dwi volume has 7 volumes.'
+ '[' 80 -eq 7 ']'
+ RPE=pairs
+ '[' 1 == 0 ']'
+ echo 'The RPE file has an odd number of volumes. Only the b0 volumes were extracted.'
+ dwiextract -bzero raw2.mif raw2.mif -nthreads 8 -quiet -force
++ mrinfo -size raw2.mif
++ grep -oE '[^[:space:]]+$'
+ ob0=6
+ echo 'This should be an even number: 6'
+ echo 'RPE assigned as: pairs'
+ '[' pairs == all ']'
+ echo 'creating dwimask (dwi2mask) from raw1.mif ...'
+ dwi2mask raw1.mif b0_dwi_brain_mask.mif -nthreads 8 -quiet -force
+ echo 'Identifying correct gradient orientation...'
+ dwigradcheck raw1.mif -grad raw1.b -mask b0_dwi_brain_mask.mif -export_grad_mrtrix cor1.b -scratch ./tmp -nthreads 8 -quiet -force
dwigradcheck: Testing gradient table alterations (0 of 48)... [===================================================]
Mean length     Axis flipped    Axis permutations    Axis basis
52.95         none                (0, 1, 2)           image
52.18         none                (0, 1, 2)           scanner
31.79            0                (0, 2, 1)           scanner
31.73         none                (0, 2, 1)           image
31.21         none                (0, 2, 1)           scanner
30.73            0                (0, 2, 1)           image
28.89            1                (0, 1, 2)           scanner
28.79            1                (0, 1, 2)           image
28.77         none                (1, 0, 2)           image
28.68            2                (1, 0, 2)           scanner
28.51            2                (1, 0, 2)           image
28.43         none                (1, 0, 2)           scanner
28.13            0                (0, 1, 2)           scanner
27.90            0                (0, 1, 2)           image
27.70            2                (0, 1, 2)           scanner
27.69            2                (0, 1, 2)           image
27.47         none                (2, 1, 0)           scanner
27.04            1                (2, 1, 0)           image
26.98            1                (2, 1, 0)           scanner
26.68         none                (2, 1, 0)           image
24.32            1                (0, 2, 1)           image
24.26            1                (0, 2, 1)           scanner
24.11            2                (1, 2, 0)           image
23.99            2                (1, 2, 0)           scanner
23.99         none                (1, 2, 0)           scanner
23.70            2                (0, 2, 1)           image
23.69            2                (0, 2, 1)           scanner
23.51         none                (1, 2, 0)           image
23.07            1                (1, 0, 2)           scanner
23.02            1                (1, 2, 0)           scanner
22.90            0                (1, 2, 0)           image
22.84            0                (1, 0, 2)           image
22.80         none                (2, 0, 1)           image
22.69            1                (1, 0, 2)           image
22.61            0                (2, 0, 1)           scanner
22.58            0                (1, 0, 2)           scanner
22.48            1                (2, 0, 1)           scanner
22.20            2                (2, 0, 1)           image
22.19            0                (1, 2, 0)           scanner
22.09            1                (2, 0, 1)           image
22.05            2                (2, 0, 1)           scanner
22.04         none                (2, 0, 1)           scanner
21.91            0                (2, 1, 0)           scanner
21.85            0                (2, 0, 1)           image
21.78            1                (1, 2, 0)           image
21.43            2                (2, 1, 0)           scanner
21.41            2                (2, 1, 0)           image
21.39            0                (2, 1, 0)           image
+ mrconvert raw1.mif -grad cor1.b dwi.mif -nthreads 8 -quiet -force
+ '[' -e raw2.mif ']'
+ dwi2mask raw2.mif rpe_b0_dwi_brain_mask.mif -nthreads 8 -quiet -force
+ cp raw2.b cor2.b
+ mrconvert raw2.mif -grad cor2.b rpe_dwi.mif -nthreads 8 -quiet -force
mrconvert: [ERROR] number of studies in base image (6) does not match number of rows in diffusion gradient table (7)
mrconvert: [ERROR] error importing diffusion gradient table for image "raw2.mif"

real    1m54.012s
user    7m40.547s
sys 0m6.043s
soichih commented 3 years ago

@bcmcpher Thanks for the patch (https://github.com/brainlife/app-mrtrix3-preproc/commit/11e79fb807af19accf0f407c4eb467e6b63d676f)

Unfortunately, I am now seeing the following error message. 

+ mrconvert raw1.mif -grad cor1.b dwi.mif -nthreads 8 -quiet -force
+ '[' -e raw2.mif ']'
+ dwi2mask raw2.mif rpe_b0_dwi_brain_mask.mif -nthreads 8 -quiet -force
+ mrconvert raw2.mif -clear dw_scheme rpe_dwi.mif -nthreads 8 -quiet -force
+ '[' true == true ']'
+ echo 'Performing PCA denoising (dwidenoise)...'
+ dwidenoise -extent 5,5,5 -noise fpe_noise.mif -estimator Exp2 dwi.mif dwi_denoise.mif -nthreads 8 -quiet -force
+ '[' -e rpe_dwi.mif ']'
+ dwidenoise -extent 5,5,5 -noise rpe_noise.mif -estimator Exp2 rpe_dwi.mif rpe_dwi_denoise.mif -nthreads 8 -quiet -force
+ mrcalc fpe_noise.mif rpe_noise.mif -add 2 -divide noise.mif -nthreads 8 -quiet -force
+ '[' true == true ']'
+ difm=dwi_denoise
+ '[' true == true ']'
+ echo 'Performing Gibbs ringing correction (mrdegibbs)...'
+ mrdegibbs -nshifts 20 -minW 1 -maxW 3 dwi_denoise.mif dwi_denoise_degibbs.mif -nthreads 8 -quiet -force
+ '[' -e rpe_dwi_denoise.mif ']'
+ mrdegibbs -nshifts 20 -minW 1 -maxW 3 rpe_dwi_denoise.mif rpe_dwi_denoise_degibbs.mif -nthreads 8 -quiet -force
+ difm=dwi_denoise_degibbs
+ '[' true == true ']'
+ echo 'Performing Eddy correction (dwifslpreproc)... rpe:pairs'
+ '[' pairs == none ']'
+ '[' pairs == pairs ']'
+ dwiextract -bzero dwi_denoise_degibbs.mif fpe_b0.mif -nthreads 8 -quiet -force
+ dwiextract -bzero rpe_dwi_denoise_degibbs.mif rpe_b0.mif -nthreads 8 -quiet -force
dwiextract: [ERROR] no valid diffusion gradient table found
dwiextract: [ERROR] error importing diffusion gradient table for image "rpe_dwi_denoise_degibbs.mif"

Do you know what this error message means?

soichih commented 3 years ago

FYI.. I've confirmed that mrconvert used to not check for volume count not matching with gradient vector count (on 3.0_RC2) but it does this check for 3.0.2 release.