[FieldTrip] please help! sourcemodel and grid don't match
Julian Keil
julian.keil at gmail.com
Wed Mar 30 15:55:53 CEST 2016
Hi Hassan,
have you double-checked that all your input units are in the same in all elements (i.e. electrodes, volume, mri, grid all have the same .unit and .pos values in the same ballpark)?
You can use ft_convert_units to fix possible differences.
Good luck,
Julian
Am 30.03.2016 um 03:19 schrieb Hassan Aleem:
> Hi all,
>
> When I plot the sourcemodel and the headmodel, the sourcemodel is many magnitudes larger than the headmodel. Not sure if it has anything to do with my code so attaching that below. Also attaching the jpgs.
>
> %% read the parameters
>
> %% Parameter 2.2
>
> mrifile='Siemens_AF8-MPRAGE_20141022_1.nii';
> mritype='MNI';
> realign='realigned';
> volumereslicing='volumeresliced';
> brainthresh=0.5;
> scalpthresh=0.10;
> scalpsmooth='scalp_smooth';
> scalp_smooth.parameter=9;
> scalpvertices=1000;
> skullvertices=1000;
> brainvertices=1000;
> meshmethod='projectmesh';
> headmodel_method='dipoli';
> subjnum='915';
> time='pre';
> timetype='hassan';
> timewindow=[0 .1];
> gridtype='template10';
> beamformer='lcmv';
> sourcemodeltype= 'templatebased';
> identifier='Parameter2-2';
> ===================================
> %%read in mri
>
> switch mritype
> case 'MNI'
> [mri] = ft_read_mri(mrifile);
> mri.coordsys ='MNI';
> case 'Analyze'
> [mri] = ft_read_mri(mrifile);
> mri.coordsys ='Analyze';
> case 'SPM'
> [mri] = ft_read_mri(mrifile);
> mri.coordsys = 'SPM';
> end
> save(strcat(mrifile,'_',identifier,'_','mri.mat'),'mri')
> cfg=[];
> ft_sourceplot(cfg,mri)
> ==============================
>
> %% realigning
> switch realign
> case 'realigned'
> cfg=[];
> cfg.method='fiducial';
> cfg.fiducial.nas = [6.60688 6.30230 -2.94229];
> cfg.fiducial.lpa = [6.04082 7.65872 0.350950];
> cfg.fiducial.rpa = [4.41106 8.71481 3.50199];
> mri_realigned=ft_volumerealign(cfg,mri);
> ft_sourceplot(cfg,mri_realigned)
> save(strcat(mrifile,'_',identifier,'_','mri_realigned.mat'),'mri_realigned')
> case 'notrealigned'
> %do nothing
> end
>
> =========================
> %% volume reslicing
> switch volumereslicing
> case 'volumeresliced'
> cfg = [];
> cfg.resolution = 1;
> mri_resliced = ft_volumereslice(cfg, mri);
> save(strcat(mrifile,'_',identifier,'_','mri_resliced.mat'),'mri_resliced')
> case 'notresliced'
> %do nothing
> end
> ft_sourceplot(cfg,mri_resliced)
> ==============================
>
> %% segmentation parameters
> switch scalpsmooth
> case 'noscalpsmooth'
> cfg = [];
> cfg.brainthreshold = brainthresh;
> cfg.scalpthreshold = scalpthresh;
> cfg.downsample = 1; %no downsampling
> cfg.output = {'brain' 'skull' 'scalp'};
> seg = ft_volumesegment(cfg, mri);
> case 'scalp_smooth'
> cfg = [];
> cfg.brainthreshold = brainthresh;
> cfg.scalpthreshold = scalpthresh;
> cfg.scalpsmooth = scalp_smooth.parameter;
> cfg.downsample = 1; %no downsampling
> cfg.output = {'brain' 'skull' 'scalp'};
> seg = ft_volumesegment(cfg, mri);
> end
>
> seg.anatomy = mri.anatomy;
> cfg = [];
> cfg.funparameter = 'trishells';
> cfg.funcolormap = [1 0 0;0 1 0;0 0 1];
> cfg.maskparameter=cfg.funparameter;
> ft_sourceplot(cfg,seg)
>
> save(strcat(mrifile,scalpsmooth,'_',identifier,'_','seg.mat'),'seg')
> =================================
> %% prepare mesh
>
> cfg = [];
> cfg.method= meshmethod;
> cfg.tissue = {'scalp', 'skull', 'brain'};
> cfg.numvertices = [scalpvertices skullvertices brainvertices];
> bnd = ft_prepare_mesh(cfg, seg);
> % i enlarged / reduced the boundary mesh a little bit, because they were causing errors, probably due to overlapping
> bnd(1).pnt = bnd(1).pnt.*1.001;
> bnd(2).pnt = bnd(2).pnt.*1;
> bnd(3).pnt = bnd(3).pnt.*0.999;
>
> figure;
> ft_plot_mesh(bnd(1))%scalp
> figure;
> ft_plot_mesh(bnd(2))%skull
> figure;
> ft_plot_mesh(bnd(3))%brain
>
> %% save mesh
> save(strcat(mrifile,cfg.method,num2str(scalpvertices),num2str(skullvertices),num2str(brainvertices),'_',identifier,'_','bnd.mat'),'bnd')
> =============================================
>
> %% build headmodel
>
> switch headmodel_method
> case 'dipoli'
> cfg = [];
> cfg.method = 'dipoli';
> case 'bemcp'
> cfg = [];
> cfg.method = 'bemcp';
> end
>
> headmodel = ft_prepare_headmodel(cfg, bnd);
> headmodel = ft_convert_units(headmodel,'cm');
>
> save(strcat(mrifile,cfg.method,'_',identifier,'_','headmodel.mat'),'headmodel');
>
> ====================================
> %% adding electrodes to headmodel
> elec=ft_read_sens('Hydrocel_GSN_128_1.0_TRIM_mod.sfp'); %load electrode file
> % visualize head surface (scalp)
> figure;
> ft_plot_mesh(headmodel.bnd(1), 'edgecolor','none','facealpha',0.8,'facecolor',[0.6 0.6 0.8]);
> hold on;
> % plot electrodes
> ft_plot_sens(elec,'style', 'sk')
> =================================
> %most likely they are not aligned, so for aligning electrodes
> cfg = [];
> cfg.method = 'interactive';%can try the 'fiducial' method but is very
> %likely not to fix it
> cfg.elec = elec;
> % nas=mri_realigned.cfg.fiducial.nas;
> % lpa=mri_realigned.cfg.fiducial.lpa;
> % rpa=mri_realigned.cfg.fiducial.rpa;
> cfg.headshape = headmodel.bnd(1);%scalp
> elec_aligned = ft_electroderealign(cfg);%use the gui to manually align the electrodes
>
> save(strcat(mrifile,'_',identifier,'_','elec_aligned.mat'),'elec_aligned')
> ==================================
> %% load EEG data
> load(strcat('subj',subjnum,'_',time,'_FT_bothStim_rejChan500hz0.1hz-30hz_fromContinuous.mat'), 'dataM6w')
> load(strcat('subj',subjnum,'_',time,'_FT_bothStim_rejChan500hz0.1hz-30hz_fromContinuous.mat'), 'dataM0')
> load(strcat('subj',subjnum,'_',time,'_FT_bothStim_rejChan500hz0.1hz-30hz_fromContinuous.mat'), 'dataM3W')
> load(strcat('subj',subjnum,'_',time,'_FT_bothStim_rejChan500hz0.1hz-30hz_fromContinuous.mat'), 'dataM3Bw')
> load(strcat('subj',subjnum,'_',time,'_FT_bothStim_rejChan500hz0.1hz-30hz_fromContinuous.mat'), 'dataM3Bb')
> load(strcat('subj',subjnum,'_',time,'_FT_bothStim_rejChan500hz0.1hz-30hz_fromContinuous.mat'), 'dataM6b')
> load(strcat('subj',subjnum,'_',time,'_FT_bothStim_rejChan500hz0.1hz-30hz_fromContinuous.mat'), 'dataM6')
> ====================================
> %% timelock data
>
> cfg = [];
> cfg.refef= 'yes'; %added at 2-5
>
>
> same = dataM0;
> same.trial = [dataM0.trial dataM3W.trial];
> same.time = [dataM0.time dataM3W.time];
>
> diff = dataM6;
> diff.trial = [dataM3Bb.trial dataM6.trial];
> diff.time = [dataM3Bb.time dataM6.time];
>
> switch timetype
> case 'simple'
> timelock_m0 = ft_timelockanalysis(cfg, dataM0);
> timelock_m6 = ft_timelockanalysis(cfg, dataM6);
> timelock_m3w = ft_timelockanalysis(cfg, dataM3W);
> timelock_m3b = ft_timelockanalysis(cfg, dataM3Bb);
> timelock_diff = ft_timelockanalysis(cfg, diff);
> timelock_same = ft_timelockanalysis(cfg, same);
> case 'hassan'
> cfg.vartrllength = 0; %was changed to 0 from 2
> cfg.covariance = 'yes';
> cfg.covariancewindow = timewindow;
> cfg.keeptrials = 'yes';
> timelock_m0 = ft_timelockanalysis(cfg, dataM0);
> timelock_m6 = ft_timelockanalysis(cfg, dataM6);
> timelock_m3w = ft_timelockanalysis(cfg, dataM3W);
> timelock_m3b = ft_timelockanalysis(cfg, dataM3Bb);
> timelock_diff = ft_timelockanalysis(cfg, diff);
> timelock_same = ft_timelockanalysis(cfg, same);
>
> end
>
> save(strcat(mrifile,timetype,time,'_',identifier,'_','timelock_datam0.mat'),'timelock_m0')
> save(strcat(mrifile,timetype,time,'_',identifier,'_','timelock_datam6.mat'),'timelock_m6')
> save(strcat(mrifile,timetype,time,'_',identifier,'_','timelock_datam3w.mat'),'timelock_m3w')
> save(strcat(mrifile,timetype,time,'_',identifier,'_','timelock_data3b.mat'),'timelock_m3b')
> save(strcat(mrifile,timetype,time,'_',identifier,'_','timelock_diff.mat'),'timelock_diff')
> save(strcat(mrifile,timetype,time,'_',identifier,'_','timelock_same.mat'),'timelock_same')
>
> ============================
> %% make sourcemodel
>
> switch sourcemodeltype
> case 'templatebased'
> template=load('standard_sourcemodel3d10mm.mat');
> elec=ft_read_sens('Hydrocel_GSN_128_1.0_TRIM_mod.sfp'); %load electrode file
> cfg = [];
> cfg.headmodel = headmodel;
> cfg.mri=mri;
> cfg.grid.warpmni = 'yes';
> cfg.grid.template = template.sourcemodel;
> cfg.grid.nonlinear = 'yes';
> cfg.resolution = 2; % resolution of the template grid in mm
> cfg.moveinward = 1; % actually uses vol mesh
> cfg.grid.unit ='mm';
> cfg.inwardshift = 0; % needs to be expressed to work with moveinward
> cfg.elec = elec;
> % cfg.headshape=('cortex_5124.surf.gii');
> mysourcemodel = ft_prepare_sourcemodel(cfg);
>
> case 'individual'
> load('standard_sourcemodel3d5mm.mat');
> cfg = [];
> cfg.grid.warpmni = 'yes';
> cfg.resolution = 5; % resolution of the template grid in mm
> cfg.grid.nonlinear = 'yes'; % use non linear normalization
> cfg.mri = mri; % use subject's mri in mni coordinates
> cfg.grid.unit = 'mm';
> mysourcemodel = ft_prepare_sourcemodel(cfg);
> case 'template'
> mysourcemodel=load('standard_sourcemodel3d10mm.mat');%doensn't work
> end
>
>
>
> save(strcat(mrifile,'resolution',num2str(cfg.resolution),'_',subjnum,'_',identifier,'_','sourcemodel.mat'),'mysourcemodel')
> =======================
> %% Plot the single subject'e head model and grid positions
>
> figure
> hold on
> ft_plot_vol(headmodel,'edgecolor','none','facealpha',0.3);
> ft_plot_mesh(mysourcemodel.pos(mysourcemodel.inside,:));
> ========================
> %% Double check that source model is aligned to the segmented volume and to the electrodes
> HERE IS WHERE I SEE THE ERROR
> figure
> hold on
> plot3(mysourcemodel.pos(mysourcemodel.inside,1),mysourcemodel.pos(mysourcemodel.inside,2),mysourcemodel.pos(mysourcemodel.inside,3),'.');
> ft_plot_vol(headmodel,'facecolor','brain','edgecolor','none','facealpha',0.3);
> ft_plot_sens(elec_aligned,'label','label','style','y*','Markersize',10);
> ========================================
> %% make leadfield
>
>
> switch gridtype
> case 'individual'
> cfg = [];
> cfg.elec = elec;
> cfg.headmodel = headmodel;
> cfg.normalize = 'yes';
> cfg.reducerank= 3;
> cfg.grid.unit = 'cm';% same unit as above, i.e. in cm
> [grid] = ft_prepare_leadfield(cfg);
> case 'templatebased'
> load('standard_sourcemodel3d10mm.mat')
> cfg = [];
> cfg.elec = elec;
> cfg.template = sourcemodel.cfg.grid;
> cfg.grid=mysourcemodel.cfg.grid;
> cfg.headmodel = headmodel;
> cfg.reducerank = 3; % default is 3 for EEG, 2 for MEG
> cfg.channel = 'all';
> cfg.normalize= 'yes';
> cfg.grid.unit = 'cm';
> cfg.grid.tight = 'yes';
> [grid] = ft_prepare_leadfield(cfg);
> case 'template10'
> load('standard_sourcemodel3d10mm.mat')
> cfg.normalize='yes';
> [grid]=sourcemodel.cfg.grid;
> end
>
>
> save(strcat(mrifile,gridtype,'_',identifier,'_','leadfield.mat'),'grid')
> ===========================
> %% Interpolate
> cfg = [];
> cfg.interpmethod = 'linear';
> cfg.parameter = 'pow';
> sourceinterp_diff = ft_sourceinterpolate(cfg, sourceanalysis_diff, mri);
> sourceinterp_same = ft_sourceinterpolate(cfg, sourceanalysis_same, mri);
> sourceinterp_m0 = ft_sourceinterpolate(cfg, sourceanalysis_m0, mri);
> sourceinterp_m6 = ft_sourceinterpolate(cfg, sourceanalysis_m6, mri);
> sourceinterp_m3w = ft_sourceinterpolate(cfg, sourceanalysis_m3w, mri);
> sourceinterp_m3b = ft_sourceinterpolate(cfg, sourceanalysis_m3b, mri);
>
>
> % save(strcat(mrifile,'_',identifier,'_','sourceinterp_diff.mat'),'sourceinterp_diff')
> % save(strcat(mrifile,'_',identifier,'_','sourceinterp_m0.mat'),'sourceinterp_m0')
> % save(strcat(mrifile,'_',identifier,'_','sourceinterp_m6.mat'),'sourceinterp_m6')
> % save(strcat(mrifile,'_',identifier,'_','sourceinterp_m3w.mat'),'sourceinterp_m3w')
> % save(strcat(mrifile,'_',identifier,'_','sourceinterp_m3b.mat'),'sourceinterp_m3b')
> % save(strcat(mrifile,'_',identifier,'_','sourceinterp_same.mat'),'sourceinterp_same')
> =================================
>
> %%beamforming
>
> cfg = [];
> cfg.method = 'slice';
> % cfg.slicedim= 1;
> cfg.funcolormap = 'jet';
> % cfg.method = 'ortho';
> % cfg.funcolorlim = 'zeromax';
> % cfg.funcolorlim = [0.01 0.2];
> cfg.funparameter= 'pow';
>
>
> test2=sourceinterp_diff;
> test2.pow=((sourceinterp_diff.pow-sourceinterp_same.pow)./(sourceinterp_diff.pow+sourceinterp_same.pow));
> test=sourceinterp_diff;
> test.pow=sourceinterp_diff.pow-sourceinterp_same.pow;
> ft_sourceplot(cfg,test)
> ft_sourceplot(cfg,test2)
> ft_sourceplot(cfg,sourceinterp_same);
> ft_sourceplot(cfg,sourceinterp_diff);
> ft_sourceplot(cfg,sourceinterp_m0);
> ft_sourceplot(cfg,sourceinterp_m6);
> ft_sourceplot(cfg,sourceinterp_m3w);
> ft_sourceplot(cfg,sourceinterp_m3b);
> <figforfilegrid.jpg><figforemail.jpg><combinedforemail.jpg>_______________________________________________
> fieldtrip mailing list
> fieldtrip at donders.ru.nl
> http://mailman.science.ru.nl/mailman/listinfo/fieldtrip
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.science.ru.nl/pipermail/fieldtrip/attachments/20160330/b905c277/attachment-0002.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: signature.asc
Type: application/pgp-signature
Size: 495 bytes
Desc: Message signed with OpenPGP using GPGMail
URL: <http://mailman.science.ru.nl/pipermail/fieldtrip/attachments/20160330/b905c277/attachment-0002.sig>
More information about the fieldtrip
mailing list