source: misclib/fmri_smooth.m

Last change on this file was 164, checked in by tjniskan, 8 years ago
  • Added option for using external regressors in misclib/fmri_analysis.m

M misclib/fmri_smooth.m
M misclib/fmri_analysis.m
M aedes_revision.m

File size: 2.7 KB
Line 
1function smooth_data_out = fmri_smooth(data,fwhm_sz,voxsize,out_file)
2
3if nargin==4
4  writeSmoothedData = true;
5else
6  writeSmoothedData = false;
7end
8
9if nargin<3 || isempty(voxsize)
10  fwhmInPixels = true;
11else
12  fwhmInPixels = false;
13end
14
15% Check if data is an Aedes structure
16if isstruct(data)
17  data = data.FTDATA;
18elseif ischar(data)
19  data = aedes_data_read(data);
20  data = data.FTDATA;
21end
22
23UseImFilter = true;
24
25smooth_data = zeros(size(data));
26
27% Check data dimensions
28if any(ndims(data)==[4,3]) && size(data,3)~=1
29  use3Dkernel = true;
30else
31  use3Dkernel = false;
32end
33
34% Calculate standard deviations using FWHM
35if fwhmInPixels
36        stds = fwhm_sz/sqrt(8*log(2));
37else
38  stds = (fwhm_sz/sqrt(8*log(2)))./voxsize;
39end
40
41% Calculate kernel size using STDs
42kernel_sz = round(6*stds);
43if kernel_sz(3)>floor((size(data,3)/2))
44  kernel_sz(3)=floor((size(data,3)-1)/2);
45end
46
47% Construct the smoothing kernel
48if use3Dkernel
49  [x,y,z] = meshgrid(-kernel_sz(2):kernel_sz(2),...
50        -kernel_sz(1):kernel_sz(1),...
51        -kernel_sz(3):kernel_sz(3));
52  s_kernel = exp(-(x).^2/(2*(stds(1)).^2)...
53        -(y).^2/(2*(stds(2)).^2)...
54        -(z).^2/(2*(stds(3)).^2));
55  s_kernel = s_kernel/sum(s_kernel(:));
56else
57  [x,y] = meshgrid(-kernel_sz(2):kernel_sz(2),...
58        -kernel_sz(1):kernel_sz(1));
59  s_kernel = exp(-(x).^2/(2*(stds(1)).^2)...
60        -(y).^2/(2*(stds(2)).^2));
61  s_kernel = s_kernel/sum(s_kernel(:));
62end
63
64% Smooth the image data -----------------------
65nVols = size(data,4);
66if ~UseImFilter
67  tmp_sz = size(data(:,:,:,1));
68  %tmp_sz(3)=tmp_sz(3)+fwhm_sz(3);
69  %ind_hi=ceil(tmp_sz/2)+floor(size(s_kernel)/2);
70  %ind_lo=ind_hi-(size(s_kernel)-1);
71  %ind_lo(ind_lo<1)=1;
72  %tmp_kernel = zeros(tmp_sz);
73  %tmp_kernel((ind_lo(1):ind_hi(1))+1,(ind_lo(2):ind_hi(2))+1,...
74  %  (ind_lo(3):ind_hi(3)))=s_kernel;
75        tmp_kernel = s_kernel;
76        tmp_kernel(tmp_sz(1),tmp_sz(2),tmp_sz(3))=0; % Pad with zeros
77end
78for ii=1:nVols
79  if ii==1
80    fprintf(1,'Smoothing volume %d/%d',ii,nVols);
81    bsz = length(sprintf('%d/%d',ii,nVols));
82  else
83    fprintf(1,repmat('\b',1,bsz));
84    bsz = length(sprintf('%d/%d',ii,nVols));
85    fprintf(1,'%d/%d',ii,nVols);
86  end
87  if UseImFilter
88    tmp_data = data(:,:,:,ii);
89    smooth_data(:,:,:,ii) = imfilter(tmp_data,s_kernel);
90  else
91    tmp_data = double(data(:,:,:,ii));
92    %pad_sz=size(tmp_data)+fwhm_sz*2;
93    %tmp_data(pad_sz(1),pad_sz(2),pad_sz(3))=0;
94    %tmp_data = cat(3,tmp_data(:,:,1),tmp_data,tmp_data(:,:,end));
95    tmp_smooth = fftshift(ifftn(fftn(tmp_data).*fftn(tmp_kernel)),3);
96    smooth_data(:,:,:,ii) = tmp_smooth;%tmp_smooth(:,:,1:size(data,3));
97  end
98end
99fprintf(1,'\n');
100
101if nargout>0
102  smooth_data_out = smooth_data;
103end
104
105% Write output files
106if writeSmoothedData
107  aedes_write_nifti(smooth_data,out_file)
108end
Note: See TracBrowser for help on using the repository browser.

Powered by Trac 1.0.9.Copyright © Juha-Pekka Niskanen 2008