Homec4science
Diffusion Slepians (master) Browse

Slepi/Code/Utilities/DPARSF_V2.3_130615/Subfunctions/y_VoxelSpecificHeadMotion.m
22f22bf83060master

Last Change
Raw File

y_VoxelSpecificHeadMotion.m
View Options

function [MeanTDvox, MeanFDvox, Header_Out] = y_VoxelSpecificHeadMotion(RealignmentParameterFile,ReferenceImage,OutputDir,GZFlag)
%function [MeanTDvox, MeanFDvox, Header_Out] = y_VoxelSpecificHeadMotion(RealignmentParameterFile,ReferenceImage,OutputDir,GZFlag)
% Generate the voxel specific head motion
% Input:
% RealignmentParameterFile - The realignment parameter file for a given participant generated by SPM. E.g., rp***.txt
% ReferenceImage - The reference image for realignment (usually the first time point (one-pass) or the mean image after an initial motion correction (two-pass))
% OutputDir - The output directory
% GZFlag - 0: write .nii files; 1: write .nii.gz files.
% Output:
% HMvox_X_4DVolume.nii(.gz) - The voxel-specific translation in X
% HMvox_Y_4DVolume.nii(.gz) - The voxel-specific translation in y
% HMvox_Z_4DVolume.nii(.gz) - The voxel-specific translation in z
% TDvox_4DVolume.nii(.gz) - The voxel-specific total displacement relative to the reference image (TDvox)
% FDvox_4DVolume.nii(.gz) - The voxel-specific frame-wise displacement (relative to the previous time point) (FDvox)
% MeanTDvox.nii(.gz) - Temporal mean of TDvox
% MeanFDvox.nii(.gz - Temporal mean of FDvox
%
% MeanTDvox - Temporal mean of TDvox, output as matrix and could be written to file.
% MeanFDvox - Temporal mean of FDvox, output as matrix and could be written to file.
% Header_Out - The NIfTI Header
%-----------------------------------------------------------
% Written by YAN Chao-Gan 120216.
% The Nathan Kline Institute for Psychiatric Research, 140 Old Orangeburg Road, Orangeburg, NY 10962, USA
% Child Mind Institute, 445 Park Avenue, New York, NY 10022, USA
% The Phyllis Green and Randolph Cowen Institute for Pediatric Neuroscience, New York University Child Study Center, New York, NY 10016, USA
% ycg.yan@gmail.com
if ~exist('GZFlag','var')
GZFlag = 0;
end
fprintf(['\nGenerating voxel specific head motion for %s...'], RealignmentParameterFile);
RP=load(RealignmentParameterFile);
nTimePoint=size(RP,1);
sinq1=sin(RP(:,4));
sinq2=sin(RP(:,5));
sinq3=sin(RP(:,6));
cosq1=cos(RP(:,4));
cosq2=cos(RP(:,5));
cosq3=cos(RP(:,6));
[Data0 Head0]=rest_ReadNiftiImage(ReferenceImage);
[nDim1,nDim2,nDim3]=size(Data0);
[D1,D2,D3] = ndgrid(1:nDim1,1:nDim2,1:nDim3);
D=zeros(4,nDim1,nDim2,nDim3);
D(1,:,:,:)=D1;
D(2,:,:,:)=D2;
D(3,:,:,:)=D3;
D(4,:,:,:)=ones(nDim1,nDim2,nDim3);
HMvox_X = zeros(nDim1,nDim2,nDim3,nTimePoint);
HMvox_Y = zeros(nDim1,nDim2,nDim3,nTimePoint);
HMvox_Z = zeros(nDim1,nDim2,nDim3,nTimePoint);
for t=1:nTimePoint
M1=[1 0 0 0;...
0 cosq1(t) sinq1(t) 0;...
0 -sinq1(t) cosq1(t) 0;...
0 0 0 1;];
M2=[cosq2(t) 0 sinq2(t) 0;...
0 1 0 0;...
-sinq2(t) 0 cosq2(t) 0;...
0 0 0 1;];
M3=[cosq3(t) sinq3(t) 0 0;...
-sinq3(t) cosq3(t) 0 0;...
0 0 1 0;...
0 0 0 1;];
MT=[1 0 0 RP(t,1);...
0 1 0 RP(t,2);...
0 0 1 RP(t,3);...
0 0 0 1;];
M_RigidBodyTransform=MT*M1*M2*M3;
Disp = gmdmp(inv(M_RigidBodyTransform)*Head0.mat, 2, D, 1) - gmdmp(Head0.mat, 2, D, 1);
%gmdmp.m is a program for General Multi Dimensional Matrix Product written by Wynton Moore
HMvox_X(:,:,:,t) = squeeze(Disp(1,:,:,:));
HMvox_Y(:,:,:,t) = squeeze(Disp(2,:,:,:));
HMvox_Z(:,:,:,t) = squeeze(Disp(3,:,:,:));
fprintf('.');
end
Header_Out = Head0;
Header_Out.pinfo = [1;0;0];
Header_Out.dt =[16,0];
Header_Out.dim = [nDim1 nDim2 nDim3 nTimePoint];
if ~(7==exist(OutputDir,'dir'))
mkdir(OutputDir);
end
HMvox_X=single(HMvox_X);
HMvox_Y=single(HMvox_Y);
HMvox_Z=single(HMvox_Z);
y_Write4DNIfTI(HMvox_X,Header_Out,[OutputDir,filesep,'HMvox_X_4DVolume.nii']);
y_Write4DNIfTI(HMvox_Y,Header_Out,[OutputDir,filesep,'HMvox_Y_4DVolume.nii']);
y_Write4DNIfTI(HMvox_Z,Header_Out,[OutputDir,filesep,'HMvox_Z_4DVolume.nii']);
TDvox = sqrt(HMvox_X.^2 + HMvox_Y.^2 + HMvox_Z.^2);
y_Write4DNIfTI(TDvox,Header_Out,[OutputDir,filesep,'TDvox_4DVolume.nii']);
MeanTDvox = mean(TDvox,4);
y_Write4DNIfTI(MeanTDvox,Header_Out,[OutputDir,filesep,'MeanTDvox.nii']);
clear TDvox
FDvox = sqrt((HMvox_X(:,:,:,2:end) - HMvox_X(:,:,:,1:end-1)).^2 + (HMvox_Y(:,:,:,2:end) - HMvox_Y(:,:,:,1:end-1)).^2 + (HMvox_Z(:,:,:,2:end) - HMvox_Z(:,:,:,1:end-1)).^2);
FDvox = cat(4, zeros(nDim1,nDim2,nDim3), FDvox);
y_Write4DNIfTI(FDvox,Header_Out,[OutputDir,filesep,'FDvox_4DVolume.nii']);
MeanFDvox=mean(FDvox,4);
y_Write4DNIfTI(MeanFDvox,Header_Out,[OutputDir,filesep,'MeanFDvox.nii']);
if GZFlag
gzip([OutputDir,filesep,'HMvox_X_4DVolume.nii']);
delete([OutputDir,filesep,'HMvox_X_4DVolume.nii']);
gzip([OutputDir,filesep,'HMvox_Y_4DVolume.nii']);
delete([OutputDir,filesep,'HMvox_Y_4DVolume.nii']);
gzip([OutputDir,filesep,'HMvox_Z_4DVolume.nii']);
delete([OutputDir,filesep,'HMvox_Z_4DVolume.nii']);
gzip([OutputDir,filesep,'TDvox_4DVolume.nii']);
delete([OutputDir,filesep,'TDvox_4DVolume.nii']);
gzip([OutputDir,filesep,'FDvox_4DVolume.nii']);
delete([OutputDir,filesep,'FDvox_4DVolume.nii']);
gzip([OutputDir,filesep,'MeanTDvox.nii']);
delete([OutputDir,filesep,'MeanTDvox.nii']);
gzip([OutputDir,filesep,'MeanFDvox.nii']);
delete([OutputDir,filesep,'MeanFDvox.nii']);
end
fprintf(['\nGenerating voxel specific head motion for %s finished.\n'], RealignmentParameterFile);
c4science · Help