surf2vol_interp

(Toolbox/surf2vol_interp.m in BrainStorm 2.0 (Alpha))

Function Synopsis

Cvol = surf2vol_interp(subjectimage,Faces,Vertices,Csurf);

Help Text

SURF2VOL_INTERP - (overwrite succinct one line summary here)

Cross-Reference Information

get_user_directory                 C:\BrainStorm_2001\Toolbox\get_user_directory.m
tri_interp                         C:\BrainStorm_2001\Toolbox\tri_interp.m

Listing of function C:\BrainStorm_2001\Toolbox\surf2vol_interp.m

%SURF2VOL_INTERP - (overwrite succinct one line summary here)

%<autobegin> ---------------------- 26-May-2004 11:34:35 -----------------------
% --------- Automatically Generated Comments Block Using AUTO_COMMENTS ---------
%
% CATEGORY: Archive
%
% At Check-in: $Author: Mosher $  $Revision: 11 $  $Date: 5/26/04 10:02a $
%
% This software is part of BrainStorm Toolbox Version 2.0 (Alpha) 24-May-2004
% 
% Principal Investigators and Developers:
% ** Richard M. Leahy, PhD, Signal & Image Processing Institute,
%    University of Southern California, Los Angeles, CA
% ** John C. Mosher, PhD, Biophysics Group,
%    Los Alamos National Laboratory, Los Alamos, NM
% ** Sylvain Baillet, PhD, Cognitive Neuroscience & Brain Imaging Laboratory,
%    CNRS, Hopital de la Salpetriere, Paris, France
% 
% See BrainStorm website at http://neuroimage.usc.edu for further information.
% 
% Copyright (c) 2004 BrainStorm by the University of Southern California
% This software distributed  under the terms of the GNU General Public License
% as published by the Free Software Foundation. Further details on the GPL
% license can be found at http://www.gnu.org/copyleft/gpl.html .
% 
% FOR RESEARCH PURPOSES ONLY. THE SOFTWARE IS PROVIDED "AS IS," AND THE
% UNIVERSITY OF SOUTHERN CALIFORNIA AND ITS COLLABORATORS DO NOT MAKE ANY
% WARRANTY, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO WARRANTIES OF
% MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE, NOR DO THEY ASSUME ANY
% LIABILITY OR RESPONSIBILITY FOR THE USE OF THIS SOFTWARE.
%<autoend> ------------------------ 26-May-2004 11:34:35 -----------------------


subjectimage = 'I:\Data\Subjects\epilepsie\01Jan23\fres_subjectimage.mat';
load('Rfres_white-CTF_0.11475Lfres_white-CTF_0.11733_tess.mat')
Faces = Faces{1};
Vertices = Vertices{1};
load('I:\Data\Studies\Epilepsie\01Jan23\run01-f.ds\BrainStorm\run01-f_data_NativeFormat_MNE_MEGresults_1757.mat','ImageGridAmp')
Csurf = ImageGridAmp(:,13);
clear ImageGridAmp;


function Cvol = surf2vol_interp(subjectimage,Faces,Vertices,Csurf);
% SURF2VOL_INTERP : Voxel base interpolation of cortical current density values
% 
% Voxel base interpolation of cortical current density values by computing the multi-linear interpolated 
% current values in the MR voxel volume from the original cortical tessellation 
% 
% subjectimage : file name of the BrainStorm's MRI file - PCS needs to be defined first in the MRI tool
% Faces        : Faces array of cortical envelope
% Vertices     : Vertice location array of cortical envelope
% Csurf        : Current denisity values as defined on the original cortical envelope
%
% Cvol         : Current density values interpolated into the MR volume

% /---Script Author--------------------------------------\
% |                                                      |
% | *** Sylvain Baillet Ph.D.                            |
% | Cognitive Neuroscience & Brain Imaging Laboratory    |
% | CNRS UPR640 - LENA                                   | 
% | Hopital de la Salpetriere, Paris, France             |
% | sylvain.baillet@chups.jussieu.fr                     |
% |                                                      |
% \------------------------------------------------------/
%  
% Date of creation: February 2003
% 
% Script History ------------------------------------------------------------
% ---------------------------------------------------------------------------


% As those obtained from MNE or non-linear imaging methods
MAINWINDOW = gcbf;
set(0,'CurrentFigure',MAINWINDOW )
mri_scan = get(gcbf,'Userdata');

image_size = findobj(gcbf,'Tag','image_size');
SLICE_thickness = findobj(gcbf,'Tag','slice_thickness');
FOV = findobj(gcbf,'Tag','fov');
NSLICE = findobj(gcbf,'Tag','nslice');
data_format = findobj(gcbf,'Tag','data_format');  
file_mri= findobj(gcbf,'Tag','file_mri');  
scan_type= findobj(gcbf,'Tag','scan_type');  



mri_scan.mript = Vertices';
clear Vertices;
load(subjectimage)

% See the whole tessellation
%GridAmp = 10*ones(size(mri_scan.mript)); 
%ImageGridAmp = GridAmp ;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

mri_scan.GridAmp = abs(Csurf); 
tmax = 1;% See integrated activity at each vertex along time

mri_scan.M = max(abs(mri_scan.GridAmp));

mri_scan.climm = max(abs(Cube(:)));

mri_scan.mri(mri_scan.mri >= mri_scan.climm) = mri_scan.climm;
mri_scan.GridAmp = mri_scan.GridAmp*mri_scan.climm/mri_scan.M;

mri_scan.thres = .2*mri_scan.climm; % Threshold for visualization of current density

mri_scan.mript = 1000*mri_scan.mript; % CHEAT

%--------

% Transform PCS coordinates in the MRI volume reference system
if ~exist('PCS','var')
    errordlg('Please Define a PCS First')
    return
end

% center = mri_scan.center;
% rot = mri_scan.rot;
% mri_scan.mript = (inv(rot)*mri_scan.mript')';      
% 
% mri_scan.mript (:,1) = mri_scan.mript (:,1) + center(1);
% mri_scan.mript (:,2) = mri_scan.mript (:,2) + center(2);
% mri_scan.mript (:,3) = mri_scan.mript (:,3) + center(3);
% icoord(:,1) = round((get(s1,'Max')-mri_scan.mript (:,2))/mri_scan.pix_siz(2)+1);
% icoord(:,2) = round(mri_scan.mript (:,1)/mri_scan.pix_siz(1)+1);
% icoord(:,3) = round((get(s2,'Max')-mri_scan.mript (:,3))/mri_scan.pix_siz(3)+1);

icoord = ClearPCS2MRI(mri_scan.mript,R,T);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%ind = sub2ind(size(mri_scan.mri),icoord(:,1),icoord(:,2),icoord(:,3));
%[uind,ii] = unique(ind);

%%%% Interpolation of the Cortical Current Density Image
if ~exist('interptoto_mat.mat','file')
    [H] = tri_interp(size(Cube),icoord, Faces)   ;
    disp('Saving Interpolation Matrix')
    User = get_user_directory;
    save interptoto_mat H
else
    load interp_mat
end

CInterp = sparse(H*mri_scan.GridAmp(:,tmax));
mri_scan.csources = {[1:size(H,2)]};
clear H mri_scan.GridAmp

ind =find(CInterp>=mri_scan.thres);

mri_scan.indmodif= ind; % indices of the voxels with modified grey levels
mri_scan.gl_orig = mri_scan.mri(ind);

%  mri_scan.mri(ind) = 1.01*mri_scan.climm + nonzeros(CInterp(ind));
mri_scan.mri(ind) = double(mri_scan.mri(ind)) + nonzeros(CInterp(ind));
minn = double(min(mri_scan.mri(ind)));
mri_scan.mri(ind) =  double(mri_scan.mri(ind)) - minn + mri_scan.climm;

%%%----- Superimpose to regular MRI
%  mri_scan.mri = mri_scan.mri + iiiblob;

[M tmp] = max(mri_scan.mri(ind));
[is ic ia] = ind2sub(size(mri_scan.mri),ind(tmp));

ax = mri_scan.ax;
sag = mri_scan.sag;
cor = mri_scan.cor;

set(sag ,'CData',(flipud(squeeze(mri_scan.mri(is,:,:))')),'erasemode','none');
set(get(sag,'Parent'),'CLim',[0 2*mri_scan.climm])
set(cor ,'CData',(fliplr(flipud(squeeze(mri_scan.mri(:,ic,:))'))),'erasemode','none');
set(get(cor,'Parent'),'CLim',[0 2*mri_scan.climm])
set(ax ,'CData',(mri_scan.mri(:,:,ia)),'erasemode','none');    
set(get(ax,'Parent'),'CLim',[0  2*mri_scan.climm])

split(1:mri_scan.climm,:) = bone(mri_scan.climm);
split(mri_scan.climm+1:2*mri_scan.climm,:) = hot(mri_scan.climm);
mri_scan.split = split;

set(MAINWINDOW,'Userdata',mri_scan);


ax = findobj(MAINWINDOW,'Tag','colorbar');
axes(ax)
colormapp = split( mri_scan.climm+1:end,:);
n = size(colormapp,1);
image([0, 1], [0 10],(n+1:2*n)');
set(ax,'Ydir','normal','Tag','colorbar')
set(ax,'xtick',[],'ycolor',frontcolor,'xcolor',frontcolor)
colormap(split)