source: misclib/fmri_filter.m @ 61

Last change on this file since 61 was 50, checked in by tjniskan, 11 years ago
  • Fixed a typo in an2_smiswrite.m
  • Added some fmri-specific functions under misclib. These do not

affect Aedes at all.

M an2_smiswrite.m
A misclib/fmri_smooth.m
A misclib/fmri_filter.m
A misclib/fmri_corr.m
M an2_revision.m

File size: 2.1 KB
Line 
1function filtered_data_out = fmri_filter(data,TR,varargin)
2
3% Defaults
4hipass = [];
5lowpass = [];
6detrending = false;
7reserve_mean = true;
8out_file = '';
9
10% Parse varargin
11for ii=1:2:length(varargin)
12  param = varargin{ii};
13  value = varargin{ii+1};
14  switch lower(param)
15        case 'hipass'
16          hipass = value;
17        case 'lowpass'
18          lowpass = value;
19        case 'detrending'
20          if strcmpi(value,'on')
21                detrending = true;
22          end
23        case 'reservemean'
24          if strcmpi(value,'off')
25                reserve_mean = false;
26          end
27        case 'filename'
28          out_file = value;
29        otherwise
30          error('Unknown parameter "%s"',param)
31  end
32end
33
34if isstruct(data)
35  data = data.FTDATA;
36elseif ischar(data)
37  data = an2_data_read(data);
38  data = data.FTDATA;
39end
40filtered_data = double(data);
41
42if isempty(hipass) && isempty(lowpass) && ~detrend
43  return
44end
45
46% Construct hi-pass filter
47if ~isempty(hipass)
48  [B_hipass,A_hipass]=butter(4,hipass/((1/TR)/2),'high');
49end
50
51% Construct low-pass filter
52if ~isempty(lowpass)
53  [B_lowpass,A_lowpass]=butter(15,lowpass/((1/TR)/2),'low');
54end
55
56
57% Filter image data -----------------------------
58if ~isempty(lowpass) || ~isempty(hipass) || detrending
59 
60  do_lp = true;
61  do_hp = true;
62  if isempty(lowpass)
63        do_lp = false;
64  end
65 
66  if isempty(hipass)
67        do_hp = false;
68  end
69 
70  % Waitbar
71  wbh=an2_wbar(0,'Filtering image data...');
72 
73  counter=1;
74  sz = ones(1,4);
75  sz(1:length(size(data)))=size(data);
76  nVox = prod(sz(1:3));
77  for ii=1:size(data,1)
78        for kk=1:size(data,2)
79          for zz=1:size(data,3)
80                tmp = squeeze(double(data(ii,kk,zz,:)));
81                mean_tmp = mean(tmp);
82                if do_hp
83                  tmp = filtfilt(B_hipass,A_hipass,tmp);
84                  if reserve_mean
85                        tmp = tmp+mean_tmp;
86                  end
87                end
88                if do_lp
89                  tmp = filtfilt(B_lowpass,A_lowpass,tmp);
90                end
91                if detrending
92                  tmp = detrend(tmp,'linear');
93                  if reserve_mean
94                        tmp = tmp+mean_tmp;
95                  end
96                end
97                filtered_data(ii,kk,zz,:) = tmp;
98                an2_wbar(counter/nVox,wbh)
99                counter=counter+1;
100          end
101        end
102  end
103  close(wbh)
104end
105
106if nargout>0
107  filtered_data_out = filtered_data;
108end
109
110% Write output to file
111if ~isempty(out_file)
112  an2_write_nifti(filtered_data,out_file);
113end
114
115
116
117
118
Note: See TracBrowser for help on using the repository browser.

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