parent_path = "E:\Projects\CBIrep_Imaging\subjects"; cd (parent_path) addpath("..\scripts\toolbox\DVARS-master\") addpath("..\scripts\toolbox\DVARS-master\Nifti_Util\") addpath("..\scripts") addpath("..\scripts\toolbox\spm12") addpath("..\scripts\toolbox\fMRI-Quality-Checker-master") addpath("..\scripts\toolbox\spm12\matlabbatch") addpath("..\scripts\toolbox\nifti_utils-master\") %% Preps % px = dir('sub*'); BIDS_path = "E:\Projects\CBIrep_Imaging\subjects"; % if ~exist(BIDS_path,'dir') % mkdir(BIDS_path) % end % physio_files = "F:\Projects\StrokePT_tDCS_rsfmri\physio"; % cd (physio_files) % phx = dir("SB*"); % timeX = {'tpA','tpB','tpC','tpD'}; % master = zeros300,5}; % noPhysio = zeros(300); % master_tag = {'subject','group','timepoint','dataset','physio','volumes','bad_volumes','meanDeltaDvar'}; tick = 1; tack = 1; % Check for registration quality (user input)? quali = "off"; %% level1: group folders cd (parent_path) gx = dir("sub*"); for g = 4 group = gx(g).name; cd (fullfile(gx(g).folder,gx(g).name)) %% level2: timepoints tx = dir('tp*'); for t = 1:length(tx) time = tx(t).name; cd (fullfile(tx(t).folder,tx(t).name)) session = strcat('ses-',num2str(find(strcmp(time,timeX)))); tix = dir("SB*") %% level3: individual datasets for m = 1:length(tix) cd (fullfile(parent_path,group,time)) dataset = tix(m).name; %% find corresponding Physio dataset physioID = find(~cellfun('isempty',regexp({phx.name},dataset,'once'))); if ~isempty(physioID) physiofolder = fullfile(physio_files,phx(physioID).name); master{tick,5} = "y"; else noPhysio{tack} = dataset; tack = tack +1; master{tick,5} = "n"; end %% sjx = strfind(dataset,'_'); subject = dataset(1:sjx(1)-1); master {tick,1} = subject; master {tick,2} = group; master {tick,3} = num2str(t); master {tick,4} = dataset; disp(strcat("running dataset:",dataset)); cd(dataset) %% copy data to BIDS folder_from = pwd; folder_to = fullfile(BIDS_path,strcat("sub-",subject(3:end)),session); % % % % % % % % if ~exist(folder_to,"dir") % % % % mkdir(folder_to) % % % % end % % % % % copy MRI data % % % % copyfile(folder_from,folder_to) % % % % % % % % % copy physio data % % % % if exist("physiofolder","var") % % % % physio_to = fullfile(folder_to,"physio",phx(physioID).name); % % % % if ~exist(physio_to,"dir") % % % % mkdir(physio_to) % % % % end % % % % copyfile(physiofolder,physio_to); % % % % end %% cd (folder_to) % movefile("fMRI","func") % movefile("T2w","anat") % a) get structural if exist(fullfile(folder_to,"anat/"))==7 && exist (fullfile(folder_to,"func/"))==7 cd (fullfile(folder_to,"anat/")) struct_orig = dir("*.1.nii.gz"); struct_name = struct_orig(length(struct_orig)).name; struct_fullfile = fullfile(struct_orig(length(struct_orig)).folder,struct_orig(length(struct_orig)).name); struct_reg = dir("*BiasBet_AnnorsfMRI.nii.gz"); if ~isempty(struct_reg) structReg_fullfile = fullfile(struct_reg.folder,struct_reg.name); struct_BET = dir("*BiasBet.nii.gz"); structBET_fullfile = fullfile(struct_BET.folder,struct_BET.name); end % copyfile(struct_fullfile,fullfile(folder_to,"anat/",strcat(struct_name(1:end-7),'_T2w.nii.gz'))) struct_orig = dir("*_T2w.nii.gz"); struct_name = struct_orig(length(struct_orig)).name; struct_fullfile = fullfile(struct_orig(length(struct_orig)).folder,struct_orig(length(struct_orig)).name); % b) get fMRI cd (fullfile(folder_to,"func/")) fmri_orig = dir("*.1.nii.gz"); clear ('fx_rep') for jx = 1:length(fmri_orig) fx_scan = niftiinfo(fmri_orig(jx).name); fx_rep(jx) = fx_scan.ImageSize(4); end fx_id = find(fx_rep>100); fx = length(fx_id); % % % % % % % % if fx>1 % % % % fmri_conc = niftiread(fmri_orig(fx_id(1)).name); % % % % info = niftiinfo(fmri_orig(fx_id(1)).name); % % % % info.ImageSize = [96,96,16,(fx*105-((fx-1)*5))]; % % % % for imgx = 2:fx % % % % fmri_add = niftiread(fmri_orig(fx_id(1)).name); % % % % fmri_conc = cat(4,fmri_conc,fmri_add(:,:,:,6:end)); % % % % end % % % % conc_fmri = strcat(dataset,'.conc.nii'); % % % % niftiwrite(fmri_conc,conc_fmri,info); % % % % gzip(conc_fmri) % % % % fMRI_fullfile = strcat(fullfile(folder_to,'func',conc_fmri),".gz"); % % % % fMRI_file = dir('*conc.nii.gz'); % % % % elseif fx==1 % % % % fMRI_fullfile = fullfile(fmri_orig(fx_id).folder,fmri_orig(fx_id).name); % % % % fMRI_file = fmri_orig(fx_id); % % % % end fMRI_fullfile = fullfile(fmri_orig(fx_id).folder,fmri_orig(fx_id).name); fMRI_file = fmri_orig(fx_id); fMRI_name = fMRI_file.name; % copyfile(fMRI_fullfile,fullfile(fMRI_file.folder,strcat(fMRI_name(1:end-7),'_bold.nii.gz'))) info = niftiinfo(fMRI_fullfile); master{tick,6} = info.ImageSize(4); fmri_smooth = dir("*SmoothBet.nii.gz"); if ~isempty(fmri_smooth) fmriSmooth_fullfile = fullfile(fmri_smooth.folder,fmri_smooth.name); end % if ~exist("analyses","dir") % mkdir("analyses") % end % cd ("analyses") % % % % %% get DVARS % % % % V1 = load_untouch_nii(fMRI_file.name); % % % % V2 = V1.img; % % % % X0 = size(V2,1); Y0 = size(V2,2); Z0 = size(V2,3); T0 = size(V2,4); % % % % I0 = prod([X0,Y0,Z0]); % % % % Y = reshape(V2,[I0,T0]); clear V2 V1; % % % % [DVARS,DVARS_Stat]=DVARSCalc(Y,'scale',1/100,'TransPower',1/3,'RDVARS','verbose',1); % % % % [V,DSE_Stat]=DSEvars(Y,'scale',1/100); % % % % % fMRIDiag_plot(V,DVARS_Stat) % % % % % saveas(gcf,strcat(fmri_name(1:end-7),'_',"DVARS_plot.png")) % % % % % pause(2) % % % % % close all % % % % % % % % bad_points = find(DVARS_Stat.pvals<0.001); % % % % bad_points (bad_points<5)=[]; % % % % perc_bad = length(bad_points)/(info.ImageSize(4)-5)*100; % % % % disp("amount of bad points: ") % % % % disp(num2str(perc_bad)) % % % % disp("%") % % % % % % % % % master{tick,6} = fmri_name; % % % % % bad_points_list{tick,2} = bad_points; % % % % master{tick,8} = mean(DSE_Stat.DeltapDvar(6:end),'omitnan'); % % % % % bad_points_list{tick,5} = mean(DSE_Stat.DeltapSvar(6:end)); % % % % master{tick,7} = perc_bad; % c) get Regression File cd (fullfile(folder_to,"func/")) ABAreg_orig = dir("*SmoothBet_AnnoSplit_rsfMRI.nii.gz"); if ~isempty(ABAreg_orig) ABAreg_fullfile = fullfile(ABAreg_orig.folder,ABAreg_orig.name); end if exist("regr","dir") cd regr\ SFRGR_orig = dir("*_SFRGR.nii.gz"); if ~isempty(SFRGR_orig) SFRGR_fullfile = fullfile(SFRGR_orig.folder,SFRGR_orig.name); else SFRGR_fullfile = []; end end %% plot raw images to check registration if ~exist("quali",'var') quali = "off" end if quali == "on" fMRI_imag = niftiread(fmriSmooth_fullfile); fMRI_imag = imrotate(fMRI_imag,-90); if ~isempty(ABAreg_orig) fMRI_Reg = niftiread(ABAreg_fullfile); fMRI_Reg = imrotate(fMRI_Reg,-90); else fMRI_Reg = zeros(96,96,16); end T2raw = niftiread(structBET_fullfile); T2raw = imrotate(T2raw,-90); if ~isempty(struct_reg) T2Reg = niftiread(structReg_fullfile); T2Reg = imrotate(T2Reg,-90); else T2Reg = zeros(256,256,28); end Loc = [8 13]; figure('units','normalized','outerposition',[0 0 1 1]) for ix = 1:2 % Raw T2w subplot(2,4,1+(4*(ix-1))) imagesc(T2raw(:,:,floor(Loc(ix)/16*48))) colormap("hot") % T2 ARA Reg subplot(2,4,2+(4*(ix-1))) imagesc(T2Reg(:,:,floor(Loc(ix)/16*48))) % colormap("hot") % Raw fMRI subplot(2,4,3+(4*(ix-1))) imagesc(fMRI_imag(:,:,Loc(ix))) % colormap("gray") % fMRI ARA Reg subplot(2,4,4+(4*(ix-1))) imagesc(fMRI_Reg(:,:,Loc(ix))) % colormap("hot") end title(strcat(subject,'-',time)) pause close %% user input to evaluate quality RegQuali{tick,1}=subject; RegQuali{tick,2}=time; % get input of image quality checkT2 = input("is T2 ok? Y/N [Y]"); if isempty(checkT2) checkT2 = "Y"; else checkT2 = "N"; end RegQuali{tick,3} = checkT2; checkfMRI = input("is fMRI ok? Y/N [Y]"); if isempty(checkfMRI) checkfMRI = "Y"; else checkfMRI = "N"; end disp('') RegQuali{tick,4} = checkfMRI; end %% tick=tick +1; end end end end % % % a) get structural % % if exist(fullfile(parent_path,subject,"MRI",time,"T2w/"))==7 && exist (fullfile(parent_path,subject,"MRI",time,"fMRI/"))==7 % % cd (fullfile(parent_path,subject,"MRI",time,"T2w/")) % % struct_orig = dir("*.1.nii.gz"); % % struct_name = struct_orig(length(struct_orig)).name; % % struct_fullfile = fullfile(struct_orig(length(struct_orig)).folder,struct_orig(length(struct_orig)).name); % % % % struct_reg = dir("*BiasBet_AnnorsfMRI.nii.gz"); % % if ~isempty(struct_reg) % % structReg_fullfile = fullfile(struct_reg.folder,struct_reg.name); % % end % % struct_BET = dir("*BiasBet.nii.gz"); % % structBET_fullfile = fullfile(struct_BET.folder,struct_BET.name); % % % % % % % b) get fMRI % % cd (fullfile(parent_path,subject,"MRI",time,"fMRI/")) % % fmri_orig = dir("*.1.nii.gz"); % % fmri_name = fmri_orig(length(fmri_orig)).name; % % fMRI_fullfile = fullfile(fmri_orig(length(fmri_orig)).folder,fmri_orig(length(fmri_orig)).name); % % fmri_smooth = dir("*SmoothBet.nii.gz"); % % fmriSmooth_fullfile = fullfile(fmri_smooth.folder,fmri_smooth.name); % % % % if ~exist("analyses","dir") % % mkdir("analyses") % % end % % cd ("analyses") % % % % % c) get Regression File % % cd (fullfile(parent_path,subject,"MRI",time,"fMRI/")) % % ABAreg_orig = dir("*SmoothBet_AnnoSplit_rsfMRI.nii.gz"); % % if ~isempty(ABAreg_orig) % % ABAreg_fullfile = fullfile(ABAreg_orig.folder,ABAreg_orig.name); % % end % % cd regr\ % % SFRGR_orig = dir("*_SFRGR.nii.gz"); % % if ~isempty(SFRGR_orig) % % SFRGR_fullfile = fullfile(SFRGR_orig.folder,SFRGR_orig.name); % % else % % SFRGR_fullfile = []; % % end % % % plot raw images to check registration % % % % % % fMRI_imag = niftiread(fmriSmooth_fullfile); % % fMRI_imag = imrotate(fMRI_imag,-90); % % % % if ~isempty(ABAreg_orig) % % fMRI_Reg = niftiread(ABAreg_fullfile); % % fMRI_Reg = imrotate(fMRI_Reg,-90); % % else % % fMRI_Reg = zeros(96,96,16); % % end % % % % T2raw = niftiread(structBET_fullfile); % % T2raw = imrotate(T2raw,-90); % % if ~isempty(struct_reg) % % T2Reg = niftiread(structReg_fullfile); % % T2Reg = imrotate(T2Reg,-90); % % else % % T2Reg = zeros(256,256,28); % % end % % % % Loc = [8 13]; % % %% % % figure('units','normalized','outerposition',[0 0 1 1]) % % for ix = 1:2 % % % Raw T2w % % subplot(2,4,1+(4*(ix-1))) % % imagesc(T2raw(:,:,floor(Loc(ix)/16*28))) % % colormap("hot") % % % T2 ARA Reg % % subplot(2,4,2+(4*(ix-1))) % % imagesc(T2Reg(:,:,floor(Loc(ix)/16*28))) % % % colormap("hot") % % % Raw fMRI % % subplot(2,4,3+(4*(ix-1))) % % imagesc(fMRI_imag(:,:,Loc(ix))) % % % colormap("gray") % % % fMRI ARA Reg % % subplot(2,4,4+(4*(ix-1))) % % imagesc(fMRI_Reg(:,:,Loc(ix))) % % % colormap("hot") % % end % % title(strcat(subject,'-',time)) % % pause % % close % % %% % % RegQuali{tick,1}=subject; % % RegQuali{tick,2}=time; % % % % % get input of image quality % % checkT2 = input("is T2 ok? Y/N [Y]"); % % if isempty(checkT2) % % checkT2 = "Y"; % % else % % checkT2 = "N"; % % end % % RegQuali{tick,3} = checkT2; % % % % checkfMRI = input("is fMRI ok? Y/N [Y]"); % % if isempty(checkfMRI) % % checkfMRI = "Y"; % % else % % checkfMRI = "N"; % % end % % disp('') % % RegQuali{tick,4} = checkfMRI; % % % % tick = tick + 1;