Description of vb_combine_meg

0001 function [bexp, bexp_ext, refmg] = ...
0002   vb_combine_meg_files(in_files, out_file, ext_spec)
0003 % Combine .meg.mat files (for VBMEG version 0.8)
0004 %
0005 % --- Syntax
0006 % [bexp, bexp_ext, refmg] = vb_combine_meg_files(in_files, out_file, ext_spec)
0007 %
0008 % --- Input
0009 % in_files (cell array): .meg.mat files to be combined
0010 % out_file (String)    : Combined .meg.mat file is saved to this file.
0011 % ext_spec (struct) : extra specification
0012 %          .saveman : <optional> <<struct>>
0013 %                   :   .data_dir : directory (relative path) for external file
0014 %                   :   .precision : usually 'float64'
0015 %
0016 % --- Output
0017 %   concatenated data as folows:
0018 %   bexp     [Nchannel x Nsample x Ntrial]
0019 %   bexp_ext [Nchannel x Nsample x Ntrial]
0020 %   refmeg   [Nchannel x Nsample x Ntrial]
0021 %
0022 %   !!!NOTICE!!! If data is stored externally, these parameters will be empty
0023 %
0024 % --- Note
0025 % [preconditions]
0026 %  * following data must correspond with each other
0027 %  * Measurement
0028 %  * CoordType
0029 %  * MEGinfo
0030 %     .Nchannel
0031 %     .Nsample
0032 %     .device
0033 %     .SampleFreq
0034 %     .Pretrigger
0035 %    (.sensor_weight)
0036 %    (.MEGch_id)
0037 %    (.MEGch_name)
0038 %
0039 % [postconditions]
0040 %
0041 %  * data format after combining is as follows:
0042 %    (Sensor position will be averaged across .meg.mat files.
0043 %     Accordingly, it is assumed that MEG data in combined files
0044 %     were recorded without head movement.)
0045 %    pick
0046 %    Qpick
0047 %    ref_pick
0048 %    ref_Qpick
0049 %
0050 %    (apply the values of the first input file)
0051 %    CoordType
0052 %    Measurement
0053 %
0054 %    (be concatenated in "trial" dimension)
0055 %    bexp
0056 %    bexp_ext
0057 %    refmg
0058 %
0059 %    (the fields of MEGinfo if it exist)
0060 %    MEGinfo (standard format)
0061 %      .Measurement   (apply the value of the first file)
0062 %      .device        (apply the value of the first file)
0063 %      .Nchannel      (apply the value of the first file)
0064 %      .Nsample       (apply the value of the first file)
0065 %      .Nrepeat       (summation of all input data)
0066 %      .SampleFreq    (apply the value of the first file)
0067 %      .Pretrigger    (apply the value of the first file)
0068 %      .sensor_weight (apply the value of the first file)
0069 %      .MEGch_id      (apply the value of the first file)
0070 %      .MEGch_name    (apply the value of the first file)
0071 %      .ActiveChannel (apply the value of the first file)
0072 %      .ActiveTrial   (apply the value of the first file)
0073 %      .Vcenter       (average)
0074 %      .Vradius       (average)
0075 %      .MEG_ID        (set empty)
0076 %      .MRI_ID        (apply the value of the first file)
0077 %      .Trial
0078 %      .ChannelInfo
0079 %      .ExtraChannelInfo
0080 %      .saveman
0081 %
0082 % --- History
0083 % 2007-01-12 (Taku Yoshioka) initial version for VBMEG 0.5
0084 % 2009-09-02 (Sako) modified according to the latest format and added return
0085 % 2010-08-31 (Sako) modified combined file could be used as a input file of
0086 %                   vb_load_meg_data file with loadspec.
0087 % 2011-02-24 (Sako) when making new file, sampling data is stored outer
0088 % 2011-06-01 (Sako) modified according to the new standard format
0089 %
0090 % Copyright (C) 2011, ATR All Rights Reserved.
0091 % License : New BSD License(see VBMEG_LICENSE.txt)
0092 
0093 % --- CHECK ARGUMENTS --- %
0094 if ~exist('in_files', 'var'), in_files = []; end
0095 if ~exist('out_file', 'var'), out_file = []; end
0096 if ~exist('ext_spec', 'var'), ext_spec = []; end
0097 [in_files, out_file, info, measurement, coord_type, saveman, d_path] = ...
0098   inner_check_arguments(in_files, out_file, ext_spec);
0099 
0100 % [pick,Qpick] = vb_load_sensor(in_files{1});
0101 
0102 % --- initialize base variables
0103 pick_sum      = [];
0104 Qpick_sum     = [];
0105 ref_pick_sum  = [];
0106 ref_Qpick_sum = [];
0107 
0108 vcenter_sum   = [];
0109 vradius_sum   = 0;
0110 
0111 n_pick      = 0;
0112 n_Qpick     = 0;
0113 n_ref_pick  = 0;
0114 n_ref_Qpick = 0;
0115 
0116 n_vcenter = 0;
0117 n_vradius = 0;
0118 
0119 n_trial   = 0;
0120 
0121 n_file = length(in_files);
0122 
0123 % --- when you want to store data outside
0124 if ~isempty(saveman)
0125   % ----- common information to store externally
0126   const = vb_define_yokogawa;
0127   data_ext = const.EXT_CHANNEL_BIN;
0128 
0129   prec = vb_saveman_get_precision(saveman);
0130 end
0131 
0132 bexp = [];
0133 bexp_ext = [];
0134 refmg = [];
0135 
0136 % --- combine file by file
0137 for i_file = 1:n_file
0138   load(in_files{i_file}, 'pick', 'Qpick', 'ref_pick', 'ref_Qpick');
0139   
0140   % --- sensor coordinates
0141   if exist('pick', 'var')
0142     if isempty(pick_sum)
0143       pick_sum = zeros(size(pick,1), size(pick,2));
0144     end
0145     pick_sum = pick_sum + pick;
0146     n_pick = n_pick + 1;
0147   end
0148 
0149   if exist('Qpick', 'var')
0150     if isempty(Qpick_sum)
0151       Qpick_sum = zeros(size(Qpick,1), size(Qpick,2));
0152     end
0153     Qpick_sum = Qpick_sum + Qpick;
0154     n_Qpick = n_Qpick + 1;
0155   end
0156 
0157   if exist('ref_pick', 'var')
0158     if isempty(ref_pick_sum)
0159       ref_pick_sum = zeros(size(ref_pick,1), size(ref_pick,2));
0160     end
0161     ref_pick_sum = ref_pick_sum + ref_pick;
0162     n_ref_pick = n_ref_pick + 1;
0163   end
0164 
0165   if exist('ref_Qpick', 'var')
0166     if isempty(ref_Qpick_sum)
0167       ref_Qpick_sum = zeros(size(ref_Qpick,1), size(ref_Qpick,2));
0168     end
0169     ref_Qpick_sum = ref_Qpick_sum + ref_Qpick;
0170     n_ref_Qpick = n_ref_Qpick + 1;
0171   end
0172   
0173   if isfield(info{i_file}, 'Vcenter')
0174     if isempty(vcenter_sum)
0175       vcenter_sum = ...
0176         zeros(size(info{i_file}.Vcenter,1), size(info{i_file}.Vcenter,2));
0177     end
0178     vcenter_sum = vcenter_sum + info{i_file}.Vcenter;
0179     n_vcenter = n_vcenter + 1;
0180   end
0181   
0182   if isfield(info{i_file}, 'Vradius')
0183     vradius_sum = vradius_sum + info{i_file}.Vradius;
0184     n_vradius = n_vradius + 1;
0185   end
0186   
0187   % --- data
0188   loadspec.ActiveChannel = 0;
0189   loadspec.ChannelType = 'MEG';
0190   cur_bexp = vb_load_meg_data(in_files{i_file}, loadspec);
0191   
0192   loadspec.ChannelType = 'EXTRA';
0193   cur_bexp_ext = vb_load_meg_data(in_files{i_file}, loadspec);
0194 
0195   loadspec.ChannelType = 'REFERENCE';
0196   cur_refmg = vb_load_meg_data(in_files{i_file}, loadspec);
0197   
0198   % --- get trial info
0199   if ~isempty(cur_bexp)
0200     trial_begin = n_trial + 1;
0201     trial_end   = trial_begin + size(cur_bexp, 3) - 1;
0202     n_trial = trial_end;
0203     
0204   elseif ~isempty(cur_bexp_ext)
0205     trial_begin = n_trial + 1;
0206     trial_end   = trial_begin + size(cur_bexp_ext, 3) - 1;
0207     n_trial = trial_end;
0208 
0209   elseif ~isempty(cur_refmg)
0210     trial_begin = n_trial + 1;
0211     trial_end   = trial_begin + size(cur_refmg, 3) - 1;
0212     n_trial = trial_end;
0213   end
0214 
0215   if ~isempty(saveman)
0216     meginfo = info{1};
0217   end
0218   
0219   if ~isempty(cur_bexp)
0220     if ~isempty(saveman)
0221       meg_labels = vb_meginfo_get_channel_label_meg(meginfo);
0222 
0223       n_ch = length(meg_labels);
0224       for i_ch = 1:n_ch
0225         vb_util_make_external_data_file(d_path, meg_labels{i_ch}, ...
0226           data_ext, prec, cur_bexp(i_ch,:,:), 'ab' );
0227       end
0228 
0229     else
0230       bexp(:,:, (trial_begin:trial_end)) = cur_bexp;
0231     end
0232   end
0233 
0234   if ~isempty(cur_bexp_ext)
0235     if ~isempty(saveman)
0236       ext_labels = vb_meginfo_get_channel_label_extra(meginfo);
0237 
0238       n_ch = length(ext_labels);
0239       for i_ch = 1:n_ch
0240         vb_util_make_external_data_file(d_path, ext_labels{i_ch}, ...
0241           data_ext, prec, cur_bexp_ext(i_ch,:,:), 'ab' );
0242       end
0243 
0244     else
0245       bexp_ext(:,:, (trial_begin:trial_end)) = cur_bexp_ext;
0246     end
0247   end
0248 
0249   if ~isempty(cur_refmg)
0250     if ~isempty(saveman)
0251       ref_labels = vb_meginfo_get_channel_label_refmg(meginfo);
0252       
0253       n_ch = length(ref_labels);
0254       for i_ch = 1:n_ch
0255         vb_util_make_external_data_file(d_path, ref_labels{i_ch}, ...
0256           data_ext, prec, cur_refmg(i_ch,:,:), 'ab' );
0257       end
0258 
0259     else
0260       refmg(:,:, (trial_begin:trial_end)) = cur_refmg;
0261     end
0262   end
0263 end
0264 
0265 % --- sensor cordinates are by averaging
0266 if n_pick > 0
0267   pick = pick_sum / n_pick;
0268 end
0269 
0270 if n_Qpick > 0
0271   Qpick = Qpick_sum / n_Qpick;
0272 end
0273 
0274 if n_ref_pick > 0
0275   ref_pick = ref_pick_sum / n_ref_pick;
0276 end
0277 
0278 if n_ref_Qpick > 0
0279   ref_Qpick = ref_Qpick_sum / n_ref_Qpick;
0280 end
0281 
0282 % --- new MEGinfo
0283 if isfield(info{1}, 'Measurement')
0284   MEGinfo.Measurement = info{1}.Measurement;
0285 else
0286   MEGinfo.Measurement = measurement{1};
0287 end
0288 MEGinfo.device     = info{1}.device;
0289 MEGinfo.Nchannel   = info{1}.Nchannel;
0290 MEGinfo.Nsample    = info{1}.Nsample;
0291 MEGinfo.Nrepeat    = n_trial;
0292 MEGinfo.SampleFreq = info{1}.SampleFreq;
0293 MEGinfo.Pretrigger = info{1}.Pretrigger;
0294 MEGinfo.MEG_ID     = '';
0295 MEGinfo.MRI_ID     = info{1}.MRI_ID;
0296 
0297 if ~isempty(saveman)
0298   MEGinfo.saveman = saveman;
0299 end
0300 
0301 if isfield(info{1}, 'sensor_weight')
0302   MEGinfo.sensor_weight = info{1}.sensor_weight;
0303 end
0304 
0305 if isfield(info{1}, 'MEGch_id')
0306   MEGinfo.MEGch_id = info{1}.MEGch_id;
0307 end
0308 
0309 if isfield(info{1}, 'MEGch_name')
0310   MEGinfo.MEGch_name = info{1}.MEGch_name;
0311 end
0312 
0313 if n_vcenter > 0
0314   MEGinfo.Vcenter = vcenter_sum / n_vcenter;
0315 end
0316 
0317 if n_vradius > 0
0318   MEGinfo.Vradius = vradius_sum / n_vradius;
0319 end
0320 
0321 % --- for not minimum version
0322 MEGinfo.Trial = [];
0323 MEGinfo.ActiveTrial = [];
0324 prev_n_trial = 0;
0325 prev_n_sample = 0;
0326 
0327 for i_file = 1:n_file
0328   cur_n_trial = length(info{i_file}.Trial);
0329   cur_n_sample = info{i_file}.Nsample;
0330   
0331   for i_tr = 1:cur_n_trial
0332     info{i_file}.Trial(i_tr).sample = ...
0333       info{i_file}.Trial(i_tr).sample + prev_n_sample;
0334     info{i_file}.Trial(i_tr).number = ...
0335       info{i_file}.Trial(i_tr).number + prev_n_trial;
0336   end
0337   MEGinfo.Trial = [MEGinfo.Trial; info{i_file}.Trial];
0338   MEGinfo.ActiveTrial = [MEGinfo.ActiveTrial; info{i_file}.ActiveTrial];
0339   prev_n_trial = prev_n_trial + cur_n_trial;
0340   prev_n_sample = prev_n_sample + cur_n_sample;
0341 end
0342 
0343 if isfield(info{1}, 'ChannelInfo')
0344   MEGinfo.ChannelInfo = info{1}.ChannelInfo;
0345 end
0346 
0347 if isfield(info{1}, 'ActiveChannel')
0348   MEGinfo.ActiveChannel = info{1}.ActiveChannel;
0349 end
0350 
0351 if isfield(info{1}, 'ExtraChannelInfo')
0352   MEGinfo.ExtraChannelInfo = info{1}.ExtraChannelInfo;
0353 end
0354 
0355 % --- top header
0356 Measurement = measurement{1};
0357 CoordType   = coord_type{1};
0358 
0359 % ----- PositionFile
0360 target_key = 'PositionFile';
0361 PositionFile = '';
0362 
0363 for  i_file = 1:n_file
0364   [state, guide_def] = ...
0365     vb_util_check_variable_in_matfile(in_files{i_file}, target_key);
0366   if state == guide_def.VALID
0367     load(in_files{i_file}, target_key);
0368     break;
0369   end
0370 end
0371 
0372 if ~isempty(out_file)
0373   if ~isempty(saveman)
0374     if (n_ref_pick > 0) && (n_ref_Qpick > 0)
0375       vb_fsave(out_file, 'Measurement', ...
0376         'pick','Qpick', 'ref_pick', 'ref_Qpick', 'CoordType', 'MEGinfo', ...
0377         'PositionFile');
0378     else
0379       vb_fsave(out_file, 'Measurement', ...
0380         'pick','Qpick', 'CoordType', 'MEGinfo', ...
0381         'PositionFile');
0382     end
0383   else
0384     if (n_ref_pick > 0) && (n_ref_Qpick > 0)
0385       vb_fsave(out_file, 'Measurement', 'bexp', 'bexp_ext', 'refmg', ...
0386         'pick','Qpick', 'ref_pick', 'ref_Qpick', 'CoordType', 'MEGinfo', ...
0387         'PositionFile');
0388     else
0389       vb_fsave(out_file, 'Measurement', 'bexp', 'bexp_ext', 'refmg', ...
0390         'pick','Qpick', 'CoordType', 'MEGinfo', ...
0391         'PositionFile');
0392     end
0393   end
0394 end
0395 return;
0396 %
0397 % --- END OF MAIN PROCEDURE -------------------------------------------------- %
0398 
0399 % --- INNER FUNCTIONS -------------------------------------------------------- %
0400 %
0401 % --- inner_check_arguments()
0402 %
0403 function [in_files, out_file, info, measurement, coord_type, ...
0404   saveman, real_dir] = ...
0405   inner_check_arguments(in_files, out_file, ext_spec)
0406 func_ = mfilename;
0407 
0408 if isempty(in_files)
0409   error('(%s) in_files is a required parameter', func_);
0410 end
0411 
0412 if isempty(out_file)
0413   % require no actions
0414 end
0415 
0416 if ~iscell(in_files)
0417   in_files = {in_files};
0418 end
0419 
0420 % --- check in_files
0421 n_files = length(in_files);
0422 
0423 for i_file = 1:n_files
0424   if exist(in_files{i_file}, 'file') ~= 2
0425     error('(%s) cannot find in_file{%d} (%s)', ...
0426       func_, i_file, in_files{i_file});
0427   end
0428   
0429   % --- collect information to check preconditions
0430   info{i_file} = vb_load_measurement_info(in_files{i_file});
0431 
0432   load(in_files{i_file}, 'Measurement', 'CoordType');
0433   measurement{i_file} = Measurement;
0434   coord_type{i_file} = CoordType;
0435 end
0436 
0437 % --- preconditions to check
0438 %  * following data must correspond with each other
0439 %  * Measurement
0440 %  * CoordType
0441 %  * MEGinfo
0442 %     .Nchannel
0443 %     .Nsample
0444 %     .device
0445 %     .SampleFreq
0446 %     .Pretrigger
0447 %    (.sensor_weight)
0448 %    (.MEGch_id)
0449 %    (.MEGch_name)
0450 
0451 % minimum set of MEGinfo
0452 n_channel_base   = info{1}.Nchannel;
0453 n_sample_base    = info{1}.Nsample;
0454 sample_freq_base = info{1}.SampleFreq;
0455 pretrigger_base  = info{1}.Pretrigger;
0456 device_base      = info{1}.device;
0457 
0458 % extra info of MEGinfo
0459 if isfield(info{1}, 'sensor_weight')
0460   sensor_weight_base = info{1}.sensor_weight;
0461 else
0462   sensor_weight_base = [];
0463 end
0464 
0465 if isfield(info{1}, 'MEGch_id')
0466   MEGch_id_base = info{1}.MEGch_id;
0467 else
0468   MEGch_id_base = [];
0469 end
0470 
0471 if isfield(info{1}, 'MEGch_name')
0472   MEGch_name_base = info{1}.MEGch_name;
0473 else
0474   MEGch_name_base = {};
0475 end
0476 
0477 if isfield(info{1}, 'ActiveChannel')
0478   ActiveChannel_base = info{1}.ActiveChannel;
0479 else
0480   ActiveChannel_base = {};
0481 end
0482 
0483 measurement_base = measurement{1};
0484 coord_type_base  = coord_type{1};
0485 
0486 for i_file = 2:n_files
0487 
0488   % --- basic fields
0489   % ----- Nchannel
0490   if info{i_file}.Nchannel ~= n_channel_base
0491     error('(%s - %s) all numbers of channel must be the same', ...
0492       func_, in_files{i_file});
0493   end
0494 
0495   % ----- Nsample
0496   if info{i_file}.Nsample ~= n_sample_base
0497     error('(%s - %s) all numbers of sample must be the same', ...
0498       func_, in_files{i_file});
0499   end
0500 
0501   % ----- SampleFreq
0502   if info{i_file}.SampleFreq ~= sample_freq_base
0503     error('(%s - %s) all sampling frequencies must be the same', ...
0504       func_, in_files{i_file});
0505   end
0506 
0507   % ----- Pretrigger
0508   if info{i_file}.Pretrigger ~= pretrigger_base
0509     error('(%s - %s) all pretriggers must be the same', ...
0510       func_, in_files{i_file});
0511   end
0512 
0513   % ----- device
0514   if ~strcmp(info{i_file}.device, device_base)
0515     error('(%s - %s) each device must be the same', ...
0516       func_, in_files{i_file});
0517   end
0518 
0519   % ----- Measurement
0520   if ~strcmp(measurement{i_file}, measurement_base)
0521     error('(%s - %s) each Measurement must be the same', ...
0522       func_, in_files{i_file});
0523   end
0524   
0525   % --- extra fields
0526   % ----- sensor_weight
0527   if isfield(info{i_file}, 'sensor_weight')
0528     if ~isequal(info{i_file}.sensor_weight, sensor_weight_base)
0529       error('(%s - %s) each sensor_weight must be the same', ...
0530         func_, in_files{i_file});
0531     end
0532   end
0533 
0534   % ----- MEGch_id
0535   if isfield(info{i_file}, 'MEGch_id')
0536     if ~isequal(info{i_file}.MEGch_id, MEGch_id_base)
0537       error('(%s - %s) each MEGch_id must be the same', ...
0538         func_, in_files{i_file});
0539     end
0540   end
0541 
0542   % ----- MEGch_name
0543   if isfield(info{i_file}, 'MEGch_name')
0544     if ~isequal(info{i_file}.MEGch_name, MEGch_name_base)
0545       error('(%s - %s) each MEGch_name must be the same', ...
0546         func_, in_files{i_file});
0547     end
0548   end
0549 
0550   % ----- CoordType
0551   if ~strcmp(coord_type{i_file}, coord_type_base)
0552     error('(%s - %s) each CoordType must be the same', ...
0553       func_, in_files{i_file});
0554   end
0555   
0556   % ----- ActiveChannel
0557   if isfield(info{i_file}, 'ActiveChannel')
0558     if ~isequal(info{i_file}.ActiveChannel, ActiveChannel_base)
0559       error('(%s - %s) each ActiveChannel must be the same', ...
0560         func_, in_files{i_file});
0561     end
0562   end
0563 
0564 end
0565 
0566 saveman = [];
0567 real_dir = [];
0568 
0569 if ~isempty(out_file)
0570   if isempty(ext_spec) ...
0571       || ~isfield(ext_spec, 'saveman') ...
0572       || ~vb_saveman_get_switch(ext_spec.saveman)
0573     [ext_spec.saveman.data_dir, real_dir] = vb_device_make_data_dir(out_file);
0574     ext_spec.saveman.precision = 'float64';
0575   else
0576     f_path = vb_get_file_parts(out_file);
0577     real_dir = [f_path '/' ext_spec.saveman.data_dir];
0578   end
0579   saveman = ext_spec.saveman;
0580 end
0581 return;
0582 %
0583 % --- end of inner_check_arguments()
0584 %
0585 % --- END OF INNER FUNCTIONS ------------------------------------------------- %
0586 
0587 % --- END OF FILE --- %
vb_combine_meg_files

PURPOSE

SYNOPSIS

DESCRIPTION

CROSS-REFERENCE INFORMATION

SUBFUNCTIONS

SOURCE CODE
