% Script generated by Brainstorm (20-Jul-2018)
% Auteur : Edouard Delaire, Melanie Pelegrini-Issac

%Definition of protocol 
%protocol_name='FNIRS_nback';
protocol_name='Protocol_NIRS';

nirs_data_folder= 'C:\Users\Debashis\Documents\MATLAB\n-back_NatureDataset\NIRS';
% eeg_data_folder='/NAS/home/edelaire/Documents/data/shin/eeg';

conditions = {'dsr','nback','wg'};
n_subject = 2;

%% No modification needed after this line

if ~brainstorm('status')
    brainstorm
end
bst_report('Start');

% Verfication of repository
if ~exist(nirs_data_folder) 
    bst_report('Error', [] , [] , ['Data folder ' nirs_data_folder ' not found']);
    ReportFile = bst_report('Save', []);
    bst_report('Open',ReportFile);
    return;     
end    

if isempty(bst_get('Protocol', protocol_name))
    gui_brainstorm('CreateProtocol', protocol_name, 1, 0); % UseDefaultAnat=1, UseDefaultChannel=0
end

% Set template for default anatomy
nst_bst_set_template_anatomy('Colin27_4NIRS_Jan19');


for i_subject = 1:n_subject
    subject_name =  sprintf('VP%03d',i_subject);
    subject_foler= fullfile( nirs_data_folder,subject_name);
    if ~exist(subject_foler) 
        bst_report('Error', [] , [] , ['Data folder ' subject_foler ' not found']);
        ReportFile = bst_report('Save', []);
        bst_report('Open',ReportFile);
    return;    
    end
    
   [sSubject, iSubject] = db_add_subject(subject_name, [], 1, 0);
   
   for i_cond = 1:length(conditions)
       cnt=load(fullfile(subject_foler, sprintf('cnt_%s.mat',conditions{i_cond})));
       mnt=load(fullfile(subject_foler, sprintf('mnt_%s.mat',conditions{i_cond})));
       mrk=load(fullfile(subject_foler, sprintf('mrk_%s.mat',conditions{i_cond})));

       % note: access to the field using cnt.(sprintf('cnt_%s',conditions{i_cond}))
       
        
        % Step1: Create channel file structure
        ChannelMat = db_template('channelmat');
        ChannelMat.Comment = 'NIRS-BRS channels';
        ChannelMat.Nirs.Wavelengths = [];
        montage = mnt.(sprintf('mnt_%s',conditions{i_cond}));

        nChannels = length(montage.clab);
        % NIRS channels
        for iChan = 1:nChannels
            % This assume measure are raw; need to change for Hbo,HbR,HbT  
            chromophores={'HbO','HbR','HbT'};
            for iC = 1:length(chromophores)
                ChannelMat.Channel(3*(iChan-1) + iC).Name = sprintf('S%dD%d%s',iChan,iChan,chromophores{iC});
                ChannelMat.Channel(3*(iChan-1) + iC).Group = chromophores{iC};
                ChannelMat.Channel(3*(iChan-1) + iC).Type = 'NIRS';
                ChannelMat.Channel(3*(iChan-1) + iC).Comment = montage.clab{iChan};
                ChannelMat.Channel(3*(iChan-1) + iC).Weight = 1;
                ChannelMat.Channel(3*(iChan-1) + iC).Loc(:,1) = montage.pos_3d(:, iChan);
                ChannelMat.Channel(3*(iChan-1) + iC).Loc(:,2) = montage.pos_3d(:, iChan);
                ChannelMat.Channel(3*(iChan-1) + iC).Orient  = [];
            end
        end
        
        iStudy = db_add_condition(subject_name, [conditions{i_cond} '_NIRS']);
        sStudy = bst_get('Study', iStudy);
        
        [tmp, iChannelStudy] = bst_get('ChannelForStudy', iStudy);
        db_set_channel(iChannelStudy, ChannelMat, 0, 0);
        
        % Step2: Load data file structure
        data = cnt.(sprintf('cnt_%s',conditions{i_cond}));
        F=zeros(3*nChannels, size(data.oxy.x,1));
        % Create matrix data
        for iChan = 1:nChannels
            F(3*(iChan-1) + 1,:) = data.oxy.x(:,iChan)';
            F(3*(iChan-1) + 2,:) = data.deoxy.x(:,iChan)';
            F(3*(iChan-1) + 3,:) = data.oxy.x(:,iChan)' + data.deoxy.x(:,iChan)';
        end    
        time= 0:1/data.oxy.fs:(size(data.oxy.x,1)-1)/data.oxy.fs;
        
        % step3: Load Events
        events = mrk.(sprintf('mrk_%s',conditions{i_cond}));
        bst_events = repmat(db_template('event'), 1, length(events.className));
        
        for iEvt = 1:length(events.className)

            bst_events(iEvt).label      = events.className{iEvt} ;
            evtTime = events.time( events.y(iEvt,:) == 1)/1000; % convert from ms to s
            bst_events(iEvt).times      = evtTime;
            bst_events(iEvt).epochs     = ones(1, size(evtTime,2));
            bst_events(iEvt).channels   = cell(1, size(evtTime,2));
            bst_events(iEvt).notes      = cell(1, size(evtTime,2));
            bst_events(iEvt).color      = rand(1,3);
            bst_events(iEvt).reactTimes = [];
        end   
        
        sDataOut = db_template('data');
        sDataOut.F            = F;
        sDataOut.Comment      = 'NIRS Hb';
        %sDataOut.ChannelFlag  = sDataIn.ChannelFlag;
        sDataOut.ChannelFlag  = ones(size(F, 1), 1);
        sDataOut.Time         =  time;
        sDataOut.DataType     = 'recordings'; 
        sDataOut.nAvg         = 1;
        sDataOut.Events       = bst_events;
        sDataOut.History      = [];
        sDataOut              = bst_history('add', sDataOut, 'process', 'Importing data');
        sDataOut.DisplayUnits = strrep(data.oxy.yUnit,'/L','.l-1');
        
        % Generate a new file name in the same folder
        OutputFile = bst_process('GetNewFilename', bst_fileparts(sStudy.FileName), 'data_HB');
        sDataOut.FileName = file_short(OutputFile);
        bst_save(OutputFile, sDataOut, 'v7');
        % Register in database
        db_add_data(iStudy, OutputFile, sDataOut); 
   end
   
%     %% Part II. Import EEG
%     
%     subject_foler= fullfile( eeg_data_folder,subject_name);
%     for i_cond = 1:length(conditions)
%        cnt=load(fullfile(subject_foler, sprintf('cnt_%s.mat',conditions{i_cond})));
%        mnt=load(fullfile(subject_foler, sprintf('mnt_%s.mat',conditions{i_cond})));
%        mrk=load(fullfile(subject_foler, sprintf('mrk_%s.mat',conditions{i_cond})));
% 
%        % note: access to the field using cnt.(sprintf('cnt_%s',conditions{i_cond}))
%        
%         
%         % Step1: Create channel file structure
%         ChannelMat = db_template('channelmat');
%         ChannelMat.Comment = 'BrainAmp channels';
%         montage = mnt.(sprintf('mnt_%s',conditions{i_cond}));
% 
%         nChannels = length(montage.clab);
%         % eeg channels
%         for iChan = 1:nChannels
%             ChannelMat.Channel(iChan).Name = montage.clab{iChan};
%             ChannelMat.Channel(iChan).Group = [];
%             ChannelMat.Channel(iChan).Type = 'EEG';
%             ChannelMat.Channel(iChan).Comment = montage.clab{iChan};
%             ChannelMat.Channel(iChan).Weight = 1;
%             ChannelMat.Channel(iChan).Loc = montage.pos_3d(:, iChan);
%             ChannelMat.Channel(iChan).Orient  = [];
%         end
%         
%         iStudy = db_add_condition(subject_name, [conditions{i_cond} '_EEG']);
%         sStudy = bst_get('Study', iStudy);
%         
%         [tmp, iChannelStudy] = bst_get('ChannelForStudy', iStudy);
%         db_set_channel(iChannelStudy, ChannelMat, 0, 0);
%         
%         % Step2: Load data file structure
%         data = cnt.(sprintf('cnt_%s',conditions{i_cond}));
%         F=zeros(nChannels, size(data.x,1));
%         % Create matrix data
%         for iChan = 1:nChannels
%             F(iChan,:) = data.x(:,iChan)';
%         end    
%         time= 0:1/data.fs:(size(data.x,1)-1)/data.fs;
%         
%         % step3: Load Events
%         events = mrk.(sprintf('mrk_%s',conditions{i_cond}));
%         bst_events = repmat(db_template('event'), 1, length(events.className));
%         
%         for iEvt = 1:length(events.className)
% 
%             bst_events(iEvt).label      = events.className{iEvt} ;
%             evtTime = events.time( events.y(iEvt,:) == 1)/1000; % convert from ms to s
%             bst_events(iEvt).times      = evtTime;
%             bst_events(iEvt).epochs     = ones(1, size(evtTime,2));
%             bst_events(iEvt).channels   = cell(1, size(evtTime,2));
%             bst_events(iEvt).notes      = cell(1, size(evtTime,2));
%             bst_events(iEvt).color      = rand(1,3);
%             bst_events(iEvt).reactTimes = [];
%         end   
%         
%         sDataOut = db_template('data');
%         sDataOut.F            = F;
%         sDataOut.Comment      = 'EEG_raw';
%         %sDataOut.ChannelFlag  = sDataIn.ChannelFlag;
%         sDataOut.ChannelFlag  = ones(size(F, 1), 1);
%         sDataOut.Time         =  time;
%         sDataOut.DataType     = 'recordings'; 
%         sDataOut.nAvg         = 1;
%         sDataOut.Events       = bst_events;
%         sDataOut.History      = [];
%         sDataOut              = bst_history('add', sDataOut, 'process', 'Importing data');
%         sDataOut.DisplayUnits = data.yUnit;
%         
%         % Generate a new file name in the same folder
%         OutputFile = bst_process('GetNewFilename', bst_fileparts(sStudy.FileName), 'data_HB');
        sDataOut.FileName = file_short(OutputFile);
        bst_save(OutputFile, sDataOut, 'v7');
        % Register in database
%         db_add_data(iStudy, OutputFile, sDataOut); 
%    end
end  
db_reload_database('current');


