test_vb_calc_current

 Test vb_calc_current.m, a sample program of VBMEG.

 [syntax]
 Z = test_vb_calc_current(proj_root,test_parm)

 [input]
 proj_root: <<string>> VBMEG project root directory. 
 test_parm: <<struct>> Test condition variables. 
 --- fields of test_parm
  cortex : <<string>> 'bv' for BrainVoyager, 'fs' for FreeSurfer. 
  data   : <<string>> 'meg' for MEG data, 'eeg' for EEG data. 
  method : <<string>> 'sarvas' for Sarvas, 'bem' for BEM. 
  prior  : <<string>> 'sparse' for sparse (noninformative), 'fmri' for
           fMRI prior. 
  tempwin: <<string>> 'single' for single time window, 'multiple' for
           multiple time windows. 
 ---

 [output]
 Z: <<matrix>> Z-current. 
 Temporal sum (51-150 ms after stimulus onset) of estimated current is
 superimposed on the inflated cortical surface. 

 [history]
 2010-06-30 Taku Yoshioka
 2010-09-29 Taku Yoshioka
  Compatibility check wich 'bexp' usage of vb_calc_current.m added. 

 Copyright (C) 2011, ATR All Rights Reserved.
 License : New BSD License(see VBMEG_LICENSE.txt)

0001 function Z = test_vb_calc_current(proj_root,test_parm)
0002 % Test vb_calc_current.m, a sample program of VBMEG.
0003 %
0004 % [syntax]
0005 % Z = test_vb_calc_current(proj_root,test_parm)
0006 %
0007 % [input]
0008 % proj_root: <<string>> VBMEG project root directory.
0009 % test_parm: <<struct>> Test condition variables.
0010 % --- fields of test_parm
0011 %  cortex : <<string>> 'bv' for BrainVoyager, 'fs' for FreeSurfer.
0012 %  data   : <<string>> 'meg' for MEG data, 'eeg' for EEG data.
0013 %  method : <<string>> 'sarvas' for Sarvas, 'bem' for BEM.
0014 %  prior  : <<string>> 'sparse' for sparse (noninformative), 'fmri' for
0015 %           fMRI prior.
0016 %  tempwin: <<string>> 'single' for single time window, 'multiple' for
0017 %           multiple time windows.
0018 % ---
0019 %
0020 % [output]
0021 % Z: <<matrix>> Z-current.
0022 % Temporal sum (51-150 ms after stimulus onset) of estimated current is
0023 % superimposed on the inflated cortical surface.
0024 %
0025 % [history]
0026 % 2010-06-30 Taku Yoshioka
0027 % 2010-09-29 Taku Yoshioka
0028 %  Compatibility check wich 'bexp' usage of vb_calc_current.m added.
0029 %
0030 % Copyright (C) 2011, ATR All Rights Reserved.
0031 % License : New BSD License(see VBMEG_LICENSE.txt)
0032 
0033 vb_disp('--- Start: test_vb_calc_current');
0034 resultdir = './vbmeg_result/';
0035 
0036 %
0037 % Filename of inverse filter
0038 %
0039 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0040 filterfile ...
0041     = [resultdir 'sbj_test_UR_' test_parm.cortex '_' test_parm.data '_' ...
0042        test_parm.method '_' test_parm.prior '_' test_parm.tempwin ...
0043        '.vbfilt.mat'];
0044 
0045 %
0046 % Current calculation parameters
0047 %
0048 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0049 curr_parm.overlap_mode = 0;
0050 curr_parm.ix_area = [];
0051 curr_parm.trial_average = false;
0052 curr_parm.ix_trial = 1;
0053 curr_parm.tsubsmpl = [];
0054 
0055 %
0056 % Load inverse filter variables
0057 %
0058 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0059 load(filterfile,'VBfilt','bayes_parm')
0060 
0061 %
0062 % Execute vb_calc_current
0063 %
0064 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0065 [Z,W,Jinfo] = vb_calc_current(proj_root,VBfilt,bayes_parm,curr_parm);
0066 ix_act = Jinfo.ix_act;
0067 ix_act_ex = Jinfo.ix_act_ex;
0068 Ztmp = Z(:,551:650,:);
0069 
0070 %
0071 % Check compatibility for usage with MEG data array 'bexp'
0072 %
0073 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0074 VBfilt2 = VBfilt;
0075 bexp{1} = vb_load_meg_data(bayes_parm.megfile{1});
0076 vb_struct2vars(bayes_parm,{'twin_meg','Tperiod','Tnext'});
0077 bexp{1} = bexp{1}(:,twin_meg(1):twin_meg(2),:);
0078 VBfilt2.Ntrial(1) = size(bexp{1},3);
0079 VBfilt2.Tsample = twin_meg(2)-twin_meg(1)+1;
0080 VBfilt2.Twindow = vb_calc_timewindow(VBfilt2.Tsample,Tperiod,Tnext);
0081 Z2 = vb_calc_current(proj_root,VBfilt2,bexp,curr_parm);
0082 d = sqrt(sum(sum((Z-Z2).^2))/sum(sum(Z.^2)));
0083 vb_disp_nonl(['Difference with ''bexp'' usage: ' sprintf('%1.3f\n',d)]);
0084 
0085 %
0086 % Calculate and plot temporal sum of cortical current
0087 %
0088 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0089 brainfile ...
0090     = [proj_root resultdir 'sbj_test_' test_parm.cortex '.brain.mat'];
0091 [V,F,xx,inf_C] = vb_load_cortex(brainfile,'Inflate');
0092 
0093 J = zeros(size(V,1),1);
0094 J(ix_act_ex) = sqrt(sum((W*Ztmp).*(W*Ztmp),2));
0095 vb_disp(['Cortical model vertices: ' num2str(size(V,1))]);
0096 vb_disp(['J-current vertices     : ' num2str(length(ix_act_ex))]);
0097 vb_disp(['Z-current vertices     : ' num2str(length(ix_act))]);
0098 vb_disp(['Sample time points     : ' num2str(size(Z,2))]);
0099 
0100 %
0101 % Plot timecourse of Z-current
0102 %
0103 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0104 figure; 
0105 subplot(3,2,[1 2]);
0106 plot(Z');
0107 
0108 %
0109 % Plot cortical current map on inflated model
0110 %
0111 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0112 subplot(3,2,3);
0113 plot_parm = vb_set_plot_parm;
0114 plot_parm.LRflag = 'L';
0115 vb_plot_cortex(plot_parm,V,F,inf_C,J);
0116 view([-90 0]);
0117 axis equal; 
0118 axis off;
0119 axis tight;
0120 
0121 subplot(3,2,4);
0122 plot_parm = vb_set_plot_parm;
0123 plot_parm.LRflag = 'R';
0124 vb_plot_cortex(plot_parm,V,F,inf_C,J);
0125 view([90 0]);
0126 axis equal; 
0127 axis off;
0128 axis tight;
0129 
0130 subplot(3,2,5);
0131 plot_parm = vb_set_plot_parm;
0132 plot_parm.LRflag = 'L';
0133 vb_plot_cortex(plot_parm,V,F,inf_C,J);
0134 view([90 0]);
0135 axis equal; 
0136 axis off;
0137 axis tight;
0138 
0139 subplot(3,2,6);
0140 plot_parm = vb_set_plot_parm;
0141 plot_parm.LRflag = 'R';
0142 vb_plot_cortex(plot_parm,V,F,inf_C,J);
0143 view([-90 0]);
0144 axis equal; 
0145 axis off;
0146 axis tight;
0147 
0148 vb_disp('--- End  : test_vb_calc_current');
0149 
0150 return;