function [hdr] = ft_read_header(filename, varargin)

% FT_READ_HEADER reads header information from a variety of EEG, MEG and other time

% series data files and represents the header information in a common

% data-independent structure. The supported formats are listed below.

%

% Use as

% hdr = ft_read_header(filename, ...)

%

% Additional options should be specified in key-value pairs and can be

% 'headerformat' = string

% 'fallback' = can be empty or 'biosig' (default = [])

% 'checkmaxfilter' = boolean, whether to check that maxfilter has been correctly applied (default = true)

% 'chanindx' = list with channel indices in case of different sampling frequencies (only for EDF)

% 'chantype' = string or cell-array with strings, channel types to be read (only for NeuroOmega and BlackRock)

% 'coordsys' = string, 'head' or 'dewar' (default = 'head')

% 'headerformat' = name of a MATLAB function that takes the filename as input (default is automatic)

% 'password' = password structure for encrypted data set (for dhn_med10, mayo_mef30, mayo_mef21)

% 'readbids' = string, 'yes', no', or 'ifmakessense', whether to read information from the BIDS sidecar files (default = 'ifmakessense')

% 'splitlabel' = string, 'yes' or 'no', split the channel labels on the '-' and return the first part (default = 'yes' for CTF and FieldLine, 'no' for others)

%

% This returns a header structure with the following fields

% hdr.Fs = sampling frequency

% hdr.nChans = number of channels

% hdr.nSamples = number of samples per trial

% hdr.nSamplesPre = number of pre-trigger samples in each trial

% hdr.nTrials = number of trials

% hdr.label = Nx1 cell-array with the label of each channel

% hdr.chantype = Nx1 cell-array with the channel type, see FT_CHANTYPE

% hdr.chanunit = Nx1 cell-array with the physical units, see FT_CHANUNIT

%

% For continuously recorded data, this will return nSamplesPre=0 and nTrials=1.

%

% For some data formats that are recorded on animal electrophysiology

% systems (e.g. Neuralynx, Plexon), the following optional fields are

% returned, which allows for relating the timing of spike and LFP data

% hdr.FirstTimeStamp number, represented as 32-bit or 64-bit unsigned integer

% hdr.TimeStampPerSample number, represented in double precision

%

% Depending on the file format, additional header information can be

% returned in the hdr.orig subfield.

%

% To use an external reading function, you can specify an external function as the

% 'headerformat' option. This function should take the filename as input argument.

% Please check the code of this function for details, and search for BIDS_TSV as

% example.

%

% The following MEG dataformats are supported

% CTF (*.ds, *.res4, *.meg4)

% Neuromag/Elekta/Megin (*.fif)

% BTi/4D (*.m4d, *.pdf, *.xyz)

% Yokogawa/Ricoh (*.ave, *.con, *.raw)

% NetMEG (*.nc)

% ITAB - Chieti (*.mhd)

% Tristan Babysquid (*.fif)

% York Instruments (*.meghdf5)

%

% The following EEG dataformats are supported

% ANT - Advanced Neuro Technology, EEProbe (*.avr, *.eeg, *.cnt)

% BCI2000 (*.dat)

% Biosemi (*.bdf)

% BrainVision (*.eeg, *.seg, *.dat, *.vhdr, *.vmrk)

% CED - Cambridge Electronic Design (*.smr)

% EGI - Electrical Geodesics, Inc. (*.egis, *.ave, *.gave, *.ses, *.raw, *.sbin, *.mff)

% GTec (*.mat, *.hdf5)

% Generic data formats (*.edf, *.gdf)

% Megis/BESA (*.avr, *.swf, *.besa)

% NeuroScan (*.eeg, *.cnt, *.avg)

% Nexstim (*.nxe)

% TMSi (*.Poly5)

% Mega Neurone (directory)

% Natus/Nicolet/Nervus (.e files)

% Nihon Kohden (*.m00, *.EEG)

% Bitalino OpenSignals (*.txt)

% OpenBCI (*.txt)

%

% The following spike and LFP dataformats are supported

% Neuralynx (*.ncs, *.nse, *.nts, *.nev, *.nrd, *.dma, *.log)

% Plextor (*.nex, *.plx, *.ddt)

% CED - Cambridge Electronic Design (*.smr)

% MPI - Max Planck Institute (*.dap)

% Neurosim (neurosim_spikes, neurosim_signals, neurosim_ds)

% Windaq (*.wdq)

% NeuroOmega (*.mat transformed from *.mpx)

% Neurodata Without Borders: Neurophysiology (*.nwb)

%

% The following NIRS dataformats are supported

% Artinis - Artinis Medical Systems B.V. (*.oxy3, *.oxy4, *.oxy5, *.oxyproj)

% BUCN - Birkbeck college, London (*.txt)

% SNIRF - Society for functional near-infrared spectroscopy (*.snirf)

%

% The following Eyetracker dataformats are supported

% EyeLink - SR Research (*.asc)

% SensoMotoric Instruments - (*.txt)

% Tobii - (*.tsv)

%

% See also FT_READ_DATA, FT_READ_EVENT, FT_WRITE_DATA, FT_WRITE_EVENT,

% FT_CHANTYPE, FT_CHANUNIT

% Copyright (C) 2003-2024 Robert Oostenveld

%

% This file is part of FieldTrip, see http://www.fieldtriptoolbox.org

% for the documentation and details.

%

% FieldTrip is free software: you can redistribute it and/or modify

% it under the terms of the GNU General Public License as published by

% the Free Software Foundation, either version 3 of the License, or

% (at your option) any later version.

%

% FieldTrip is distributed in the hope that it will be useful,

% but WITHOUT ANY WARRANTY; without even the implied warranty of

% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the

% GNU General Public License for more details.

%

% You should have received a copy of the GNU General Public License

% along with FieldTrip. If not, see <http://www.gnu.org/licenses/>.

%

% $Id$

% TODO channel renaming should be made a general option (see bham_bdf)

persistent cacheheader % for caching the full header

persistent cachechunk % for caching the res4 chunk when doing realtime analysis on the CTF scanner

persistent db_blob % for fcdc_mysql

global FIFF % for fiff-files

if isempty(db_blob)

db_blob = false;

end

if iscell(filename)

% use recursion to read the header from multiple files

ft_warning('concatenating header from %d files', numel(filename));

hdr = cell(size(filename));

for i=1:numel(filename)

hdr{i} = ft_read_header(filename{i}, varargin{:});

end

allhdr = cat(1, hdr{:});

if numel(unique([allhdr.label]))==sum([allhdr.nChans])

% each file has different channels, concatenate along the channel dimension

for i=1:numel(filename)

assert(isequal(hdr{i}.Fs, hdr{1}.Fs), 'sampling rates are not consistent over files');

assert(isequal(hdr{i}.nSamples, hdr{1}.nSamples), 'number of samples is not consistent over files');

assert(isequal(hdr{i}.nTrials, hdr{1}.nTrials), 'number of trials is not consistent over files');

end

combined = hdr{1}; % copy the first header as the general one

combined.label = [allhdr.label];

combined.chanunit = [allhdr.chanunit];

combined.chantype = [allhdr.chantype];

combined.nChans = sum([allhdr.nChans]);

combined.orig = hdr; % store the original header details of all files

else

% each file has the same channels, concatenate along the time dimension

ntrl = nan(size(filename));

nsmp = nan(size(filename));

for i=1:numel(filename)

assert(isequal(hdr{i}.Fs, hdr{1}.Fs), 'sampling rates are not consistent over files');

assert(isequal(hdr{i}.label, hdr{1}.label), 'channels are not consistent over files');

ntrl(i) = hdr{i}.nTrials;

nsmp(i) = hdr{i}.nSamples;

end

% the subsequent code concatenates the files over time, i.e. each file has the same channels

combined = hdr{1}; % copy the first header as the general one

combined.orig = hdr; % store the original header details of all files

if all(ntrl==1)

% each file is a continuous recording

combined.nTrials = ntrl(1);

combined.nSamples = sum(nsmp);

elseif all(nsmp==nsmp(1))

% each file holds segments of the same length

combined.nTrials = sum(ntrl);

combined.nSamples = nsmp(1);

else

ft_error('cannot concatenate files');

end

end % concatenate over channels or over time

% return the header of the concatenated datafiles

hdr = combined;

return

end

% get the options

headerformat = ft_getopt(varargin, 'headerformat');

retry = ft_getopt(varargin, 'retry', false); % the default is not to retry reading the header

chanindx = ft_getopt(varargin, 'chanindx'); % this is used for EDF with different sampling rates

coordsys = ft_getopt(varargin, 'coordsys', 'head'); % this is used for ctf and neuromag_mne, it can be head or dewar

coilaccuracy = ft_getopt(varargin, 'coilaccuracy'); % empty, or a number between 0-2

coildeffile = ft_getopt(varargin, 'coildeffile'); % empty, or a filename

chantype = ft_getopt(varargin, 'chantype', {});

password = ft_getopt(varargin, 'password', struct([]));

readbids = ft_getopt(varargin, 'readbids', 'ifmakessense');

splitlabel = ft_getopt(varargin, 'splitlabel'); % default for CTF and FieldLine is dealt with below

% this should be a cell array

if ~iscell(chantype); chantype = {chantype}; end

% optionally get the data from the URL and make a temporary local copy

filename = fetch_url(filename);

if isempty(headerformat)

% only do the autodetection if the format was not specified

headerformat = ft_filetype(filename);

end

if iscell(headerformat)

% this happens for datasets specified as cell-array for concatenation

headerformat = headerformat{1};

end

if strcmp(headerformat, 'compressed')

% the file is compressed, unzip on the fly

inflated = true;

filename = inflate_file(filename);

headerformat = ft_filetype(filename);

else

inflated = false;

end

% for backward compatibility with https://github.com/fieldtrip/fieldtrip/issues/1585

if islogical(readbids)

% it should be either yes/no/ifmakessense

if readbids

readbids = 'yes';

else

readbids = 'no';

end

end

realtime = any(strcmp(headerformat, {'fcdc_buffer', 'ctf_shm', 'fcdc_mysql'}));

% The checkUniqueLabels flag is used for the realtime buffer in case

% it contains fMRI data. It prevents 1000000 voxel names to be checked

% for uniqueness. fMRI users will probably never use channel names

% for anything.

if realtime

% skip the rest of the initial checks to increase the speed for realtime operation

checkUniqueLabels = false;

% the cache and fallback option should always be false for realtime processing

cache = false;

fallback = false;

else

% check whether the file or directory exists

if ~exist(filename, 'file')

ft_error('file or directory ''%s'' does not exist', filename);

end

checkUniqueLabels = true;

% get the rest of the options, this is skipped for realtime operation

cache = ft_getopt(varargin, 'cache');

fallback = ft_getopt(varargin, 'fallback');

checkmaxfilter = ft_getopt(varargin, 'checkmaxfilter', true);

if isempty(cache)

if any(strcmp(headerformat, {'bci2000_dat', 'eyelink_asc', 'gtec_mat', 'gtec_hdf5', 'mega_neurone', 'nihonkohden_m00', 'smi_txt', 'biosig'}))

cache = true;

else

cache = false;

end

end

% ensure that the headerfile and datafile are defined, which are sometimes different than the name of the dataset

[filename, headerfile, datafile] = dataset2files(filename, headerformat);

if ~strcmp(filename, headerfile) && ~ft_filetype(filename, 'ctf_ds') && ~ft_filetype(filename, 'fcdc_buffer_offline') && ~ft_filetype(filename, 'fcdc_matbin')

filename = headerfile; % this function should read the headerfile, not the dataset

headerformat = ft_filetype(filename); % update the filetype

end

end % if skip initial check

% implement the caching in a data-format independent way

if cache && exist(headerfile, 'file') && ~isempty(cacheheader)

% try to get the header from cache

details = dir(headerfile);

if isequal(details, cacheheader.details)

% the header file has not been updated, fetch it from the cache

% fprintf('got header from cache\n');

hdr = rmfield(cacheheader, 'details');

switch ft_filetype(datafile)

case {'ctf_ds' 'ctf_meg4' 'ctf_old' 'read_ctf_res4'}

% for realtime analysis end-of-file-chasing the res4 does not correctly

% estimate the number of samples, so we compute it on the fly

sz = 0;

files = dir([filename filesep '*.*meg4']);

for j=1:numel(files)

sz = sz + files(j).bytes;

end

hdr.nTrials = floor((sz - 8) / (hdr.nChans*4) / hdr.nSamples);

end

return

end % if the details correspond

end % if cache

% the support for head/dewar coordinates is still limited

if strcmp(coordsys, 'dewar') && ~any(strcmp(headerformat, {'fcdc_buffer', 'ctf_ds', 'ctf_meg4', 'ctf_res4', 'neuromag_fif', 'neuromag_mne'}))

ft_error('dewar coordinates are not supported for %s', headerformat);

end

% deal with data that is organized according to BIDS

if strcmp(readbids, 'yes') || strcmp(readbids, 'ifmakessense')

[p, f, x] = fileparts(filename);

% check whether it is a BIDS dataset with json and tsv sidecar files

% data in a BIDS tsv file (like physio and stim) will be explicitly dealt with in BIDS_TSV

isbids = startsWith(f, 'sub-') && ~strcmp(x, '.tsv');

if isbids

% try to read the metadata from the BIDS sidecar files

sidecar = bids_sidecar(filename);

if ~isempty(sidecar)

data_json = ft_read_json(sidecar);

end

sidecar = bids_sidecar(filename, 'channels');

if ~isempty(sidecar)

channels_tsv = ft_read_tsv(sidecar);

end

sidecar = bids_sidecar(filename, 'electrodes');

if ~isempty(sidecar)

electrodes_tsv = ft_read_tsv(sidecar);

end

sidecar = bids_sidecar(filename, 'optodes');

if ~isempty(sidecar)

optodes_tsv = ft_read_tsv(sidecar);

end

sidecar = bids_sidecar(filename, 'coordsystem');

if ~isempty(sidecar)

coordsystem_json = ft_read_json(sidecar);

end

end

end

% start with an empty header

hdr = [];

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

% read the data with the low-level reading function

% please maintain this list in alphabetical order

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

switch headerformat

case '4d'

orig = read_4d_hdr(datafile);

hdr.Fs = orig.header_data.SampleFrequency;

hdr.nChans = orig.header_data.TotalChannels;

hdr.nSamples = orig.header_data.SlicesPerEpoch;

hdr.nSamplesPre = round(orig.header_data.FirstLatency*orig.header_data.SampleFrequency);

hdr.nTrials = orig.header_data.TotalEpochs;

%hdr.label = {orig.channel_data(:).chan_label}';

hdr.label = orig.Channel;

[hdr.grad, elec] = bti2grad(orig);

if ~isempty(elec)

hdr.elec = elec;

end

% remember original header details

hdr.orig = orig;

case {'4d_pdf', '4d_m4d', '4d_xyz'}

orig = read_bti_m4d(filename);

hdr.Fs = orig.SampleFrequency;

hdr.nChans = orig.TotalChannels;

hdr.nSamples = orig.SlicesPerEpoch;

hdr.nSamplesPre = round(orig.FirstLatency*orig.SampleFrequency);

hdr.nTrials = orig.TotalEpochs;

hdr.label = orig.ChannelOrder(:);

[hdr.grad, elec] = bti2grad(orig);

if ~isempty(elec)

hdr.elec = elec;

end

% remember original header details

hdr.orig = orig;

case 'AnyWave'

orig = read_ahdf5_hdr(datafile);

hdr.orig = orig;

hdr.Fs = orig.channels(1).samplingRate;

hdr.nChans = numel(orig.channels);

hdr.nSamples = orig.numberOfSamples;

hdr.nTrials = orig.numberOfBlocks;

hdr.nSamplesPre = 0;

hdr.label = orig.label;

hdr.reference = orig.reference(:);

hdr.chanunit = orig.unit(:);

hdr.chantype = orig.type(:);

case 'bci2000_dat'

% this requires the load_bcidat mex file to be present on the path

ft_hastoolbox('BCI2000', 1);

% this is inefficient, since it reads the complete data

[signal, states, parameters, total_samples] = load_bcidat(filename);

% convert into a FieldTrip-like header

hdr = [];

hdr.nChans = size(signal,2);

hdr.nSamples = total_samples;

hdr.nSamplesPre = 0; % it is continuous

hdr.nTrials = 1; % it is continuous

hdr.Fs = parameters.SamplingRate.NumericValue;

% there are some differences in the fields that are present in the

% *.dat files, probably due to different BCI2000 versions

if isfield(parameters, 'ChannelNames') && isfield(parameters.ChannelNames, 'Value') && ~isempty(parameters.ChannelNames.Value)

hdr.label = parameters.ChannelNames.Value;

elseif isfield(parameters, 'ChannelNames') && isfield(parameters.ChannelNames, 'Values') && ~isempty(parameters.ChannelNames.Values)

hdr.label = parameters.ChannelNames.Values;

else

% give this warning only once

ft_warning('creating fake channel names');

for i=1:hdr.nChans

hdr.label{i} = sprintf('%d', i);

end

end

% remember the original header details

hdr.orig.parameters = parameters;

% also remember the complete data upon request

if cache

hdr.orig.signal = signal;

hdr.orig.states = states;

hdr.orig.total_samples = total_samples;

end

case 'besa_besa'

if isempty(chanindx)

hdr = read_besa_besa(filename);

else

% this will keep the non-selected channels hidden from the user

hdr = read_besa_besa(filename, [], chanindx);

end

case 'besa_avr'

orig = read_besa_avr(filename);

hdr.Fs = 1000/orig.di;

hdr.nChans = size(orig.data,1);

hdr.nSamples = size(orig.data,2);

hdr.nSamplesPre = -(hdr.Fs * orig.tsb/1000); % convert from ms to samples

hdr.nTrials = 1;

if isfield(orig, 'label') && iscell(orig.label)

hdr.label = orig.label;

elseif isfield(orig, 'label') && ischar(orig.label)

hdr.label = tokenize(orig.label, ' ');

else

% give this warning only once

ft_warning('creating fake channel names');

for i=1:hdr.nChans

hdr.label{i} = sprintf('%d', i);

end

end

case 'besa_swf'

orig = read_besa_swf(filename);

hdr.Fs = 1000/orig.di;

hdr.nChans = size(orig.data,1);

hdr.nSamples = size(orig.data,2);

hdr.nSamplesPre = -(hdr.Fs * orig.tsb/1000); % convert from ms to samples

hdr.nTrials = 1;

hdr.label = orig.label;

case 'biosig'

% this requires the biosig toolbox

ft_hastoolbox('BIOSIG', 1);

hdr = read_biosig_header(filename);

case {'biosemi_bdf', 'bham_bdf'}

hdr = read_biosemi_bdf(filename);

if any(diff(hdr.orig.SampleRate))

ft_error('channels with different sampling rate not supported');

end

if ~ft_senstype(hdr, 'ext1020')

% assign the channel type and units for the known channels

hdr.chantype = repmat({'unknown'}, size(hdr.label));

hdr.chanunit = repmat({'unknown'}, size(hdr.label));

chan = ~cellfun(@isempty, regexp(hdr.label, '^[A-D]\d*$'));

hdr.chantype(chan) = {'eeg'};

hdr.chanunit(chan) = {'uV'};

end

if ft_filetype(filename, 'bham_bdf')

% TODO channel renaming should be made a general option

% this is for the Biosemi system used at the University of Birmingham

labelold = { 'A1' 'A2' 'A3' 'A4' 'A5' 'A6' 'A7' 'A8' 'A9' 'A10' 'A11' 'A12' 'A13' 'A14' 'A15' 'A16' 'A17' 'A18' 'A19' 'A20' 'A21' 'A22' 'A23' 'A24' 'A25' 'A26' 'A27' 'A28' 'A29' 'A30' 'A31' 'A32' 'B1' 'B2' 'B3' 'B4' 'B5' 'B6' 'B7' 'B8' 'B9' 'B10' 'B11' 'B12' 'B13' 'B14' 'B15' 'B16' 'B17' 'B18' 'B19' 'B20' 'B21' 'B22' 'B23' 'B24' 'B25' 'B26' 'B27' 'B28' 'B29' 'B30' 'B31' 'B32' 'C1' 'C2' 'C3' 'C4' 'C5' 'C6' 'C7' 'C8' 'C9' 'C10' 'C11' 'C12' 'C13' 'C14' 'C15' 'C16' 'C17' 'C18' 'C19' 'C20' 'C21' 'C22' 'C23' 'C24' 'C25' 'C26' 'C27' 'C28' 'C29' 'C30' 'C31' 'C32' 'D1' 'D2' 'D3' 'D4' 'D5' 'D6' 'D7' 'D8' 'D9' 'D10' 'D11' 'D12' 'D13' 'D14' 'D15' 'D16' 'D17' 'D18' 'D19' 'D20' 'D21' 'D22' 'D23' 'D24' 'D25' 'D26' 'D27' 'D28' 'D29' 'D30' 'D31' 'D32' 'EXG1' 'EXG2' 'EXG3' 'EXG4' 'EXG5' 'EXG6' 'EXG7' 'EXG8' 'Status'};

labelnew = { 'P9' 'PPO9h' 'PO7' 'PPO5h' 'PPO3h' 'PO5h' 'POO9h' 'PO9' 'I1' 'OI1h' 'O1' 'POO1' 'PO3h' 'PPO1h' 'PPO2h' 'POz' 'Oz' 'Iz' 'I2' 'OI2h' 'O2' 'POO2' 'PO4h' 'PPO4h' 'PO6h' 'POO10h' 'PO10' 'PO8' 'PPO6h' 'PPO10h' 'P10' 'P8' 'TPP9h' 'TP7' 'TTP7h' 'CP5' 'TPP7h' 'P7' 'P5' 'CPP5h' 'CCP5h' 'CP3' 'P3' 'CPP3h' 'CCP3h' 'CP1' 'P1' 'Pz' 'CPP1h' 'CPz' 'CPP2h' 'P2' 'CPP4h' 'CP2' 'CCP4h' 'CP4' 'P4' 'P6' 'CPP6h' 'CCP6h' 'CP6' 'TPP8h' 'TP8' 'TPP10h' 'T7' 'FTT7h' 'FT7' 'FC5' 'FCC5h' 'C5' 'C3' 'FCC3h' 'FC3' 'FC1' 'C1' 'CCP1h' 'Cz' 'FCC1h' 'FCz' 'FFC1h' 'Fz' 'FFC2h' 'FC2' 'FCC2h' 'CCP2h' 'C2' 'C4' 'FCC4h' 'FC4' 'FC6' 'FCC6h' 'C6' 'TTP8h' 'T8' 'FTT8h' 'FT8' 'FT9' 'FFT9h' 'F7' 'FFT7h' 'FFC5h' 'F5' 'AFF7h' 'AF7' 'AF5h' 'AFF5h' 'F3' 'FFC3h' 'F1' 'AF3h' 'Fp1' 'Fpz' 'Fp2' 'AFz' 'AF4h' 'F2' 'FFC4h' 'F4' 'AFF6h' 'AF6h' 'AF8' 'AFF8h' 'F6' 'FFC6h' 'FFT8h' 'F8' 'FFT10h' 'FT10'};

% rename the channel labels

for i=1:length(labelnew)

chan = strcmp(labelold(i), hdr.label);

hdr.label(chan) = labelnew(chan);

end

end

case {'biosemi_old'}

% this uses the openbdf and readbdf functions that were copied from EEGLAB

orig = openbdf(filename);

if any(orig.Head.SampleRate~=orig.Head.SampleRate(1))

ft_error('channels with different sampling rate not supported');

end

hdr.Fs = orig.Head.SampleRate(1);

hdr.nChans = orig.Head.NS;

hdr.label = cellstr(orig.Head.Label);

% it is continuous data, therefore append all records in one trial

hdr.nSamples = orig.Head.NRec * orig.Head.Dur * orig.Head.SampleRate(1);

hdr.nSamplesPre = 0;

hdr.nTrials = 1;

hdr.orig = orig;

% close the file between separate read operations

fclose(orig.Head.FILE.FID);

case 'blackrock_nev'

% read header for nsX file associated with NEV file

% use ft_read_event to read event information in .nev file

ft_hastoolbox('NPMK', 1);

% ensure that the filename contains a full path specification,

% otherwise the low-level function fails

[p,n] = fileparts(filename);

if isempty(p)

filename = which(filename);

[p,n] = fileparts(filename);

end

NEV = openNEV(filename,'noread','nosave');

% searching for associated nsX file in same folder

files=dir(strcat(fullfile(p,n),'.ns*'));

if isempty(files)

ft_error('no .ns* file associated to %s in %s',n,p);

end

%searching for nsX file with same sampling freq that NEV

for i=1:numel(files)

nsX_hdr = ft_read_header(fullfile(p,files(i).name),'chantype',chantype);

if nsX_hdr.Fs == NEV.MetaTags.SampleRes

hdr = nsX_hdr;

break

end

end

if isempty(hdr)

ft_error('no .ns* file with same sampling frequency as %s (%i)',n,NEV.MetaTags.SampleRes);

end

case 'blackrock_nsx'

ft_hastoolbox('NPMK', 1);

% ensure that the filename contains a full path specification,

% otherwise the low-level function fails

p = fileparts(filename);

if isempty(p)

filename = which(filename);

end

orig = openNSx(filename, 'noread');

channelstype=regexp({orig.ElectrodesInfo.Label},'[A-Za-z]+','match','once');

chaninfo=table({orig.ElectrodesInfo.ElectrodeID}',...

transpose(deblank({orig.ElectrodesInfo.Label})),[channelstype]',...

{orig.ElectrodesInfo.ConnectorBank}',{orig.ElectrodesInfo.ConnectorPin}',...

transpose(deblank({orig.ElectrodesInfo.AnalogUnits})),...

'VariableNames',{'id' 'label' 'chantype' 'bank' 'pin' 'unit'});

if isempty(chantype)

chantype = unique(channelstype,'stable');

end

% selecting channel according to chantype

orig_label=deblank({orig.ElectrodesInfo.Label});

orig_unit=deblank({orig.ElectrodesInfo.AnalogUnits});

channels={}; channelstype={}; channelsunit={}; skipfactor=[];

for c=1:length(chantype)

chantype_split=strsplit(chantype{c},':');

if numel(chantype_split) == 2

chantype{c}=chantype_split{1};

skipfactor=[skipfactor,str2double(chantype_split{2})];

elseif numel(chantype_split) > 2

ft_error('Use : to specify skipfactor, e.g. analog:10')

end

chan_sel=~cellfun(@isempty,regexp(orig_label,chantype{c}));

if sum(chan_sel)==0

if ~strcmp(chantype{c},'chaninfo')

ft_error('unknown chantype %s, available channels are %s',chantype{c},strjoin(orig_label));

end

else

channels=[channels, orig_label(chan_sel)];

channelsunit=[channelsunit, orig_unit(chan_sel)];

channelstype=[channelstype, repmat(chantype(c), 1, sum(chan_sel))];

end

end

skipfactor=unique(skipfactor);

if isempty(skipfactor)

skipfactor=1;

elseif length(skipfactor)>1

ft_error('inconsistent skip factors across channels');

end

% If no channel selected issue error specifying available chantypes

if isempty(channels)

ft_error('No channel selected. Availabe chantypes are: %s', strjoin(unique(chaninfo.chantype)));

end

hdr.Fs = orig.MetaTags.SamplingFreq/skipfactor;

hdr.nChans = length(channels);

hdr.nSamples = floor(orig.MetaTags.DataPoints/skipfactor);

hdr.nSamplesPre = 0;

hdr.nTrials = 1; %?

hdr.label = deblank(channels)';

hdr.chantype = channelstype;

hdr.chanunit = channelsunit;

hdr.orig = orig;

hdr.orig.chaninfo = chaninfo;

hdr.orig.skipfactor = skipfactor;

case {'brainvision_vhdr', 'brainvision_seg', 'brainvision_eeg', 'brainvision_dat'}

orig = read_brainvision_vhdr(filename);

hdr.Fs = orig.Fs;

hdr.nChans = orig.NumberOfChannels;

hdr.label = orig.label;

hdr.nSamples = orig.nSamples;

hdr.nSamplesPre = orig.nSamplesPre;

hdr.nTrials = orig.nTrials;

hdr.orig = orig;

% assign the channel type and units for the known channels

hdr.chantype = repmat({'eeg'}, size(hdr.label));

hdr.chanunit = repmat({'uV'}, size(hdr.label));

case 'bucn_nirs'

orig = read_bucn_nirshdr(filename);

hdr = rmfield(orig, 'time');

hdr.orig = orig;

case 'ced_son'

% check that the required low-level toolbox is available

ft_hastoolbox('neuroshare', 1);

% use the reading function supplied by Gijs van Elswijk

orig = read_ced_son(filename,'readevents','no','readdata','no');

orig = orig.header;

% In Spike2, channels can have different sampling rates, units, length

% etc. etc. Here, channels need to have to same properties.

if length(unique([orig.samplerate]))>1,

ft_error('channels with different sampling rates are not supported');

else

hdr.Fs = orig(1).samplerate;

end

hdr.nChans = length(orig);

% nsamples of the channel with least samples

hdr.nSamples = min([orig.nsamples]);

hdr.nSamplesPre = 0;

% only continuous data supported

if sum(strcmpi({orig.mode},'continuous')) < hdr.nChans,

ft_error('not all channels contain continuous data');

else

hdr.nTrials = 1;

end

hdr.label = {orig.label};

case 'combined_ds'

hdr = read_combined_ds(filename);

case {'ctf_ds', 'ctf_meg4', 'ctf_res4'}

% check the presence of the required low-level toolbox

ft_hastoolbox('ctf', 1);

% default for CTF is to remove the site-specific numbers from each channel name, e.g. 'MZC01-1706' becomes 'MZC01'

splitlabel = ft_getopt(varargin, 'splitlabel', true);

orig = readCTFds(filename);

if isempty(orig)

% this is to deal with data from the 64 channel system and the error

% readCTFds: .meg4 file header=MEG4CPT Valid header options: MEG41CP MEG42CP

ft_error('could not read CTF with this implementation, please try again with the ''ctf_old'' file format');

end

hdr.Fs = orig.res4.sample_rate;

hdr.nChans = orig.res4.no_channels;

hdr.nSamples = orig.res4.no_samples;

hdr.nSamplesPre = orig.res4.preTrigPts;

hdr.nTrials = orig.res4.no_trials;

hdr.label = cellstr(orig.res4.chanNames);

% read the balance coefficients, these are used to compute the synthetic gradients

try

coeftype = cellstr(char(orig.res4.scrr(:).coefType));

catch

coeftype = {};

end

try

[alphaMEG,MEGlist,Refindex] = getCTFBalanceCoefs(orig,'NONE', 'T');

orig.BalanceCoefs.none.alphaMEG = alphaMEG;

orig.BalanceCoefs.none.MEGlist = MEGlist;

orig.BalanceCoefs.none.Refindex = Refindex;

catch

ft_warning('cannot read balancing coefficients for NONE');

end

if any(~cellfun(@isempty,strfind(coeftype, 'G1BR')))

try

[alphaMEG,MEGlist,Refindex] = getCTFBalanceCoefs(orig,'G1BR', 'T');

orig.BalanceCoefs.G1BR.alphaMEG = alphaMEG;

orig.BalanceCoefs.G1BR.MEGlist = MEGlist;

orig.BalanceCoefs.G1BR.Refindex = Refindex;

catch

ft_warning('cannot read balancing coefficients for G1BR');

end

end

if any(~cellfun(@isempty,strfind(coeftype, 'G2BR')))

try

[alphaMEG,MEGlist,Refindex] = getCTFBalanceCoefs(orig, 'G2BR', 'T');

orig.BalanceCoefs.G2BR.alphaMEG = alphaMEG;

orig.BalanceCoefs.G2BR.MEGlist = MEGlist;

orig.BalanceCoefs.G2BR.Refindex = Refindex;

catch

ft_warning('cannot read balancing coefficients for G2BR');

end

end

if any(~cellfun(@isempty,strfind(coeftype, 'G3BR')))

try

[alphaMEG,MEGlist,Refindex] = getCTFBalanceCoefs(orig, 'G3BR', 'T');

orig.BalanceCoefs.G3BR.alphaMEG = alphaMEG;

orig.BalanceCoefs.G3BR.MEGlist = MEGlist;

orig.BalanceCoefs.G3BR.Refindex = Refindex;

catch

ft_warning('cannot read balancing coefficients for G3BR');

end

end

if any(~cellfun(@isempty,strfind(coeftype, 'G3AR')))

try

[alphaMEG,MEGlist,Refindex] = getCTFBalanceCoefs(orig, 'G3AR', 'T');

orig.BalanceCoefs.G3AR.alphaMEG = alphaMEG;

orig.BalanceCoefs.G3AR.MEGlist = MEGlist;

orig.BalanceCoefs.G3AR.Refindex = Refindex;

catch

% May not want a warning here if these are not commonly used.

% Already get a (fprintf) warning from getCTFBalanceCoefs.m

% ft_warning('cannot read balancing coefficients for G3AR');

end

end

% add a gradiometer structure for forward and inverse modelling

try

[grad, elec] = ctf2grad(orig, strcmp(coordsys, 'dewar'), coilaccuracy, coildeffile);

if ~isempty(grad)

hdr.grad = grad;

end

if ~isempty(elec)

hdr.elec = elec;

end

catch

% this fails if the res4 file is not correctly closed, e.g. during realtime processing

tmp = lasterror;

disp(tmp.message);

ft_warning('could not construct gradiometer definition from the header');

end

% for realtime analysis EOF chasing the res4 does not correctly

% estimate the number of samples, so we compute it on the fly from the

% meg4 file sizes.

sz = 0;

files = dir([filename '/*.*meg4']);

for j=1:numel(files)

sz = sz + files(j).bytes;

end

hdr.nTrials = floor((sz - 8) / (hdr.nChans*4) / hdr.nSamples);

% add the original header details

hdr.orig = orig;

case {'ctf_old', 'read_ctf_res4'}

% read it using the open-source MATLAB code that originates from CTF and that was modified by the FCDC

orig = read_ctf_res4(headerfile);

hdr.Fs = orig.Fs;

hdr.nChans = orig.nChans;

hdr.nSamples = orig.nSamples;

hdr.nSamplesPre = orig.nSamplesPre;

hdr.nTrials = orig.nTrials;

hdr.label = orig.label;

% add a gradiometer structure for forward and inverse modelling

try

hdr.grad = ctf2grad(orig);

catch

% this fails if the res4 file is not correctly closed, e.g. during realtime processing

tmp = lasterror;

disp(tmp.message);

ft_warning('could not construct gradiometer definition from the header');

end

% add the original header details

hdr.orig = orig;

case 'ctf_read_res4'

% check that the required low-level toolbos ix available

ft_hastoolbox('eegsf', 1);

% read it using the CTF importer from the NIH and Daren Weber

orig = ctf_read_res4(fileparts(headerfile), 0);

% convert the header into a structure that FieldTrip understands

hdr = [];

hdr.Fs = orig.setup.sample_rate;

hdr.nChans = length(orig.sensor.info);

hdr.nSamples = orig.setup.number_samples;

hdr.nSamplesPre = orig.setup.pretrigger_samples;

hdr.nTrials = orig.setup.number_trials;

for i=1:length(orig.sensor.info)

hdr.label{i} = orig.sensor.info(i).label;

end

hdr.label = hdr.label(:);

% add a gradiometer structure for forward and inverse modelling

try

hdr.grad = ctf2grad(orig);

catch

% this fails if the res4 file is not correctly closed, e.g. during realtime processing

tmp = lasterror;

disp(tmp.message);

ft_warning('could not construct gradiometer definition from the header');

end

% add the original header details

hdr.orig = orig;

case 'ctf_shm'

% read the header information from shared memory

hdr = read_shm_header(filename);

case {'curry_dat', 'curry_cdt'}

orig = load_curry_data_file(filename);

hdr = [];

hdr.Fs = orig.fFrequency;

hdr.nChans = orig.nChannels;

hdr.nSamples = orig.nSamples;

hdr.nSamplesPre = sum(orig.time<0);

hdr.nTrials = orig.nTrials;

hdr.label = orig.labels(:);

hdr.orig = orig;

case 'dataq_wdq'

orig = read_wdq_header(filename);

hdr = [];

hdr.Fs = orig.fsample;

hdr.nChans = orig.nchan;

hdr.nSamples = orig.nbytesdat/(2*hdr.nChans);

hdr.nSamplesPre = 0;

hdr.nTrials = 1;

for k = 1:hdr.nChans

if isfield(orig.chanhdr(k), 'annot') && ~isempty(orig.chanhdr(k).annot)

hdr.label{k,1} = orig.chanhdr(k).annot;

else

hdr.label{k,1} = orig.chanhdr(k).label;

end

end

% add the original header details

hdr.orig = orig;

case {'deymed_ini' 'deymed_dat'}

% the header is stored in a *.ini file

orig = read_deymed_ini(headerfile);

hdr = [];

hdr.Fs = orig.Fs;

hdr.nChans = orig.nChans;

hdr.nSamples = orig.nSamples;

hdr.nSamplesPre = 0;

hdr.nTrials = 1;

hdr.label = orig.label(:);

hdr.orig = orig; % remember the original details

case 'dhn_med10'

ft_hastoolbox('mayo_mef', 1); % make sure mayo_mef exists

hdr = read_dhn_med10(filename, password);

case 'edf'

% this reader is largely similar to the bdf reader

if isempty(chanindx)

hdr = read_edf(filename);

else

% this will keep the non-selected channels hidden from the user

hdr = read_edf(filename, [], chanindx);

end

case 'eep_avr'

% check that the required low-level toolbox is available

ft_hastoolbox('eeprobe', 1);

% read the whole average and keep only header info (it is a bit silly, but the easiest to do here)

hdr = read_eep_avr(filename);

hdr.Fs = hdr.rate;

hdr.nChans = size(hdr.data,1);

hdr.nSamples = size(hdr.data,2);

hdr.nSamplesPre = hdr.xmin*hdr.rate/1000;

hdr.nTrials = 1; % it can always be interpreted as continuous data

% remove the data and variance if present

hdr = removefields(hdr, {'data', 'variance'});

case 'eep_cnt'

% check that the required low-level toolbox is available

ft_hastoolbox('eeprobe', 1);

% read the first sample from the continuous data, this will also return the header

orig = read_eep_cnt(filename, 1, 1);

hdr.Fs = orig.rate;

hdr.nSamples = orig.nsample;

hdr.nSamplesPre = 0;

hdr.label = orig.label;

hdr.nChans = orig.nchan;

hdr.nTrials = 1; % it can always be interpreted as continuous data

hdr.orig = orig; % remember the original details

case 'eeglab_set'

hdr = read_eeglabheader(filename);

case 'eeglab_erp'

hdr = read_erplabheader(filename);

case 'emotiv_mat'

% This is a MATLAB *.mat file that is created using the Emotiv MATLAB

% example code. It contains a 25xNsamples matrix and some other stuff.

orig = load(filename);

hdr.Fs = 128;

hdr.nChans = 25;

hdr.nSamples = size(orig.data_eeg,1);

hdr.nSamplesPre = 0;

hdr.nTrials = 1;

hdr.label = {'ED_COUNTER','ED_INTERPOLATED','ED_RAW_CQ','ED_AF3','ED_F7','ED_F3','ED_FC5','ED_T7','ED_P7','ED_O1','ED_O2','ED_P8','ED_T8','ED_FC6','ED_F4','ED_F8','ED_AF4','ED_GYROX','ED_GYROY','ED_TIMESTAMP','ED_ES_TIMESTAMP','ED_FUNC_ID','ED_FUNC_VALUE','ED_MARKER','ED_SYNC_SIGNAL'};

% store the complete information in hdr.orig

% ft_read_data and ft_read_event will get it from there

hdr.orig = orig;

case 'eyelink_asc'

asc = read_eyelink_asc(filename);

hdr.nChans = size(asc.dat,1);

hdr.nSamples = size(asc.dat,2);

hdr.nSamplesPre = 0;

hdr.nTrials = 1;

hdr.FirstTimeStamp = asc.dat(1,1);

hdr.TimeStampPerSample = median(diff(asc.dat(1,:)));

hdr.Fs = 1000/hdr.TimeStampPerSample; % these timestamps are in miliseconds

% give this warning only once

ft_warning('creating fake channel names');

for i=1:hdr.nChans

hdr.label{i} = sprintf('%d', i);

end

% remember all header and data details upon request

if cache

hdr.orig = asc;

else

% remember the original header details

hdr.orig = removefields(asc, 'dat');

end

case 'spmeeg_mat'

hdr = read_spmeeg_header(filename);

case 'ced_spike6mat'

hdr = read_spike6mat_header(filename);

case 'egi_egia'

[fhdr,chdr,ename,cnames,fcom,ftext] = read_egis_header(filename);

[p, f, x] = fileparts(filename);

if any(chdr(:,4)-chdr(1,4))

ft_error('Sample rate not the same for all cells.');

end

hdr.Fs = chdr(1,4); %making assumption that sample rate is same for all cells

hdr.nChans = fhdr(19);

for i = 1:hdr.nChans

% this should be consistent with ft_senslabel

hdr.label{i,1} = ['E' num2str(i)];

end

%since NetStation does not properly set the fhdr(11) field, use the number of subjects from the chdr instead

hdr.nTrials = chdr(1,2)*fhdr(18); %number of trials is numSubjects * numCells

hdr.nSamplesPre = ceil(fhdr(14)/(1000/hdr.Fs));

if any(chdr(:,3)-chdr(1,3))

ft_error('Number of samples not the same for all cells.');

end

hdr.nSamples = chdr(1,3); %making assumption that number of samples is same for all cells

% remember the original header details

hdr.orig.fhdr = fhdr;

hdr.orig.chdr = chdr;

hdr.orig.ename = ename;

hdr.orig.cnames = cnames;

hdr.orig.fcom = fcom;

hdr.orig.ftext = ftext;

case 'egi_egis'

[fhdr,chdr,ename,cnames,fcom,ftext] = read_egis_header(filename);

[p, f, x] = fileparts(filename);

if any(chdr(:,4)-chdr(1,4))

ft_error('Sample rate not the same for all cells.');

end

hdr.Fs = chdr(1,4); %making assumption that sample rate is same for all cells

hdr.nChans = fhdr(19);

for i = 1:hdr.nChans

% this should be consistent with ft_senslabel

hdr.label{i,1} = ['E' num2str(i)];

end

hdr.nTrials = sum(chdr(:,2));

hdr.nSamplesPre = ceil(fhdr(14)/(1000/hdr.Fs));

% assuming that a utility was used to insert the correct baseline

% duration into the header since it is normally absent. This slot is

% actually allocated to the age of the subject, although NetStation

% does not use it when generating an EGIS session file.