From 30a9514a18dd655f6c750c9cb11cdccbd04b7bb4 Mon Sep 17 00:00:00 2001 From: willtack Date: Thu, 20 May 2021 17:09:57 -0400 Subject: [PATCH] pw output --- app/for_testing/readme.txt | 102 ------------------ app/for_testing/run_full_analysis.sh | 33 ------ manifest.json | 4 +- run | 39 +++++-- src/full_analysis.m | 154 ++------------------------- 5 files changed, 40 insertions(+), 292 deletions(-) delete mode 100644 app/for_testing/readme.txt delete mode 100755 app/for_testing/run_full_analysis.sh diff --git a/app/for_testing/readme.txt b/app/for_testing/readme.txt deleted file mode 100644 index 698611f..0000000 --- a/app/for_testing/readme.txt +++ /dev/null @@ -1,102 +0,0 @@ -full_analysis Executable - -1. Prerequisites for Deployment - -Verify that version 9.6 (R2019a) of the MATLAB Runtime is installed. -If not, you can run the MATLAB Runtime installer. -To find its location, enter - - >>mcrinstaller - -at the MATLAB prompt. - -Alternatively, download and install the Linux version of the MATLAB Runtime for R2019a -from the following link on the MathWorks website: - - http://www.mathworks.com/products/compiler/mcr/index.html - -For more information about the MATLAB Runtime and the MATLAB Runtime installer, see -"Distribute Applications" in the MATLAB Compiler documentation -in the MathWorks Documentation Center. - -2. Files to Deploy and Package - -Files to Package for Standalone -================================ --full_analysis --run_full_analysis.sh (shell script for temporarily setting environment variables and - executing the application) - -to run the shell script, type - - ./run_full_analysis.sh - - at Linux or Mac command prompt. is the directory - where version 9.6 of the MATLAB Runtime is installed or the directory where - MATLAB is installed on the machine. is all the - arguments you want to pass to your application. For example, - - If you have version 9.6 of the MATLAB Runtime installed in - /mathworks/home/application/v96, run the shell script as: - - ./run_full_analysis.sh /mathworks/home/application/v96 - - If you have MATLAB installed in /mathworks/devel/application/matlab, - run the shell script as: - - ./run_full_analysis.sh /mathworks/devel/application/matlab --MCRInstaller.zip - Note: if end users are unable to download the MATLAB Runtime using the - instructions in the previous section, include it when building your - component by clicking the "Runtime included in package" link in the - Deployment Tool. --This readme file - - - -3. Definitions - -For information on deployment terminology, go to -http://www.mathworks.com/help and select MATLAB Compiler > -Getting Started > About Application Deployment > -Deployment Product Terms in the MathWorks Documentation -Center. - -4. Appendix - -A. Linux systems: -In the following directions, replace MR/v96 by the directory on the target machine where - MATLAB is installed, or MR by the directory where the MATLAB Runtime is installed. - -(1) Set the environment variable XAPPLRESDIR to this value: - -MR/v96/X11/app-defaults - - -(2) If the environment variable LD_LIBRARY_PATH is undefined, set it to the following: - -MR/v96/runtime/glnxa64:MR/v96/bin/glnxa64:MR/v96/sys/os/glnxa64:MR/v96/sys/opengl/lib/glnxa64 - -If it is defined, set it to the following: - -${LD_LIBRARY_PATH}:MR/v96/runtime/glnxa64:MR/v96/bin/glnxa64:MR/v96/sys/os/glnxa64:MR/v96/sys/opengl/lib/glnxa64 - - For more detailed information about setting the MATLAB Runtime paths, see Package and - Distribute in the MATLAB Compiler documentation in the MathWorks Documentation Center. - - - - NOTE: To make these changes persistent after logout on Linux - or Mac machines, modify the .cshrc file to include this - setenv command. - NOTE: The environment variable syntax utilizes forward - slashes (/), delimited by colons (:). - NOTE: When deploying standalone applications, you can - run the shell script file run_full_analysis.sh - instead of setting environment variables. See - section 2 "Files to Deploy and Package". - - - - - - diff --git a/app/for_testing/run_full_analysis.sh b/app/for_testing/run_full_analysis.sh deleted file mode 100755 index 0372147..0000000 --- a/app/for_testing/run_full_analysis.sh +++ /dev/null @@ -1,33 +0,0 @@ -#!/bin/sh -# script for execution of deployed applications -# -# Sets up the MATLAB Runtime environment for the current $ARCH and executes -# the specified command. -# -exe_name=$0 -exe_dir=`dirname "$0"` -echo "------------------------------------------" -if [ "x$1" = "x" ]; then - echo Usage: - echo $0 \ args -else - echo Setting up environment variables - MCRROOT="$1" - echo --- - LD_LIBRARY_PATH=.:${MCRROOT}/runtime/glnxa64 ; - LD_LIBRARY_PATH=${LD_LIBRARY_PATH}:${MCRROOT}/bin/glnxa64 ; - LD_LIBRARY_PATH=${LD_LIBRARY_PATH}:${MCRROOT}/sys/os/glnxa64; - LD_LIBRARY_PATH=${LD_LIBRARY_PATH}:${MCRROOT}/sys/opengl/lib/glnxa64; - export LD_LIBRARY_PATH; - echo LD_LIBRARY_PATH is ${LD_LIBRARY_PATH}; - shift 1 - args= - while [ $# -gt 0 ]; do - token=$1 - args="${args} \"${token}\"" - shift - done - eval "\"${exe_dir}/full_analysis\"" $args -fi -exit - diff --git a/manifest.json b/manifest.json index fff7446..f35a4d9 100755 --- a/manifest.json +++ b/manifest.json @@ -2,7 +2,7 @@ "name": "hasl-decode", "label": "HASL Decoder: Hadamard-encoded ASL (HASL) processing pipeline", "description": "Hadamard-encoded multi-delay ASL data(HASL) processing pipeline for measuring cerebrovascular reactivity. Decodes the Hadamard matrix into three post-labelling delays and calculates cerebral blood flow, arterial transit time and transit time-corrected cerebral blood flow.", - "version": "0.5.1_0.5.1", + "version": "0.6.2_0.6.2", "inputs": { "asl": { "base": "file", @@ -41,7 +41,7 @@ "custom": { "gear-builder": { "category": "analysis", - "image": "willtack/hasl:0.5.1" + "image": "willtack/hasl:0.6.2" }, "flywheel": { "suite": "BrainScienceCenter" diff --git a/run b/run index 354d0ec..8cd8930 100755 --- a/run +++ b/run @@ -68,20 +68,23 @@ if [[ $exit_status != 0 ]]; then fi # Rename output files -TT_MAP_ORIG=$(find $OUTPUT_DIR -type f | grep tt.nii) -TT_MAP="${OUTPUT_DIR}"/"${FILE_PREFIX}_tt.nii" +TT_FILE_ORIG=$(find $OUTPUT_DIR -type f | grep tt.nii) +TT_FILE="${OUTPUT_DIR}"/"${FILE_PREFIX}_tt.nii" TT_MNI="${OUTPUT_DIR}"/"${FILE_PREFIX}_MNI152_tt.nii" CBF_FILE="${OUTPUT_DIR}"/"${FILE_PREFIX}_cbf.nii" TTCBF_FILE="${OUTPUT_DIR}"/"${FILE_PREFIX}_ttcbf.nii" -mv $TT_MAP_ORIG $TT_MAP +PW_FILE="${OUTPUT_DIR}"/"${FILE_PREFIX}_pw.nii" +mv $TT_FILE_ORIG $TT_FILE mv $(find $OUTPUT_DIR -type f | grep -v tt | grep cbf.nii) $CBF_FILE mv $(find $OUTPUT_DIR -type f | grep ttcbf.nii) $TTCBF_FILE +mv $(find $OUTPUT_DIR -type f | grep pw.nii) $PW_FILE + # Register TT map to MNI152 via MPRAGE scan -flirt -ref "${mprageFileBet}" -in "${TT_MAP}" -omat $(pwd)/tt2struc.mat +flirt -ref "${mprageFileBet}" -in "${TT_FILE}" -omat $(pwd)/tt2struc.mat flirt -ref ${FSLDIR}/data/standard/MNI152_T1_2mm_brain -in "${mprageFileBet}" -omat $(pwd)/struc2mni_affine.mat fnirt --in="${mprageFile}" --aff=struc2mni_affine.mat --cout=struc2standard_warp --config=T1_2_MNI152_2mm -applywarp --ref=${FSLDIR}/data/standard/MNI152_T1_2mm.nii.gz --in=${TT_MAP} --warp=struc2standard_warp \ +applywarp --ref=${FSLDIR}/data/standard/MNI152_T1_2mm.nii.gz --in=${TT_FILE} --warp=struc2standard_warp \ --premat=tt2struc.mat --out=${TT_MNI} # Invert warps to bring vascular territory masks into same space as transit time map @@ -94,7 +97,7 @@ mkdir "${TERRITORIES_DIR}" # Apply warps to every mask for mask in ${FLYWHEEL_BASE}/territories/*.nii.gz; do new_mask_name="${FILE_PREFIX}_$(basename ${mask} .nii.gz)_native.nii.gz" - applywarp --ref="${TT_MAP}" --in="${mask}" --warp=struc2standard_warp_inv --postmat=struc2tt.mat --out=${TERRITORIES_DIR}/${new_mask_name} + applywarp --ref="${TT_FILE}" --in="${mask}" --warp=struc2standard_warp_inv --postmat=struc2tt.mat --out=${TERRITORIES_DIR}/${new_mask_name} fslmaths ${TERRITORIES_DIR}/${new_mask_name} -thrP 50 -bin ${TERRITORIES_DIR}/${new_mask_name} done @@ -105,11 +108,27 @@ wmSeg=$(find . -type f | grep _seg_2) for x in /flywheel/v0/*_seg*.nii.gz; do basename=$(basename $x .nii.gz) output="${TERRITORIES_DIR}"/${basename}_native.nii.gz - flirt -in $x -ref ${TT_MAP} -applyxfm -init struc2tt.mat -out $output + flirt -in $x -ref ${TT_FILE} -applyxfm -init struc2tt.mat -out $output fslmaths $output -thrP 50 -bin $output # cp $output ${TERRITORIES_DIR}/ done +# split perfusion weighted image and create deltaM/M0 images +deltaM_DIR="${OUTPUT_DIR}/deltaM" +mkdir "${deltaM_DIR}" +fslsplit $PW_FILE "${deltaM_DIR}/${FILE_PREFIX}_deltaM_" +for x in /flywheel/v0/output/deltaM/*; do + basename=$(basename $x .nii.gz) + fslmaths "${x}" -div "${m0File}" "${deltaM_DIR}/${basename}_dMM0" +done +rename 's/0007/control/' "${deltaM_DIR}/*.nii.gz" + +# warp file stuff +WARPS_DIR="${OUTPUT_DIR}/warp_files" +mkdir ${WARPS_DIR} +mv struc2standard*.nii.gz "${WARPS_DIR}/" +mv *.mat "${WARPS_DIR}/" + # Get a list of the files in the output directory outputs=$(find $OUTPUT_DIR/* -maxdepth 0 -type f) @@ -118,15 +137,13 @@ if [[ -z $outputs ]]; then echo "$CONTAINER No results found in output directory... Exiting" exit 1 else - WARPS_DIR="${OUTPUT_DIR}/warp_files" - mkdir ${WARPS_DIR} - mv struc2standard*.nii.gz "${WARPS_DIR}/" - mv *.mat "${WARPS_DIR}/" cd $OUTPUT_DIR + pushd deltaM; zip -r ../"${FILE_PREFIX}"_deltaM.zip *; popd pushd vascular_territories; zip -r ../"${FILE_PREFIX}"_vascular_territories.zip *; popd pushd warp_files; zip -r ../"${FILE_PREFIX}"_warp_files.zip *; popd rm -rf MPRAGE ASL M0 ${TERRITORIES_DIR} # remove the intermediate directories. just keep the derived maps rm -rf ${WARPS_DIR} + rm -rf ${deltaM_DIR} echo -e "$CONTAINER Wrote: `ls`" echo -e "$CONTAINER Done!" diff --git a/src/full_analysis.m b/src/full_analysis.m index 50a764d..665ecf3 100644 --- a/src/full_analysis.m +++ b/src/full_analysis.m @@ -1,8 +1,4 @@ -<<<<<<< HEAD function []=full_analysis(aslpath, m0path, mpragepath, outputdir) -======= -function []=full_analysis(aslpath, m0path, mpragepath, hyperC02, stimulus_start, outputdir) ->>>>>>> 4d0410692edfc486821256d898544b48d4c43d84 %-------------------------------------------------------------------------------------------------- % Step 1 @@ -12,14 +8,6 @@ %-------------------------------------------------------------------------------------------------- fprintf("Starting Step 1") subj_folder = outputdir; % data dir -<<<<<<< HEAD -% stimulus_start = str2double(stimulus_start); -% hyperC02 = logical(str2double(hyperC02)); -======= -stimulus_start = str2double(stimulus_start); -hyperC02 = logical(str2double(hyperC02)); ->>>>>>> 4d0410692edfc486821256d898544b48d4c43d84 -% hyperC02_bool = logical(hyperC02); [~, asl_name, asl_ext] = fileparts(aslpath); [~, m0_name, m0_ext] = fileparts(m0path); @@ -128,7 +116,6 @@ m0_path = nii_phase_extract(hasl_m0_path, 'M0', 1); % % set parameters -<<<<<<< HEAD asl_para = auxil_asl_para_init(); asl_para.LD = 3.5; @@ -142,21 +129,6 @@ % % smooth kernel pixel_size = nii_pixel_size(hasl_asl_path); flt_kernel = auxil_msk_gen_kernel_gaussian(pixel_size, 2.5); -======= -asl_para = hasl_para_init(); - -asl_para.LD = 3.5; -asl_para.PLD = 1; -asl_para.PLD_Num = 3; -asl_para.PLD_Lin = 1; -asl_para.T1b = 1.65; -asl_para = hasl_para_proc_state(asl_para); -asl_para = hasl_para_proc_ld_pld(asl_para); - -% % smooth kernel -pixel_size = nii_pixel_size(hasl_asl_path); -smooth_kernel = msk_gen_kernel_gaussian(pixel_size, 4.0); ->>>>>>> 4d0410692edfc486821256d898544b48d4c43d84 % % asl_img = nii_load_dimg(hasl_asl_path); @@ -166,13 +138,22 @@ state_num = asl_para.State_Num; loop_num = phase_num./state_num; -<<<<<<< HEAD % state = 1:phase_num; asl_img = auxil_asl_mean(asl_img, asl_para); pw_img = auxil_hasl_decode(asl_img, asl_para); pw_img = convn(pw_img, flt_kernel, 'same'); +% save perfusion-weighted image +nii = load_nii(hasl_asl_path); +nii.img = pw_img; +nii.hdr.dime.dim(5) = 8; % 7 PLDs + 1 control +save_nii(nii, fullfile(subj_folder, 'pw.nii')) + +% ninfo = niftiinfo(m0path); +% ninfo.ImageSize(4) = 7 % 7 PLDs +% niftiwrite(pw_img, fullfile(subj_folder, 'pw'), ninfo); + tt_img = auxil_asl_calc_tt(pw_img, brainmsk, asl_para); cbf_img = auxil_asl_cbf(pw_img, m0_img, brainmsk, asl_para); @@ -196,118 +177,3 @@ nii.hdr.dime.dim(5)=1; nii.img = tt_img; save_nii(nii, fullfile(subj_folder,'tt.nii')); -======= -nii = load_nii(hasl_m0_path); - - -% % % % % -if ~hyperC02 - state = 1:phase_num; - hasl_img = hasl_anymean(asl_img, asl_para, state, loop_num); - hasl_pw = hasl_decode(hasl_img, asl_para); - cbf_img = hasl_cbf(hasl_pw, m0_img, brainmsk, asl_para); - cbf_img = hasl_filter_apply(cbf_img, brainmsk, smooth_kernel); - tt_img = hasl_calc_tt(hasl_pw, brainmsk, asl_para); - tt_img = hasl_filter_apply(tt_img, brainmsk, smooth_kernel); - ttccbf_img = hasl_ttccbf(hasl_pw, m0_img, tt_img, brainmsk, asl_para); - - % % save nii for CBF, Transit Time and ttCBF - nii.hdr.dime.dim(5)=1; - nii.img = cbf_img(:,:,:,4); - save_nii(nii, fullfile(subj_folder, 'cbf.nii')); - - nii.hdr.dime.dim(5)=1; - nii.img = ttccbf_img(:,:,:,4); - save_nii(nii, fullfile(subj_folder, 'ttcbf.nii')); - - tt_img = tt_img.*1000; % *1000:change unit from 's' to 'ms' - nii.hdr.dime.dim(5)=1; - nii.img = tt_img; - save_nii(nii, fullfile(subj_folder,'tt.nii')); - -else - - normalCO2_state = 1:stimulus_start-3; - hyperCO2_state = stimulus_start:phase_num; % ASL phase numbers - - % % - normalCO2stat_img = asl_img(:,:,:,normalCO2_state); - hyperCO2stat_img = asl_img(:,:,:,hyperCO2_state); - - % % decode hasl - normalCO2_hasl_img = hasl_anymean(normalCO2stat_img, asl_para, normalCO2_state, loop_num); - normalCO2_hasl_pw = hasl_decode(normalCO2_hasl_img, asl_para); - - hyperCO2_hasl_img = hasl_anymean(hyperCO2stat_img, asl_para, hyperCO2_state, loop_num); - hyperCO2_hasl_pw = hasl_decode(hyperCO2_hasl_img, asl_para); - - % % CBF,Transit time,ttCBF calculation in different CO2 state - normalCO2cbf_img = hasl_cbf(normalCO2_hasl_pw, m0_img, brainmsk, asl_para); - normalCO2cbf_img = hasl_filter_apply(normalCO2cbf_img, brainmsk, smooth_kernel); - - hyperCO2cbf_img = hasl_cbf(hyperCO2_hasl_pw, m0_img, brainmsk, asl_para); - hyperCO2cbf_img = hasl_filter_apply(hyperCO2cbf_img, brainmsk, smooth_kernel); - - normalCO2tt_img = hasl_calc_tt(normalCO2_hasl_pw, brainmsk, asl_para); - normalCO2tt_img = hasl_filter_apply(normalCO2tt_img, brainmsk, smooth_kernel); - - hyperCO2tt_img = hasl_calc_tt(hyperCO2_hasl_pw, brainmsk, asl_para); - hyperCO2tt_img = hasl_filter_apply(hyperCO2tt_img, brainmsk, smooth_kernel); - - normalCO2ttccbf_img = hasl_ttccbf(normalCO2_hasl_pw, m0_img, normalCO2tt_img, brainmsk, asl_para); - normalCO2ttccbf_img = hasl_filter_apply(normalCO2ttccbf_img, brainmsk, smooth_kernel); - - hyperCO2ttccbf_img = hasl_ttccbf(hyperCO2_hasl_pw, m0_img, hyperCO2tt_img, brainmsk, asl_para); - hyperCO2ttccbf_img = hasl_filter_apply(hyperCO2ttccbf_img, brainmsk, smooth_kernel); - - % % save nii for normalCO2(BL) and hyperCO2(HC) state CBF, Transit Time and ttCBF - nii.hdr.dime.dim(5)=1; - nii.img = normalCO2cbf_img(:,:,:,4); - save_nii(nii, fullfile(subj_folder, 'normalCO2_cbf.nii')); - nii.hdr.dime.dim(5)=1; - nii.img = hyperCO2cbf_img(:,:,:,4); - save_nii(nii, fullfile(subj_folder, 'hyperCO2_cbf.nii')); - - nii.hdr.dime.dim(5)=1; - nii.img = normalCO2ttccbf_img(:,:,:,4); - save_nii(nii, fullfile(subj_folder, 'normalCO2_ttcbf.nii')); - nii.hdr.dime.dim(5)=1; - nii.img = hyperCO2ttccbf_img(:,:,:,4); - save_nii(nii, fullfile(subj_folder, 'hyperCO2_ttcbf.nii')); - - normalCO2tt_img = normalCO2tt_img.*1000; % *1000:change unit from 's' to 'ms' - nii.hdr.dime.dim(5)=1; - nii.img = normalCO2tt_img; - save_nii(nii, fullfile(subj_folder,'normalCO2_tt.nii')); - hyperCO2tt_img = hyperCO2tt_img.*1000; - nii.hdr.dime.dim(5)=1; - nii.img = hyperCO2tt_img; - save_nii(nii, fullfile(subj_folder,'hyperCO2_tt.nii')); - - - % % save change and change% of CBF, tt and ttCBF between BL and HC - diff_cbf_img = hyperCO2cbf_img(:,:,:,4) - normalCO2cbf_img(:,:,:,4); - nii.img = diff_cbf_img; - save_nii(nii, fullfile(subj_folder, 'diff_cbf.nii')); - - cbf_ratio_img = (diff_cbf_img./normalCO2cbf_img(:,:,:,4)).*100; - nii.img = cbf_ratio_img; - save_nii(nii, fullfile(subj_folder, 'change%_cbf.nii')); - - diff_ttccbf_img = hyperCO2ttccbf_img(:,:,:,4) - normalCO2ttccbf_img(:,:,:,4); - nii.img = diff_ttccbf_img; - save_nii(nii, fullfile(subj_folder, 'diff_ttcbf.nii')); - - ttccbf_ratio_img = (diff_ttccbf_img./normalCO2ttccbf_img(:,:,:,4)).*100; - nii.img = ttccbf_ratio_img; - save_nii(nii, fullfile(subj_folder, 'change%_ttcbf.nii')); - - diff_tt_img = normalCO2tt_img - hyperCO2tt_img; - nii.img = diff_tt_img; - save_nii(nii, fullfile(subj_folder, 'diff_tt.nii')); - - tt_ratio_img = (diff_tt_img./normalCO2tt_img).*100; - nii.img = tt_ratio_img; - save_nii(nii, fullfile(subj_folder, 'change%_tt.nii')); -end ->>>>>>> 4d0410692edfc486821256d898544b48d4c43d84