[51] | 1 | function calc_asl_cbf(DATA,ROI,AddInfo) |
---|
| 2 | % CALC_ASL_CBF - Calculate CBF maps from ASL data using T1-map data if |
---|
| 3 | % possible (Aedes plugin) |
---|
| 4 | % |
---|
| 5 | % |
---|
| 6 | % Synopsis: |
---|
| 7 | % |
---|
| 8 | % Description: |
---|
| 9 | % |
---|
| 10 | % Examples: |
---|
| 11 | % |
---|
| 12 | % See also: |
---|
| 13 | % |
---|
| 14 | |
---|
| 15 | % This function is a part of Aedes - A graphical tool for analyzing |
---|
| 16 | % medical images |
---|
| 17 | % |
---|
| 18 | % Copyright (C) 2006 Juha-Pekka Niskanen <Juha-Pekka.Niskanen@uku.fi> |
---|
| 19 | % |
---|
| 20 | % Department of Physics, Department of Neurobiology |
---|
| 21 | % University of Kuopio, FINLAND |
---|
| 22 | % |
---|
| 23 | % This program may be used under the terms of the GNU General Public |
---|
| 24 | % License version 2.0 as published by the Free Software Foundation |
---|
| 25 | % and appearing in the file LICENSE.TXT included in the packaging of |
---|
| 26 | % this program. |
---|
| 27 | % |
---|
| 28 | % This program is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE |
---|
| 29 | % WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. |
---|
| 30 | |
---|
| 31 | % Defaults |
---|
| 32 | lambda = 0.9; |
---|
| 33 | t1map = 1.2; |
---|
| 34 | |
---|
[52] | 35 | if AddInfo.isDataMixed |
---|
[51] | 36 | % Check that the number of images is OK |
---|
| 37 | if rem(length(DATA),2)~=0 || length(DATA)==1 |
---|
| 38 | h=errordlg('The number of input images has to be divisible by 2!',... |
---|
| 39 | 'Error','modal'); |
---|
| 40 | return |
---|
| 41 | end |
---|
| 42 | controlIm = zeros([size(DATA{1}.FTDATA) length(DATA)/2]); |
---|
| 43 | labelIm = zeros([size(DATA{1}.FTDATA) length(DATA)/2]); |
---|
| 44 | for ii=1:2:size(controlIm,3) |
---|
| 45 | controlIm(:,:,ii) = DATA{ii}.FTDATA; |
---|
| 46 | labelIm(:,:,ii) = DATA{ii+1}.FTDATA; |
---|
| 47 | end |
---|
| 48 | for ii=1:length(DATA) |
---|
| 49 | ParentFileName{ii} = fullfile(DATA{ii}.HDR.fpath,DATA{ii}.HDR.fname); |
---|
| 50 | end |
---|
| 51 | else |
---|
[184] | 52 | |
---|
| 53 | if strcmpi(DATA{1}.DataFormat,'bruker_reco') |
---|
| 54 | % Assume that Bruker data is 4D |
---|
| 55 | if rem(size(DATA{1}.FTDATA,4),2)~=0 || size(DATA{1}.FTDATA,4)==1 |
---|
| 56 | h=errordlg('The number of input images has to be divisible by 2!',... |
---|
| 57 | 'Error','modal'); |
---|
| 58 | return |
---|
| 59 | end |
---|
| 60 | |
---|
| 61 | controlIm = squeeze(DATA{1}.FTDATA(:,:,2,:)); |
---|
| 62 | labelIm = squeeze(DATA{1}.FTDATA(:,:,1,:)); |
---|
| 63 | ParentFileName = fullfile(DATA{1}.HDR.fpath,DATA{1}.HDR.fname); |
---|
| 64 | else |
---|
| 65 | % Check that the number of images is OK |
---|
| 66 | if rem(size(DATA{1}.FTDATA,3),2)~=0 || size(DATA{1}.FTDATA,3)==1 |
---|
| 67 | h=errordlg('The number of input images has to be divisible by 2!',... |
---|
| 68 | 'Error','modal'); |
---|
| 69 | return |
---|
| 70 | end |
---|
| 71 | if isfield(DATA{1}.PROCPAR,'ltype') && ... |
---|
| 72 | length(DATA{1}.PROCPAR.ltype)==size(DATA{1}.FTDATA,3) |
---|
| 73 | control_ind = find(strcmpi(DATA{1}.PROCPAR.ltype,'c')); |
---|
| 74 | label_ind = find(strcmpi(DATA{1}.PROCPAR.ltype,'l')); |
---|
| 75 | controlIm = DATA{1}.FTDATA(:,:,control_ind); |
---|
| 76 | labelIm = DATA{1}.FTDATA(:,:,label_ind); |
---|
| 77 | ParentFileName = fullfile(DATA{1}.HDR.fpath,DATA{1}.HDR.fname); |
---|
| 78 | else |
---|
| 79 | controlIm = DATA{1}.FTDATA(:,:,1:2:size(DATA{1}.FTDATA,3)); |
---|
| 80 | labelIm = DATA{1}.FTDATA(:,:,2:2:size(DATA{1}.FTDATA,3)); |
---|
| 81 | ParentFileName = fullfile(DATA{1}.HDR.fpath,DATA{1}.HDR.fname); |
---|
| 82 | end |
---|
| 83 | end |
---|
[51] | 84 | end |
---|
| 85 | |
---|
| 86 | % Ask for lambda |
---|
[80] | 87 | resp = aedes_inputdlg('Lambda coefficient?','Lambda coefficient?',... |
---|
[51] | 88 | num2str(lambda)); |
---|
| 89 | if isempty(resp) |
---|
| 90 | % Canceled |
---|
| 91 | return |
---|
| 92 | end |
---|
| 93 | lambda = str2num(resp{1}); |
---|
| 94 | |
---|
| 95 | % Ask for T1-map |
---|
| 96 | default_path = DATA{1}.HDR.fpath; |
---|
| 97 | if isempty(default_path) |
---|
| 98 | default_path = [pwd,filesep]; |
---|
| 99 | end |
---|
| 100 | |
---|
| 101 | resp = questdlg(['Select if you want to use a constant value for T1 ',... |
---|
| 102 | 'or an existing T1-map.'],... |
---|
| 103 | 'Use constant T1 value or T1-map?',... |
---|
| 104 | 'Use T1-map','Use constant T1','Cancel',... |
---|
| 105 | 'Use T1-map'); |
---|
| 106 | if isempty(resp) || strcmpi(resp,'Cancel') |
---|
| 107 | % Canceled |
---|
| 108 | return |
---|
[52] | 109 | elseif strcmpi(resp,'Use constant T1') |
---|
[51] | 110 | % Ask for T1 value |
---|
[80] | 111 | resp = aedes_inputdlg('Input T1 value','Input T1 value',... |
---|
[51] | 112 | num2str(t1map)); |
---|
| 113 | if isempty(resp) |
---|
| 114 | % Canceled |
---|
| 115 | return |
---|
| 116 | end |
---|
| 117 | fname = ''; |
---|
| 118 | fpath = ''; |
---|
| 119 | t1map = str2num(resp{1}); |
---|
| 120 | t1val = t1map; |
---|
| 121 | else |
---|
| 122 | [fname,fpath,findex] = uigetfile({'*.t1','T1-Files (*.t1)';... |
---|
| 123 | '*.*','All Files (*.*)'},'Select T1-map to open',... |
---|
| 124 | default_path); |
---|
| 125 | if isequal(fname,0) |
---|
| 126 | % Canceled |
---|
| 127 | return |
---|
| 128 | end |
---|
| 129 | try |
---|
[52] | 130 | tmp=load(fullfile(fpath,fname),'-mat'); |
---|
| 131 | t1map = tmp.Data; |
---|
[51] | 132 | catch |
---|
| 133 | h=errordlg({'Could not open T1-map from file',... |
---|
[52] | 134 | fullfile(fpath,fname)},'Error','modal'); |
---|
[51] | 135 | return |
---|
| 136 | end |
---|
| 137 | |
---|
| 138 | |
---|
| 139 | % Check that the T1-map size matches with ASL. If it doesn't, resize it to |
---|
| 140 | % match... |
---|
| 141 | t1map = t1map(:,:,1); |
---|
| 142 | if ~isequal(size(t1map),[size(controlIm,1),size(controlIm,2)]) |
---|
| 143 | t1map = imresize(t1map,[size(controlIm,1),size(controlIm,2)]); |
---|
| 144 | end |
---|
| 145 | t1val = []; |
---|
| 146 | end |
---|
| 147 | |
---|
[185] | 148 | % Ask if the user wants to calculate a single asl map or time series |
---|
| 149 | resp = questdlg(['Calculate an average CBF map (recommended) ',... |
---|
| 150 | 'or a CBF time series (i.e. one CBF map / image pair)?'],... |
---|
| 151 | 'Calculate an average CBF map?',... |
---|
| 152 | 'Average CBF','CBF map / image pair','Cancel',... |
---|
| 153 | 'Average CBF'); |
---|
| 154 | if isempty(resp) || strcmpi(resp,'Cancel') |
---|
| 155 | % Canceled |
---|
| 156 | return |
---|
| 157 | end |
---|
| 158 | if strcmpi(resp,'Average CBF') |
---|
| 159 | CBFts = false; |
---|
| 160 | else |
---|
| 161 | CBFts = true; |
---|
| 162 | end |
---|
| 163 | |
---|
[52] | 164 | % Convert milliseconds to seconds |
---|
[129] | 165 | %t1map=t1map./1000; |
---|
[52] | 166 | |
---|
[185] | 167 | if CBFts |
---|
| 168 | % Calculate differences and CBF maps |
---|
| 169 | asl_map = zeros(size(controlIm)); |
---|
| 170 | cbf_map = zeros(size(controlIm)); |
---|
| 171 | for ii=1:size(controlIm,3) |
---|
| 172 | asl_map(:,:,ii) = controlIm(:,:,ii)-labelIm(:,:,ii); |
---|
| 173 | cbf_map(:,:,ii) = (100.*60.*lambda.*asl_map(:,:,ii))./(t1map.*2.*controlIm(:,:,ii)); |
---|
| 174 | end |
---|
[51] | 175 | |
---|
[185] | 176 | else |
---|
| 177 | % Calculate difference fotr ASL |
---|
| 178 | mean_control = mean(controlIm,3); |
---|
| 179 | mean_label = mean(labelIm,3); |
---|
| 180 | asl_map = mean_control-mean_label; |
---|
| 181 | |
---|
| 182 | % Calculate CBF and save it into a MAT-file |
---|
| 183 | cbf_map = (100.*60.*lambda.*asl_map)./(t1map.*2.*mean_control); |
---|
| 184 | |
---|
| 185 | end |
---|
[51] | 186 | cbf_map(isinf(cbf_map))=0; |
---|
| 187 | cbf_map(isnan(cbf_map))=0; |
---|
| 188 | |
---|
| 189 | % Ask where to save the resulting CBF map |
---|
| 190 | [fn,fp,fi] = uiputfile({'*.mat','Matlab MAT-Files (*.mat)';... |
---|
[185] | 191 | '*.*','All Files (*.*)'},'Save CBF-map as...',... |
---|
| 192 | [default_path,'asl_cbf_map.mat']); |
---|
[51] | 193 | if isequal(fname,0) |
---|
[185] | 194 | % Canceled |
---|
| 195 | return |
---|
[51] | 196 | end |
---|
| 197 | |
---|
| 198 | % Construct the map structure |
---|
| 199 | Data = cbf_map; |
---|
| 200 | Param.ParentFileName = ParentFileName; |
---|
[52] | 201 | Param.T1FileName = fullfile(fpath,fname); |
---|
[51] | 202 | Param.T1value = t1val; |
---|
| 203 | Param.Lambda = lambda; |
---|
| 204 | |
---|
| 205 | % Save the map |
---|
| 206 | try |
---|
[52] | 207 | save(fullfile(fp,fn),'Data','Param','-mat'); |
---|
[51] | 208 | catch |
---|
| 209 | h=errordlg({'Could not save CBF-map to',... |
---|
[52] | 210 | fullfile(fp,fn)},'Error','modal'); |
---|
[51] | 211 | return |
---|
| 212 | end |
---|
| 213 | |
---|
[184] | 214 | % Open in a new Aedes window |
---|
| 215 | aedes(fullfile(fp,fn)) |
---|
[51] | 216 | |
---|
| 217 | |
---|
| 218 | |
---|
| 219 | |
---|
| 220 | |
---|
| 221 | |
---|
| 222 | |
---|
| 223 | |
---|
[184] | 224 | |
---|
| 225 | |
---|