[WIP] NORDIC denoising on dMRI images

Dockerfile

@@ -358,6 +358,10 @@ RUN set -uex; \
     rm itksnap.tar.gz
 ENV PATH="/opt/c3d-1.0.0-Linux-x86_64/bin:${PATH}"
 
+# Install NORDIC denoising
+RUN git clone https://github.com/SteenMoeller/NORDIC_Raw.git /opt/NORDIC_Raw && cd /opt/NORDIC_Raw && git checkout v1.1
+ENV NORDIC_Raw="/opt/NORDIC_Raw"
+
 COPY . /opt/micapipe/
 
 RUN bash -c 'cd /opt/micapipe && mv fix_settings.sh /opt/fix1.068/settings.sh && mv fsl_conf/* /opt/fsl-6.0.2/etc/flirtsch/'

functions/02_proc-dwi.sh

@@ -23,16 +23,17 @@ nocleanup=$5
 threads=$6
 tmpDir=$7
 dwi_main=$8
-dwi_rpe=$9
-dwi_processed=${10}
-rpe_all=${11}
-regAffine=${12}
-dwi_str=${13}
-b0thr=${14}
-bvalscale=${15}
-synth_reg=${16}
-dwi_upsample=${17}
-PROC=${18}
+dwi_phase=$9
+dwi_rpe=${10}
+dwi_processed=${11}
+rpe_all=${12}
+regAffine=${13}
+dwi_str=${14}
+b0thr=${15}
+bvalscale=${16}
+synth_reg=${17}
+dwi_upsample=${18}
+PROC=${19}
 here=$(pwd)
 
 #------------------------------------------------------------------------------#
@@ -51,6 +52,7 @@ bids_variables "$BIDS" "$id" "$out" "$SES"
 Info "Inputs of proc_dwi"
 Note "tmpDir        :" "$tmpDir"
 Note "dwi_main      :" "$dwi_main"
+Note "dwi_phase     :" "$dwi_phase"
 Note "dwi_rpe       :" "$dwi_rpe"
 Note "rpe_all       :" "$rpe_all"
 Note "dwi_acq       :" "$dwi_str"
@@ -80,6 +82,14 @@ if [[ "$dwi_main" != "DEFAULT" ]]; then
     IFS=',' read -ra bids_dwis <<< "$dwi_main"
     bids_dwis=("${bids_dwis[@]}")
 fi
+# Manage manual inputs: DWI phase image(s)
+if [[ "$dwi_phase" != "DEFAULT" ]]; then
+    IFS=',' read -ra bids_phase_dwis <<< "$dwi_phase"
+    bids_phase_dwis=("${bids_phase_dwis[@]}")
+    for i in "${bids_phase_dwis[@]}"; do
+      if [[ ! -f ${i} ]]; then Error "Provided dwi_phase image's path is wrong of file doesn't exist!!, check:\n\tls $dwi_phase"; exit; fi
+    done
+fi
 # Manage manual inputs: DWI reverse phase encoding
 if [[ "$dwi_rpe" != "DEFAULT" ]]; then
     IFS=',' read -ra dwi_reverse <<< "$dwi_rpe"
@@ -166,10 +176,31 @@ if [[ "$dwi_processed" == "FALSE" ]] && [[ ! -f "$dwi_corr" ]]; then
     if [ ! -f "$dwi_res" ] || [ ! -f "$dwi_dns" ]; then ((N++))
     Info "DWI denoise and concatenation"
     # Concatenate shells -if only one shell then just convert to mif and rename.
+    i=0
           for dwi in "${bids_dwis[@]}"; do
                 dwi_nom=$(echo "${dwi##*/}" | awk -F ".nii" '{print $1}')
                 bids_dwi_str=$(echo "$dwi" | awk -F . '{print $1}')
-                Do_cmd mrconvert "$dwi" -json_import "${bids_dwi_str}.json" -fslgrad "${bids_dwi_str}.bvec" "${bids_dwi_str}.bval" "${tmp}/${dwi_nom}.mif" "${bvalstr}"
+                # if phase images are present then use nordic denoising
+                if [[ -f "${bids_phase_dwis[i]}" && ("${#bids_phase_dwis[@]}" -eq "${#bids_dwis[@]}") ]]; then
+                  Info "Phase image found, running NORDIC denoising!"
+                  # Run MATLAB command with the specified arguments
+                  matlab -nodisplay -nojvm -nosplash -nodesktop -r " \
+                  try; \
+                  addpath('${NORDIC_Raw}/'); \
+                  ARG.temporal_phase=3; \
+                  ARG.phase_filter_width=3; \
+                  ARG.DIROUT = '${tmp}/'; \
+                  NIFTI_NORDIC('${dwi}', '${bids_phase_dwis[i]}', '${dwi_nom}_nordic', ARG); \
+                  end; \
+                  quit;"
+                  #>> ${ARG_DIROUT}/log_NORDIC_$(date '+%Y-%m-%d').txt
+                  Do_cmd mrconvert "${tmp}/${dwi_nom}_nordic.nii" -json_import "${bids_dwi_str}.json" -fslgrad "${bids_dwi_str}.bvec" "${bids_dwi_str}.bval" "${tmp}/${dwi_nom}.mif" "${bvalstr}"
+                  nordic_run=true
+                else
+                  Info "Phase image not found or incomplete, cannot run NORDIC denoising! Verify or add flag -dwi_phase to use phase"
+                  Do_cmd mrconvert "${dwi}" -json_import "${bids_dwi_str}.json" -fslgrad "${bids_dwi_str}.bvec" "${bids_dwi_str}.bval" "${tmp}/${dwi_nom}.mif" "${bvalstr}"
+                fi
+                i=$[$i +1]
                 Do_cmd dwiextract "${tmp}/${dwi_nom}.mif" "${tmp}/${dwi_nom}_b0.mif" -bzero
                 Do_cmd mrmath "${tmp}/${dwi_nom}_b0.mif" mean "${tmp}/${dwi_nom}_b0.nii.gz" -axis 3 -nthreads "$threads"
           done
@@ -203,12 +234,17 @@ if [[ "$dwi_processed" == "FALSE" ]] && [[ ! -f "$dwi_corr" ]]; then
           fi
 
           # Denoise DWI and calculate residuals
-          Info "DWI MP-PCA denoising and Gibbs ringing correction"
-          dwi_dns_tmp="${tmp}/MP-PCA_dwi.mif"
-          Do_cmd dwidenoise "$dwi_cat" "$dwi_dns_tmp" -nthreads "$threads"
-          mrcalc "$dwi_cat" "$dwi_dns_tmp" -subtract - -nthreads "$threads" | mrmath - mean "$dwi_resPCA" -axis 3
-          Do_cmd mrdegibbs "$dwi_dns_tmp" "$dwi_dns" -nthreads "$threads"
-          mrcalc "$dwi_dns_tmp" "$dwi_dns" -subtract - -nthreads "$threads" | mrmath - mean "$dwi_resGibss" -axis 3
+          Info "DWI MP-PCA denoising (if NORDIC not applied) and Gibbs ringing correction"
+          if [[ ${nordic_run} = true ]] ; then
+            Do_cmd mrdegibbs "$dwi_cat" "$dwi_dns" -nthreads "$threads"
+            mrcalc "$dwi_cat" "$dwi_dns" -subtract - -nthreads "$threads" | mrmath - mean "$dwi_resGibss" -axis
+          else
+            dwi_dns_tmp="${tmp}/MP-PCA_dwi.mif"
+            Do_cmd dwidenoise "$dwi_cat" "$dwi_dns_tmp" -nthreads "$threads"
+            mrcalc "$dwi_cat" "$dwi_dns_tmp" -subtract - -nthreads "$threads" | mrmath - mean "$dwi_resPCA" -axis 3
+            Do_cmd mrdegibbs "$dwi_dns_tmp" "$dwi_dns" -nthreads "$threads"
+            mrcalc "$dwi_dns_tmp" "$dwi_dns" -subtract - -nthreads "$threads" | mrmath - mean "$dwi_resGibss" -axis 3
+          fi
           ((Nsteps++))
     else
           Info "Subject ${id} has DWI in mif, denoised and concatenaded"; ((Nsteps++)); ((N++))

functions/utilities.sh

@@ -96,6 +96,7 @@ export idBIDS="${subject}${ses}"
   # BIDS Files
   bids_T1ws=($(ls "$subject_bids"/anat/*T1w.nii* 2>/dev/null))
   bids_dwis=($(ls "${subject_bids}/dwi/${subject}${ses}"*_dir-AP_*dwi.nii* 2>/dev/null))
+  bids_phase_dwis=($(ls "${subject_bids}/dwi/${subject}${ses}"*_dir-AP_*part-phase_dwi.nii* 2>/dev/null))
   bids_T1map=$(ls "$subject_bids"/anat/*mp2rage*T1map.nii* 2>/dev/null)
   bids_inv1=$(ls "$subject_bids"/anat/*inv1*T1map.nii* 2>/dev/null)
   bids_inv2=$(ls "$subject_bids"/anat/*inv2*T1map.nii* 2>/dev/null)
@@ -184,13 +185,17 @@ bids_print.variables-post() {
 bids_print.variables-dwi() {
   # This functions prints BIDS variables names and files if found
   Info "Variables for DWI processing"
-  Note "proc_dwi dir    :" "$proc_dwi"
-  Note "bids_dwis       :" "N-${#bids_dwis[@]}, $bids_dwis"
-  Note "dwi_reverse     :" "N-${#dwi_reverse[@]}, $dwi_reverse"
-
-  Note "T1 nativepro    :" "$(find "$T1nativepro" 2>/dev/null)"
-  Note "T1 5tt          :" "$(find "$T15ttgen" 2>/dev/null)"
-  Note "MNI152_mask     :" "$MNI152_mask"
+  Note "proc_dwi dir        :" "$proc_dwi"
+  for i in "${!bids_dwis[@]}"; do
+    file.exist "bids_dwis $[$i+1]/${#bids_dwis[@]}       :" ${bids_dwis[i]}
+  done
+  for i in "${!bids_dwis[@]}"; do
+    file.exist "bids_phase_dwis $[$i+1]/${#bids_dwis[@]} :" ${bids_phase_dwis[i]}
+  done
+  Note "dwi_reverse         :" "N-${#dwi_reverse[@]}, $dwi_reverse"
+  Note "T1 nativepro        :" "$(find "$T1nativepro" 2>/dev/null)"
+  Note "T1 5tt              :" "$(find "$T15ttgen" 2>/dev/null)"
+  Note "MNI152_mask         :" "$MNI152_mask"
 }
 
 bids_print.variables-func() {
@@ -261,6 +266,7 @@ bids_variables_unset() {
   unset icafixTraining
   unset bids_T1ws
   unset bids_dwis
+  unset bids_phase_dwis
   unset bids_T1map
   unset bids_inv1
   unset dwi_reverse

generate_micapipe_images.sh

@@ -68,6 +68,8 @@ generate() {
     --run-bash "wget https://sourceforge.net/projects/c3d/files/c3d/1.0.0/c3d-1.0.0-Linux-x86_64.tar.gz/download -O itksnap.tar.gz &&
                 tar -xfv itksnap.tar.gz -C /opt/" \
     --env PATH="/opt/itksnap/bin/:${PATH}" \
+    --run-bash "git clone https://github.com/SteenMoeller/NORDIC_Raw.git /opt/NORDIC_Raw && cd /opt/NORDIC_Raw && git checkout v1.1" \
+    --env NORDIC_Raw="/opt/NORDIC_Raw"
     --run-bash "wget https://www.dropbox.com/s/47lu1nojrderls1/install_R_env.sh?dl=0 -O /opt/install_R_env.sh &&
                 bash /opt/install_R_env.sh && cd /opt/afni-latest && rPkgsInstall -pkgs ALL" \
     --copy . /opt/micapipe \

micapipe

@@ -82,6 +82,8 @@ echo -e "
 \t\033[38;5;197m-proc_dwi\033[0m
 \t   \033[38;5;120m-dwi_main\033[0m            : Full path to DWI with N number of directions and b0. If used, it will overwrite the defaults.
 \t\t\t            Default = <bids>/<subject>/dwi/*_dir-AP_dwi.nii*
+\t   \033[38;5;120m-dwi_phase\033[0m           : Full path to DWI phase part with N number of directions. If used, it will overwrite the defaults.
+\t\t\t            Default = <bids>/<subject>/dwi/*_dir-AP_part-phase_dwi.nii*
 \t   \033[38;5;120m-dwi_rpe\033[0m             : DWI b0 image(s) with phase reversal encoding. If used, it will overwrite the defaults.
 \t\t\t            Default = <bids>/<subject>/dwi/<sub>_dir-PA_dwi.nii*
 \t\t\t            Set to 'FALSE' if no DWI with reverse fase encoding will be provided
@@ -455,6 +457,10 @@ do
     dwi_main="$2"
     shift;shift;
   ;;
+  -dwi_phase)
+    dwi_phase="$2"
+    shift;shift;
+  ;;
   -dwi_rpe)
     dwi_rpe="$2"
     shift;shift;
@@ -608,6 +614,14 @@ else
     for i in "${!dwis_main[@]}"; do dwis_main[i]=$(realpath ${dwis_main[$i]}); done # Full path
     dwi_main=$(IFS=','; echo "${dwis_main[*]}")
 fi
+# User provided DWI phase image(s)
+if [ -z ${dwi_phase} ]; then
+    dwi_phase=DEFAULT
+else
+    IFS=',' read -ra dwis_phase <<< $dwi_phase
+    for i in "${!dwis_phase[@]}"; do dwis_phase[i]=$(realpath ${dwis_phase[$i]}); done # Full path
+    dwi_phase=$(IFS=','; echo "${dwis_phase[*]}")
+fi
 # User provided DWI reverse phase encoded
 if [ -z ${dwi_rpe} ]; then
     dwi_rpe=DEFAULT
@@ -831,7 +845,7 @@ fi
 if [ "$procDWI" = "TRUE" ]; then
     if [[ ${dwi_str} != "DEFAULT" ]]; then log_str=_${dwi_str}; else log_str=""; fi
     log_file_str=$dir_logs/proc_dwi_$(date +'%Y-%m-%d_%H.%M.%S')${log_str}
-    COMMAND="${scriptDir}/02_proc-dwi.sh $BIDS $id $out $SES $nocleanup $threads $tmpDir $dwi_main $dwi_rpe $dwi_processed $rpe_all $regAffine $dwi_str $b0thr $bvalscale $synth_reg $dwi_upsample"
+    COMMAND="${scriptDir}/02_proc-dwi.sh $BIDS $id $out $SES $nocleanup $threads $tmpDir $dwi_main $dwi_phase $dwi_rpe $dwi_processed $rpe_all $regAffine $dwi_str $b0thr $bvalscale $synth_reg $dwi_upsample"
     # mica.q - Diffusion processing
     if [[ $micaq == "TRUE" ]]; then
         Info "MICA qsub - Diffusion processing"