Here is a function which lets you see the overlap of
different T images, and it prints a short results table to
the matlab figure window. You will have to edit the paths
to tell it where to find the canonical images on your
system, and you have to set the degrees of freedom manually
for each image, but otherwise it should work.
Antonia.
function varargout = show_T(varargin)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% this file plots T images on the canonical T1 hires
%% with a nice gui
%% optional inputs are image files for blobs to show
%% if used without inputs, figure is cleared by default
%% and user is asked for blobs file
spm_defaults
global blobs st
if(nargin>0)
action = varargin{1};
else
action = 'new';
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
switch lower(action)
%----------------------------------------------------------
case 'new'
%% set up blobs array
cols = [1 0 0; 0 1 0; 0 0 1];
for i=1:3
blobs(i).name = [];
blobs(i).colour = cols(i,:);
blobs(i).prob = 0.001;
blobs(i).k = 10;
blobs(i).df = 15;
blobs(i).tails = 1;
end
show_t('setup');
case 'input'
%% read in a blobs array
blobs = varargin{2};
show_t('setup');
case 'setup'
%% set up the gui
F = spm_figure('GetWin','Graphics');
spm_figure('Clear',F)
%% show hires
if(strcmp(spm('ver'),'SPM2'))
hname =
'c:/matlab701/toolbox/spm2/canonical/single_subj_T1.mnc';
else
hname =
'c:/matlab701/toolbox/spm99/canonical/single_subj_T1.img';
end
H = spm_vol(hname);%fullfile(clustroot,clust_stub));
hh=spm_orthviews('Image',hname,[0.05,0.25,0.9,0.9]);
%% set up hires callback
st.callback = ['show_T(''show_coords'',st.centre);'];
h = 20;
r = linspace(150,50,3);
w = [250,50,50,35,35,35,50,35];
c = [20,cumsum(w)];
%% do column headers
uicontrol(F,'Style','Text',...
'String',['Probability'],...
'HorizontalAlignment','Left',...
'Position',[c(3),r(1)+20,w(3)-5,h], ...
'BackgroundColor',[1,1,1]);
uicontrol(F,'Style','Text',...
'String',['Extent'],...
'HorizontalAlignment','Left',...
'Position',[c(4),r(1)+20,w(4)-5,h], ...
'BackgroundColor',[1,1,1]);
uicontrol(F,'Style','Text',...
'String',['DF'],...
'HorizontalAlignment','Left',...
'Position',[c(5),r(1)+20,w(5)-5,h], ...
'BackgroundColor',[1,1,1]);
uicontrol(F,'Style','Text',...
'String',['Tails'],...
'HorizontalAlignment','Left',...
'Position',[c(6),r(1)+20,w(6)-5,h], ...
'BackgroundColor',[1,1,1]);
%% make a gui
for i=1:3
uicontrol(F,'Style','Text',...
'String',[blobs(i).name],...
'Tag',['Filename',num2str(i)],...
'HorizontalAlignment','Left',...
'Position',[c(1),r(i),w(1)-5,h], ...
'BackgroundColor',[1,1,1]);
uicontrol(F,'String','New',...
'Tag',['GetFile',num2str(i)],...
'HorizontalAlignment','Left',...
'Position',[c(2),r(i),w(2)-5,h], ...
'BackgroundColor',[1,1,1],...
'Callback',['show_T(''get_file'',',num2str(i),')']);
uicontrol(F,'Style','Edit',...
'String',[num2str(blobs(i).prob)],...
'Tag',['P',num2str(i)],...
'HorizontalAlignment','Left',...
'Position',[c(3),r(i),w(3)-5,h], ...
'BackgroundColor',[1,1,1]);
uicontrol(F,'Style','Edit',...
'String',[num2str(blobs(i).k)],...
'Tag',['K',num2str(i)],...
'HorizontalAlignment','Left',...
'Position',[c(4),r(i),w(4)-5,h], ...
'BackgroundColor',[1,1,1]);
uicontrol(F,'Style','Edit',...
'String',[num2str(blobs(i).df)],...
'Tag',['DF',num2str(i)],...
'HorizontalAlignment','Left',...
'Position',[c(5),r(i),w(5)-5,h], ...
'BackgroundColor',[1,1,1]);
uicontrol(F,'Style','Edit',...
'String',[num2str(blobs(i).tails)],...
'Tag',['Tail',num2str(i)],...
'HorizontalAlignment','Left',...
'Position',[c(6),r(i),w(6)-5,h], ...
'BackgroundColor',[1,1,1]);
uicontrol(F,'Tag',['colour',num2str(i)],...
'Position',[c(7),r(i),w(7)-5,h], ...
'BackgroundColor',blobs(i).colour);
end
%% add a redraw button
uicontrol(F,'String','Redraw Blobs',...
'Position',[300,200,150,50],...
'FontSize',12,...
'BackgroundColor',[0.8 0.8 1],...
'Callback','show_T(''redraw'')');
%% add a render button
uicontrol(F,'String','Render blobs',...
'Position',[300,250,150,50],...
'FontSize',12,...
'BackgroundColor',[0.8 0.8 1],...
'Callback','show_T(''render'')');
%% add coordinates buttons
uicontrol(F,'Style','Text',...
'String','Current Crosshair Location',...
'BackgroundColor',[1 1 1],...
'Position',[300,420,150,20])
uicontrol(F,'Style','Edit',...
'String','',...
'Tag',['Label_coords'],...
'UserData','',...
'HorizontalAlignment','Right',...
'Position',[300,400,80,20], ...
'BackgroundColor',[1,1,1]);
uicontrol(F,'String','Set',...
'Position',[390,400,50,20],...
'Callback',['spm_orthviews(''Reposition'',str2num(get(findobj(''Tag'',''Label_coords''),''String'')))'])
%% redraw
% show_T('redraw')
F = spm_figure('GetWin','Graphics');
spm_orthviews('rmBlobs',1:24)
%----------------------------------------------------------
case 'get_blobs'
%% read in the blobs data to display or render
F = spm_figure('GetWin','Graphics');
%% check which blobs files are valid
todo = zeros(1,3);
for i=1:3
try
B = spm_vol(blobs(i).name);
if(length(B)>0)
todo(i) = 1;
end
catch
todo(i) = 0;
end
end
%% first read all the values in from the gui
for i=find(todo)
str =
get(findobj(F,'Tag',['P',num2str(i)]),'String');
try
p = str2num(str);
if(p>0 & p<1)
blobs(i).prob = p;
else
disp(['Invalid Probability: ',str])
spm('beep'), return;
end
catch
disp(['Invalid Probability: ',str])
spm('beep'), return;
end
str =
get(findobj(F,'Tag',['K',num2str(i)]),'String');
try
blobs(i).k = str2num(str);
catch
disp(['Invalid Extent threshold: ',str])
spm('beep'), return;
end
str =
get(findobj(F,'Tag',['DF',num2str(i)]),'String');
try
blobs(i).df = str2num(str);
catch
disp(['Invalid Probability: ',str])
spm('beep'), return;
end
str =
get(findobj(F,'Tag',['Tail',num2str(i)]),'String');
try
blobs(i).tails = str2num(str);
catch
disp(['Invalid Tails: ',str])
spm('beep'), return;
end
if ~(blobs(i).tails==1 | blobs(i).tails==2)
disp(['Invalid Tails: ',str])
spm('beep'), return;
end
end
%% load up one image
kk=find(todo);
if(length(kk)>0)
B = spm_vol(blobs(kk(1)).name);
%% set up data structure
for i=1:3
data(i).XYZ = [];
data(i).t = [];
data(i).mat = B.mat;
data(i).dim = B.dim;
data(i).colour = blobs(i).colour;
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%5
%% then draw all the images
% todo
for kk=find(todo)
%% read in image
B = spm_vol(blobs(kk).name);
BV = spm_read_vols(B);
% blobs
%% calc T threshold
t =
abs(tinv(blobs(kk).prob/blobs(kk).tails,blobs(kk).df));
%% apply it
BV3 = BV.*(BV>t);
ccx = []; ccy = []; ccz = []; cci = [];
for i=1:size(BV3,3)
[cx,cy] = find(BV3(:,:,i)>0);
if(length(cx)>0)
ind =
sub2ind(size(BV3),cx,cy,repmat(i,size(cx)));
ci = BV3(ind);
ccx = [ccx;cx];
ccy = [ccy;cy];
ccz = [ccz;ones(size(cx))*i];
cci = [cci;ci];
end
end
locn = [ccx,ccy,ccz];
lev = [cci];
clear coords;
if(length(lev)>0)
% keyboard
a = spm_clusters(locn');
b = hist(a,1:max(a));
good = find(b>blobs(kk).k);
keep = zeros(size(a));
for i=1:length(good)
keep = keep | a==good(i);
nv = sum(a==good(i));
tt = max(lev(a==good(i)));
pxyz = locn(lev==tt & a'==good(i),:);
pxyz = pxyz(1,:);
pmni = pxyz'.*2+B.mat(1:3,4);
coords(i,:) =
[pmni',nv,tt,(1-tcdf(tt,blobs(kk).df))*blobs(kk).tails];
end
blocn{kk} = locn(keep,:);
blev{kk} = lev(keep,:);
disp(['File: ',blobs(kk).name])
disp('MNI_x MNI_y MNI_z size
prob')
disp(coords)
data(kk).XYZ = blocn{kk}';
data(kk).t = blev{kk};
data(kk).colour = blobs(kk).colour;
%spm_orthviews('AddColouredBlobs',F,blocn{kk}',blev{kk}',B.mat,blobs(kk).colour)
end
end
varargout{1} = todo;
varargout{2} = data;
%----------------------------------------------------------
case 'redraw'
F = spm_figure('GetWin','Graphics');
%% clear hires
spm_orthviews('rmBlobs',1:24)
[todo,data] = show_t('get_blobs');
for kk=find(todo)
spm_orthviews('AddColouredBlobs',F,data(kk).XYZ,data(kk).t, ...
data(kk).mat,data(kk).colour)
end
if(sum(todo)==2)
kk=find(todo);
cc = data(kk(1)).colour+data(kk(2)).colour;
cc(cc>1) = 1;
overlap =
intersect(data(kk(1)).XYZ',data(kk(2)).XYZ','rows');
disp(['Files ',num2str(kk(1)),' and
',num2str(kk(2)),' overlap by ',num2str(size(overlap,1)),'
voxels'])
spm_orthviews('AddColouredBlobs',F,overlap',ones(size(overlap,1),1)',data(kk(1)).mat,cc)
end
if(sum(todo)==3)
overlap =
intersect(data(1).XYZ',data(2).XYZ','rows');
disp(['Files 1 and 2 overlap by
',num2str(size(overlap,1)),' voxels'])
cc = data(1).colour+data(2).colour;
cc(cc>1) = 1;
spm_orthviews('AddColouredBlobs',F,overlap',ones(size(overlap,1),1)',data(kk(1)).mat,cc)
overlap =
intersect(data(2).XYZ',data(3).XYZ','rows');
cc = data(2).colour+data(3).colour;
cc(cc>1) = 1;
disp(['Files 2 and 3 overlap by
',num2str(size(overlap,1)),' voxels'])
spm_orthviews('AddColouredBlobs',F,overlap',ones(size(overlap,1),1)',data(kk(1)).mat,cc)
overlap =
intersect(data(1).XYZ',data(3).XYZ','rows');
cc = data(1).colour+data(3).colour;
cc(cc>1) = 1;
disp(['Files 1 and 3 overlap by
',num2str(size(overlap,1)),' voxels'])
spm_orthviews('AddColouredBlobs',F,overlap',ones(size(overlap,1),1)',data(kk(1)).mat,cc)
overlap =
intersect(data(1).XYZ',data(3).XYZ','rows');
overlap = intersect(overlap,data(2).XYZ','rows');
cc = data(1).colour+data(2).colour+data(3).colour;
cc(cc>1) = 1;
disp(['Files 1, 2 and 3 overlap by
',num2str(size(overlap,1)),' voxels'])
spm_orthviews('AddColouredBlobs',F,overlap',ones(size(overlap,1),1)',data(kk(1)).mat,cc)
end
disp('------------------------------------------------')
spm_orthviews('Redraw')
case 'render'
F = spm_figure('GetWin','Graphics');
%% clear hires
spm_orthviews('rmBlobs',1:24)
[todo,data] = show_t('get_blobs');
n = 1;
for kk=find(todo)
dd(n) = data(kk);
n=n+1;
end
% keyboard
brt = 0.5;
rendfile =
'c:\MATLAB701\toolbox\spm99\rend\render_single_subj.mat';
spm_render(dd,brt,rendfile);
case 'get_file'
fno = varargin{2}
F = spm_figure('GetWin','Graphics');
V = spm_get(1,'spmT*.img','Select the T image');
blobs(fno).name = V;
set(findobj(F,'Tag',['Filename',num2str(fno)]),'String',blobs(fno).name);
case 'show_coords'
%% this should mainly be called by clicking on the
figure
%% and it shows the current coordinates in a text box
centre = varargin{2};
centre = centre(:);
if(~all(size(centre)==[3,1]));
disp('cannot show coords')
end
F = spm_figure('getwin');
hCoord = findobj(F,'Tag','Label_coords');
set(hCoord,'String',num2str(round(centre')));
end
---Christian Keysers <[log in to unmask]> wrote ---
> Dear SPM
>
> I often have to do further computations on random effect
analysis, such
> as for instance combining the results of two different
tasks (e.g. motor
> and visual tasks) to look for overlaps. I do this by doing
two separate
> rfx analysis, which create me their own spmT maps. Then
using image
> calc, I do i1.*(i2>3.11), which gives me the T values of
the first RFX
> but only in voxels where the other task was significant as
well. Not to
> report that in a paper, its nice to have a table. In SPM99
I handled
> that by creating a bogus contrast in the first RFX, look
at it, and then
> rename my image calc image as spmT003.img/hdr. Then
visualizing my bogus
> contrast, I saw my image calc SPMT, and could visualise a
nice results
> table. In SPM2 that doesn't work anymore, because of the
header info's
> etc which are different. Any idea how to get a table from
an image calc
> T map?
>
> Christian
>
> --
> Christian Keysers, PhD
> Associate Professor
>
> BCN Neuro-Imaging Center
> University Medical Center Groningen
> Antonius Deusinglaan 2 (room 125)
> 9713 AW Groningen
>
> Phone: +31 50 3638794
> Fax: +31 50 3638875
>
--- End of quote ---
|
