sload.m

Go to the documentation of this file.

function [signal,H] = sload(FILENAME,varargin)
% SLOAD loads signal data of various data formats
% [signal,header] = sload(FILENAME)
% [signal,header] = sload(FILENAME,CHAN)
%       read selected channels in list CHAN
%       CHAN=0 [default], reads all channels
%
% [signal,header] = sload(dir('f*.emg') ... )
% [signal,header] = sload('f*.emg' ...)
%   loads all files 'f*.emg'
%
% [signal,header] = sload(..., PropertyName1,PropertyValue1,...)
%   PropertyName(s)     PropertyValue
%   'UCAL'          -       data uncalibrated (not scaled)
%               'On'        data uncalibrated (not scaled)
%               'Off'       data calibrated (scaled) [default]
%   'OVERFLOWDETECTION' 'On'        [default] 
%               'Off'       no overflow detection 
%   'OUTPUT'        'single'    single precision data [default: 'double'] 
%   'NUMBER_OF_NAN_IN_BREAK'   N        inserts N NaN's between two concatanated segments
%                       default: N=100
%   'SampleRate'        Fs      target sampling rate (supports resampling)
%   'EOG_CORRECTION'    'On'        uses two-channel regression analysis - if possible 
%               'Off'       no correction of EOG artifacts [default]
%       'CNT32', '32bit'            force CNT 32bit format 
%   'BDF'           N       decodes BDF status channel into eventtable 
%       using Mode = N. For details see HELP BDF2BIOSIG_EVENTS
%   
% The list of supported formats is available here: 
% http://pub.ist.ac.at/~schloegl/biosig/TESTED
%
%    SLOAD loads all the data (of the selected channels)
%    at once. In case of large data files, This can be 
%    a problem. Instead, You can use 
%       HDR = sopen(FILENAME,'r');
%       [s,HDR]=sread(HDR,duration_segment1); 
%       [s,HDR]=sread(HDR,duration_segment2); 
%       ....
%       [s,HDR]=sread(HDR,duration_segmentM); 
%       HDR = sclose(HDR); 
%
% see also: SVIEW, SOPEN, SREAD, SCLOSE, SAVE2BKR, BDF2BIOSIG_EVENTS, TLOAD
%
% In order to increase the speed, install mexSLOAD.mex from biosig4c++
%
% Reference(s):
 
 
%   $Id: sload.m 3169 2012-12-04 15:39:03Z schloegl $
%   Copyright (C) 1997-2007,2008,2009,2010,2011,2012 by Alois Schloegl 
%       This is part of the BIOSIG-toolbox http://biosig.sf.net/
 
%
%    BioSig 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.
%
%    BioSig 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 BioSig.  If not, see <http://www.gnu.org/licenses/>.
 
if length(varargin)<2; 
    MODE = ''; 
else    
    MODE = varargin{2};
end;
 
CHAN = 0; 
STATE.CNT32 = 0; 
STATE.UCAL = 0; 
STATE.OVERFLOWDETECTION = 1; 
STATE.NUMBER_OF_NAN_IN_BREAK = 100;     %% default value, always include when several files are concatanated
STATE.FLAG_NUM_NAN   = 0;       %% flags whether explicit argument NUMBER_OF_NAN_IN_BREAK is used, only then NaNs are introduced between sweeps(segments) within single files
STATE.EOG_CORRECTION = 0 ; 
STATE.BDF = [] ; 
STATE.OUTPUT = 'double';
 
Fs = NaN; 
k = 1; 
while (k<=length(varargin))
    if isnumeric(varargin{k})
        if (k==1), CHAN=varargin{k}; end;
    elseif ~isempty(strfind(varargin{k},'UCAL'))
        MODE = varargin{k}; 
        if strcmpi(varargin{k+1},'on'); 
            STATE.UCAL = 1; 
            k=k+1;
        elseif strcmpi(varargin{k+1},'off'); 
            STATE.UCAL = 0; 
            k=k+1;
        else
            STATE.UCAL = 1; 
        end; 
    elseif ~isempty(strfind(varargin{k},'OVERFLOWDETECTION:OFF'))
        MODE = varargin{k}; 
        STATE.OVERFLOWDETECTION = 0;
    elseif strcmpi(varargin{k},'OVERFLOWDETECTION')
        if strcmpi(varargin{k+1},'on'); 
            STATE.OVERFLOWDETECTION = 1;
        elseif strcmpi(varargin{k+1},'off'); 
            STATE.OVERFLOWDETECTION = 0;
        elseif isnumeric(varargin{k+1})
            STATE.OVERFLOWDETECTION = varargin{k+1};
        end;            
        k=k+1;
    elseif ~isempty(strfind(varargin{k},'OUTPUT:SINGLE'))
        MODE = varargin{k}; 
        STATE.OUTPUT = 'single';
    elseif strcmpi(varargin{k},'OUTPUT')
        STATE.OUTPUT = varargin{k+1};
    elseif strcmpi(varargin{k},'EOG_CORRECTION')
        if strcmpi(varargin{k+1},'on'); 
            STATE.EOG_CORRECTION = 1;
        elseif strcmpi(varargin{k+1},'off'); 
            STATE.EOG_CORRECTION = 0;
        elseif isnumeric(varargin{k+1})
            STATE.EOG_CORRECTION = varargin{k+1};
        end;            
        k=k+1;
    elseif strcmpi(varargin{k},'Number_of_nan_in_break')
        STATE.NUMBER_OF_NAN_IN_BREAK = varargin{k+1};
        STATE.FLAG_NUM_NAN = 1; 
        k=k+1;
    elseif strcmpi(varargin{k},'SampleRate')
        Fs = varargin{k+1};
        k=k+1;
    elseif strcmpi(varargin{k},'BDF')
        STATE.BDF = varargin{k+1};
        k=k+1;
    elseif strcmpi(varargin{k},'CNT32') || strcmpi(varargin{k},'32bit')  
        STATE.CNT32 = 1;
    elseif exist(varargin{k},'file') && (k==1)
            CHAN = varargin{k}; 
    end;    
    k=k+1;
end;
 
 
%%% resolve wildcards %%%
if ischar(FILENAME) 
if any(FILENAME=='*')
        p = fileparts(FILENAME);
        f = dir(FILENAME);
        EOGix = zeros(1,length(f));
        for k = 1:length(f);
                f(k).name = fullfile(p,f(k).name);
                [p,g,e]   = fileparts(f(k).name);
                lg = length(g);
                if (lg>2) && strcmp(upper(g(lg+(-2:0))),'EOG')
                        EOGix(k) = 1;
                end
        end;
        FILENAME={f([find(EOGix),find(~EOGix)]).name};
end;        
end;
 
if nargout>0,
    signal = [];
end;
 
if (iscell(FILENAME) && (length(FILENAME)==1)),
    FILENAME = FILENAME{k};
end;    
    
if ((iscell(FILENAME) || isstruct(FILENAME))),
    signal = [];
    for k = 1:length(FILENAME),
        if iscell(FILENAME(k))
            f = FILENAME{k};
        else  
            f = FILENAME(k);
        end 
 
                [s,h] = sload(f,varargin{:});
        if k==1,
            H = h;
            signal = s;  
            H.SegLen = [0,size(s,1)];
                        %H.EVENT.POS = [H.EVENT.POS; 0];
                        %H.EVENT.TYP = [H.EVENT.TYP; hex2dec('7ffe')];
                        %if isfield(H.EVENT,'CHN');
                        %        H.EVENT.CHN = [H.EVENT.CHN; 0];
                        %end;
                        %if isfield(H.EVENT,'DUR');
                        %        H.EVENT.DUR = [H.EVENT.DUR; size(s,1)];
                        %end;
                        %if isfield(H.EVENT,'Desc');    % TFM-Excel-Beat-to-Beat
                        %        H.EVENT.Desc = [H.EVENT.Desc; {'New Segment'}; h.EVENT.Desc];
                        %end;
        else
            H.FILE(k) = h.FILE;
                        H.T0(k,1:6) = h.T0;
            if ~isnan(h.SampleRate) && (H.SampleRate ~= h.SampleRate),
                fprintf(2,'Warning SLOAD: sampling rates of multiple files differ %i!=%i.\n',H.SampleRate, h.SampleRate);
            end;
                        if size(s,2)==size(signal,2), %(H.NS == h.NS) 
                signal = [signal; repmat(NaN,STATE.NUMBER_OF_NAN_IN_BREAK,size(s,2)); s];
                H.SegLen = [H.SegLen,size(signal,1)];
            else
                error('ERROR SLOAD: incompatible channel numbers %i!=%i of multiple files\n',H.NS,h.NS);
            end;
 
                        if ~isempty(h.EVENT.POS),
                                H.EVENT.POS = [H.EVENT.POS; size(signal,1); h.EVENT.POS+size(signal,1)-size(s,1)];
                                H.EVENT.TYP = [H.EVENT.TYP; hex2dec('7ffe'); h.EVENT.TYP];
                                if isfield(H.EVENT,'CHN');
                                H.EVENT.CHN = [H.EVENT.CHN; 0; h.EVENT.CHN];
                                end;
                                if isfield(H.EVENT,'DUR');
                        H.EVENT.DUR = [H.EVENT.DUR; size(s,1); h.EVENT.DUR];
                                end;
                                if isfield(H.EVENT,'Desc'); % TFM-Excel-Beat-to-Beat
                                        H.EVENT.Desc = [H.EVENT.Desc; {'New Segment'}; h.EVENT.Desc];
                                end;
                        end;
 
                        if ~isequal(H.Label,h.Label) || ~isequal(H.PhysDimCode,h.PhysDimCode)
                warning('Labels and PhysDim of multiple files differ!\n');
                                for k2 = 1:length(H.InChanSelect),
                                    k1 = H.InChanSelect(k2); 
                                        if ~strcmp(H.Label{k1},h.Label{k1}) || (H.PhysDimCode(k1)~=h.PhysDimCode(k1)),
                                                warning('#%02i:  %s | %s | (%i)%s | (%i)%s\n',k1, H.Label{k1},h.Label{k1}, H.PhysDimCode(k1),H.PhysDim{k1},h.PhysDimCode(k1),h.PhysDim{k1});
                                        end;
                                end;
                        end;
                        if ~isequal(H.Calib,h.Calib) && H.FLAG.UCAL,
                error('SLOAD: concatanating uncalibrated data with different scaling factors does not make sense!');
                        end;
 
                        if isfield(H,'TriggerOffset'),
                                if H.TriggerOffset ~= h.TriggerOffset,
                                        fprintf(2,'Warning SLOAD: Triggeroffset %f does not fit.%f \n',H.TriggerOffset,h.TriggerOffset);
                                end;
                        end;
 
                        if isfield(h,'TRIG'), 
                                if ~isfield(H,'TRIG'),
                                        H.TRIG = [];
                                end;
                                H.TRIG = [H.TRIG(:); h.TRIG(:)+size(signal,1)-size(s,1)];
                        else h.TRIG=[];        
                        end;
                        
                        if isfield(h,'Classlabel'), 
                if isfield(h,'ArtifactSelection'),
                                    if (any(h.ArtifactSelection>1) || (length(h.ArtifactSelection) < length(h.Classlabel)))
                                            sel = zeros(size(h.Classlabel));
                                                sel(h.ArtifactSelection) = 1; 
                                        else
                                                sel = h.ArtifactSelection(:);
                                        end;
                                else 
                                    h.ArtifactSelection = repmat(logical(0),length(h.Classlabel),1);        
                                end;         
 
                                if isfield(H,'ArtifactSelection')
                                    H.ArtifactSelection = [H.ArtifactSelection; h.ArtifactSelection(:)];
                                else
                                    H.ArtifactSelection = [repmat(logical(0),length(H.Classlabel),1); h.ArtifactSelection(:)];
                                end;
                                H.Classlabel = [H.Classlabel(:); h.Classlabel(:)];
                                
                        elseif isfield(H,'Classlabel') && (length(h.TRIG)>0), 
                    warning(sprintf('field classlabel missing in %s. This might corrupt the resulting classlabels!',h.FileName));
 
                        end;
            clear s;
                end;
    end;
        
    % fprintf(1,'SLOAD: data segments are concatenated with NaNs in between.\n');
    return; 
end;
%%% end of multi-file section 
 
 
%%%% start of single file section
%%%%%%%%%% --------- EOG CORRECTION -------------- %%%%%%
if STATE.EOG_CORRECTION,
try
        h = get_regress_eog(FILENAME,'REG fft16');
catch,
    fprintf(H.FILE.stderr,'Error: SLOAD (EOG_CORRECTION): %s\n',lasterr); 
    H.FLAG.EOG_CORRECTION = 0; 
end;
end; 
 
 
FlagLoaded = 0;
if exist('mexSLOAD','file')==3,
    try
        valid_rerefmx = 1;
        if ischar(CHAN)
                HDR = sopen(CHAN,'r'); HDR=sclose(HDR); 
                if isfield(HDR.Calib),
                        ReRefMx = HDR.Calib; 
                else
                    valid_rerefmx = 0;
                end;         
        elseif all(size(CHAN)>1) || any(floor(CHAN)~=CHAN) || any(CHAN<0) || (any(CHAN==0) && (numel(CHAN)>1));
                ReRefMx = CHAN; 
                CHAN = find(any(CHAN,2));
        elseif all(CHAN>0) && all(floor(CHAN)==CHAN), 
            [tmp,ix]= sort(CHAN);
                ReRefMx = sparse(CHAN,1:length(CHAN),1);
        else
                ReRefMx = [];
            valid_rerefmx=0;
        end
 
        if STATE.EOG_CORRECTION,
            ReRefMx = h.r0*ReRefMx;
            valid_rerefmx=1;
        end
        if STATE.OVERFLOWDETECTION,
            arg1 = 'OVERFLOWDETECTION:ON';
        else
            arg1 = 'OVERFLOWDETECTION:OFF';
        end
        if STATE.UCAL,
            arg2 = {'UCAL:ON'};
        else
            arg2 = {'UCAL:OFF'};
        end
        if STATE.CNT32
            arg2 = {arg2,'CNT32'};
        end
        if ~valid_rerefmx,
            [signal,HDR] = mexSLOAD(FILENAME,0,arg1,arg2{:});
            if isfield(HDR,'ErrNum') && (HDR.ErrNum>0),
                fprintf(1,'%s\n',HDR.ErrMsg);
                H = HDR;
                x = HDR.nonsense;   % hack: triggers CATCH branch in TRY-CATCH-END;
            end  
 
            if isfield(HDR.FLAG,'ROW_BASED_CHANNELS') && HDR.FLAG.ROW_BASED_CHANNELS, signal = signal.'; end;
            FlagLoaded   = isfield(HDR,'NS');
            HDR.InChanSelect = 1:HDR.NS;
        else
            InChanSelect = find(any(ReRefMx,2));
            [signal,HDR] = mexSLOAD(FILENAME,InChanSelect,arg1,arg2{:});
            if isfield(HDR,'ErrNum') && HDR.ErrNum==3,
                %% file not found - fopen failed
                H = HDR;
                return;
            end  
            if isfield(HDR.FLAG,'ROW_BASED_CHANNELS') && HDR.FLAG.ROW_BASED_CHANNELS, signal = signal.'; end;
            FlagLoaded   = isfield(HDR,'NS');
            HDR.InChanSelect = InChanSelect(InChanSelect <= HDR.NS);
            signal = signal*ReRefMx(InChanSelect,:); %% can be sparse if just a single channel is loaded
            signal = full(signal);  %% make sure signal is not sparse 
        end;
        HDR = bdf2biosig_events(HDR, STATE.BDF);
 
        HDR.T0 = datevec(HDR.T0);
        HDR.Patient.Birthday = datevec(HDR.Patient.Birthday);
        HDR.Calib = [HDR.Off(:)';diag(HDR.Cal)];
        [HDR.FILE.Path,HDR.FILE.Name,HDR.FILE.Ext] = fileparts(FILENAME);
        HDR.FileName = FILENAME;
        HDR = leadidcodexyz(HDR);
 
        H=HDR;
        H.FLAG.EOG_CORRECTION = STATE.EOG_CORRECTION; 
        if isfield(HDR,'Patient') && isfield(HDR.Patient,'Weight') && isfield(HDR.Patient,'Height')
            %% Body Mass Index 
            HDR.Patient.BMI = HDR.Patient.Weight * HDR.Patient.Height^-2 * 1e4;
 
                %% Body Surface Area
            % DuBois D, DuBois EF. A formula to estimate the approximate surface area if height and weight be known. Arch Intern Medicine. 1916; 17:863-71.
            % Wang Y, Moss J, Thisted R. Predictors of body surface area. J Clin Anesth. 1992; 4(1):4-10.
                HDR.Patient.BSA = 0.007184 * HDR.Patient.Weight^0.425 * HDR.Patient.Height^0.725;
        end;
 
        if strcmp(H.TYPE,'GDF');
                        % Classlabels according to 
                        % http://biosig.cvs.sourceforge.net/*checkout*/biosig/biosig/doc/eventcodes.txt
 
                        % sort event table before extracting HDR.Classlabel and HDR.TRIG
            [H.EVENT.POS,ix] = sort(H.EVENT.POS);
            H.EVENT.TYP = H.EVENT.TYP(ix);
            if isfield(H.EVENT,'CHN')
                H.EVENT.CHN = H.EVENT.CHN(ix);
            end;            
            if isfield(H.EVENT,'DUR')
                H.EVENT.DUR = H.EVENT.DUR(ix);
            end;    
            if isfield(H.EVENT,'TimeStamp')
                H.EVENT.TimeStamp = H.EVENT.TimeStamp(ix);
            end;    
 
            if (length(H.EVENT.TYP)>0)
                                ix = (H.EVENT.TYP>hex2dec('0300')) & (H.EVENT.TYP<hex2dec('030d'));
                                ix = ix | ((H.EVENT.TYP>=hex2dec('0320')) & (H.EVENT.TYP<=hex2dec('037f')));
                                ix = ix | (H.EVENT.TYP==hex2dec('030f')); % unknown/undefined cue
                                H.Classlabel = mod(H.EVENT.TYP(ix),256);
                                H.Classlabel(H.Classlabel==15) = NaN; % unknown/undefined cue
                        end;
                        
                        % Trigger information and Artifact Selection 
                        ix = find(H.EVENT.TYP==hex2dec('0300')); 
                        H.TRIG = H.EVENT.POS(ix);
                        ArtifactSelection = repmat(logical(0),length(ix),1);
                        for k = 1:length(ix),
                                ix2 = find(H.EVENT.POS(ix(k))==H.EVENT.POS);
                                if any(H.EVENT.TYP(ix2)==hex2dec('03ff'))
                                        ArtifactSelection(k) = logical(1);                
                                end;
                        end;
                        if any(ArtifactSelection), % define only if necessary
                                H.ArtifactSelection = ArtifactSelection; 
                        end;
                        
            % apply channel selections to EVENT table
            if valid_rerefmx && ~isempty(H.EVENT.POS) && (isfield(H.EVENT,'CHN')),  % only if channels are selected. 
                sel = (H.EVENT.CHN(:)==0);  % memory allocation, select all general events
                for k = find(~sel'),        % select channel specific elements
                    sel(k) = any(H.EVENT.CHN(k)==InChanSelect);
                end;
                H.EVENT.POS = H.EVENT.POS(sel);
                H.EVENT.TYP = H.EVENT.TYP(sel);
                H.EVENT.DUR = H.EVENT.DUR(sel); % if EVENT.CHN available, also EVENT.DUR is defined. 
                H.EVENT.CHN = H.EVENT.CHN(sel);
                if isfield(H.EVENT,'TimeStamp')
                    H.EVENT.TimeStamp = H.EVENT.TimeStamp(sel);
                end;    
                % assigning new channel number 
                a = zeros(1,HDR.NS);
                for k = 1:length(InChanSelect),     % select channel specific elements
                    a(InChanSelect(k)) = k;     % assigning to new channel number. 
                end;
                ix = H.EVENT.CHN>0;
                H.EVENT.CHN(ix) = a(H.EVENT.CHN(ix));   % assigning new channel number
            end;    
            
        elseif strcmp(H.TYPE,'BDF');
            H.BDF.ANNONS = zeros(HDR.NRec*HDR.SPR,1);
            ix  = HDR.EVENT.TYP==hex2dec('7ffe');
            ix1 = find(~ix); 
            POS = HDR.EVENT.POS(ix1);
            TYP = HDR.EVENT.TYP(ix1);
            [s,ix2] = sort([POS; HDR.NRec*HDR.SPR+1]);
            for k = 1:length(s)-1,
                H.BDF.ANNONS(s(k):s(k+1)-1) = TYP(ix2(k));
            end;
            H.BDF.ANNONS(HDR.EVENT.POS(ix)) = H.BDF.ANNONS(HDR.EVENT.POS(ix)) | hex2dec('10000');
            
        elseif strcmp(H.TYPE,'BKR');
                    H.Classlabel = [];
                    H.TRIG = [];
                tmp=fullfile(H.FILE.Path,[H.FILE.Name,'.mat']);
                if ~exist(tmp,'file'),
                        tmp=fullfile(H.FILE.Path,[H.FILE.Name,'.MAT']);
                end
                x = [];
                if exist(tmp,'file'),
                        x = load('-mat',tmp);
                end;
 
                        if isfield(x,'header'),
                                H.MAT  = x.header;
                                if isfield(x.header,'Setup'), 
                                        if isfield(x.header.Setup,'Bits'), 
                                                H.Bits = x.header.Setup.Bits;
                                                [datatyp, limits, datatypes] = gdfdatatype(H.Bits+255);
                                                % THRESHOLD for Overflow detection
                                                if ~isfield(H,'THRESHOLD')
                                                    H.THRESHOLD = repmat(limits, H.NS, 1);
                                                end;         
                                        end;
                                end;
                                if isfield(x.header,'Result') && isfield(x.header.Result,'Classlabel'),
                                        H.Classlabel = x.header.Result.Classlabel;
                                end;
                                if isfield(x.header,'Paradigm')
                                        if isempty(H.Classlabel) && isfield(x.header.Paradigm,'Classlabel')
                                                H.Classlabel = x.header.Paradigm.Classlabel;
                                        end;
                                        H.BCI.Paradigm = x.header.Paradigm;
                                        if isfield(H.BCI.Paradigm,'TriggerOnset');
                                                H.TriggerOffset = H.BCI.Paradigm.TriggerOnset;
                                        elseif isfield(H.BCI.Paradigm,'TriggerTiming');
                                            %    H.BCI.Paradigm.TriggerTiming,
                                                H.TriggerOffset = H.BCI.Paradigm.TriggerTiming;
                                                fprintf(2,'Warning BKROPEN: Paradigm.TriggerOnset is unknown. Paradigm.TriggerTiming= %f ms is used instead\n',H.TriggerOffset);
                                        end;
                                end;
 
                                if isfield(x.header,'PhysioRec'), % R. Leeb's data 
                                       H.Label = cellstr(x.header.PhysioRec);
                                end;
                                if isfield(x.header,'BKRHeader'), % R. Scherer Data 
                                        if isfield(x.header.BKRHeader,'TO'),
                                                H.T0 = x.header.BKRHeader.TO;
                                        end;
                                        if isfield(x.header.BKRHeader,'Label'),
                                                H.Label = cellstr(x.header.BKRHeader.Label);
                                                ns = H.NS-length(H.Label);
                                                if ns == 1;
                                                        H.Label = strvcat(H.Label,'TRIGGER');
                                                elseif ns > 1;
                                                        H.Label = strvcat(H.Label,char(repmat('n.a.',ns,1)));
                                                end;
                                        end;
                                end;
                                if isfield(x.header,'Model'), % More 
                                        if isfield(x.header.Model,'AnalogInput'), 
                                                for k = 1:length(x.header.Model.AnalogInput),
                                                        H.Filter.HighPass(k) = x.header.Model.AnalogInput{k}{5};
                                                        H.Filter.LowPass(k)  = x.header.Model.AnalogInput{k}{6};
                                                        H.Filter.Notch(k)    = strcmpi(x.header.Model.AnalogInput{k}{7},'on');
 
                                                        H.MAT.Cal(k) = x.header.Model.AnalogInput{k}{3};
                                                end
                                        end;
                                end;
                                if ~isempty(strmatch('TRIGGER',H.Label))
                                        H.AS.TRIGCHAN = H.NS; %strmatch('TRIGGER',H.Label); 
                                end;
                        end;
                        if 1; ~isfield(H,'Classlabel');
                                tmp=fullfile(H.FILE.Path,[H.FILE.Name,'.par']);
                                if ~exist(tmp,'file'),
                                    tmp=fullfile(H.FILE.Path,[H.FILE.Name,'.PAR']);
                                end
                                if exist(tmp,'file'),
                                        H.Classlabel = load(tmp);
                                end;
                        end;
 
                        %%% Artifact Selection files 
                        tmp1=fullfile(H.FILE.Path,[H.FILE.Name,'.sel']);
                        if ~exist(tmp1,'file'),
                                tmp1=fullfile(H.FILE.Path,[H.FILE.Name,'.SEL']);
                        end
                        tmp2 = fullfile(H.FILE.Path,[H.FILE.Name,'_artifact.mat']);
                        SW   = (exist(tmp1,'file')>0) + 2*(exist(tmp2,'file')>0);
                        if SW == 0, 
                        elseif SW == 1,
                                if exist('OCTAVE_VERSION','builtin')
                                        H.ArtifactSelection = load('-ascii',tmp1);
                                else
                                        H.ArtifactSelection = load(tmp1);
                                end;
                        elseif SW == 2,
                                if exist('OCTAVE_VERSION','builtin')
                                        tmp = load('-mat',tmp2);
                                else
                                        tmp = load(tmp2);
                                end;
                                H.ArtifactSelection = tmp.artifact(:);
                        elseif SW == 3,
                                fprintf(H.FILE.stderr,'Warning BKROPEN: more than one ArtifactSelection files. File %s is used.\n',tmp1);
                                if exist('OCTAVE_VERSION')>5
                                        H.ArtifactSelection = load('-ascii',tmp1);
                                else
                                        H.ArtifactSelection = load(tmp1);
                                end;
                        end;
                        if isfield(H,'ArtifactSelection'),
                                if any(H.ArtifactSelection>1) || (length(H.ArtifactSelection)<length(H.Classlabel))
                                        sel = zeros(size(H.Classlabel));
                                        sel(H.ArtifactSelection) = 1;
                                        H.ArtifactSelection = sel(:);
                                end;
                                H.ArtifactSelection = H.ArtifactSelection(:);
                        end;
 
                        if isfield(H.AS,'TRIGCHAN') % & isempty(H.EVENT.POS)
                                if H.AS.TRIGCHAN <= H.NS, %size(H.data,2),
                                        H.THRESHOLD(H.AS.TRIGCHAN,1:2) = [-1-2^15,2^15]; % do not apply overflow detection for Trigger channel 
                                        data = mexSLOAD(H.FileName,H.AS.TRIGCHAN,'UCAL:ON','OVERFLOWDETECTION:OFF');
                                        TRIGon = gettrigger(data);
                                        %TRIGoff = gettrigger(-double(H.data(:,H.AS.TRIGCHAN)));
                                        if isfield(H,'TriggerOffset')
                                                TRIGon  = TRIGon - round(H.TriggerOffset/1000*H.SampleRate);
                                        %        TRIGoff = TRIGoff - round(H.TriggerOffset/1000*H.SampleRate);
                                        end;
                                end;
                                H.TRIG = TRIGon(:);
                                H.EVENT.POS = TRIGon(:); %[TRIGon(:); TRIGoff(:)]; 
                                H.EVENT.TYP = repmat(hex2dec('0300'),numel(TRIGon),1); %repmat(hex2dec('8300'),numel(TRIGoff),1)];
                        end;
                        if length(H.TRIG)~=length(H.Classlabel),
                                % hack to deal with BCI22 data
                                fprintf(2,'Warning BKROPEN: Number of triggers (%i) and number of Classlabels (%i) do not fit\n',length(H.TRIG),length(H.Classlabel));
                                H.TRIG = [];
                                H.Classlabel = [];
                                H.ArtifactSelection = [];
                        end;
            %% end of BKR
 
        elseif strcmp(H.TYPE,'BrainVision');
                    try 
                        H = bv2biosig_events(H); 
                    catch
                        %warning('bv2biosig_events not executed');
                    end; 
 
        elseif strcmp(H.TYPE,'EDF') && isfield(H.EVENT,'Desc');
            if length(H.EVENT.TYP)==length(H.EVENT.Desc)
                            [H.EVENT.CodeDesc, CodeIndex, H.EVENT.TYP] = unique(H.EVENT.Desc);
                    end;
 
            %% end of BKR
 
        elseif strcmp(H.TYPE,'BrainVision');
            HDR = sopen(fullfile(H.FILE.Path,[H.FILE.Name,'.vmrk']));
            HDR = sclose(HDR);
            H.EVENT = HDR.EVENT;
 
        elseif strcmp(H.TYPE,'PDP');
            signal = repmat(NaN,max(HDR.EVENT.POS),HDR.NS); 
                end;        
 
        if isempty(signal);
            signal = repmat(NaN,round(max(HDR.EVENT.POS)*HDR.SampleRate/HDR.EVENT.SampleRate),HDR.NS); 
            for k = 1:HDR.NS, 
                ix = find(HDR.EVENT.CHN==k);
                signal(round(HDR.EVENT.POS(ix)*HDR.SampleRate/HDR.EVENT.SampleRate),k)=HDR.EVENT.DUR(ix);
            end;
                end;        
 
            H.CHANTYP = repmat(' ',1,H.NS);
        for k=1:H.NS,
            if     ~isempty(strfind(lower(H.Label{k}),'eeg'))   H.CHANTYP(k) = 'E';
            elseif ~isempty(strfind(lower(H.Label{k}),'meg'))   H.CHANTYP(k) = 'E';
            elseif ~isempty(strfind(lower(H.Label{k}),'emg'))   H.CHANTYP(k) = 'M';
            elseif ~isempty(strfind(lower(H.Label{k}),'eog'))   H.CHANTYP(k) = 'O';
            elseif ~isempty(strfind(lower(H.Label{k}),'ecg'))   H.CHANTYP(k) = 'C';
            elseif ~isempty(strfind(lower(H.Label{k}),'air'))   H.CHANTYP(k) = 'R';
            elseif ~isempty(strfind(lower(H.Label{k}),'trig'))  H.CHANTYP(k) = 'T';
            end; 
        end;
 
    catch
        %fprintf(1,lasterr);
        fprintf(1, 'SLOAD: mexSLOAD(''%s'') failed - the slower M-function is used.\n', FILENAME);
    end;
else 
    global FLAG_HINT_mexSLOAD;
    if isempty(FLAG_HINT_mexSLOAD) 
        fprintf(1, 'Hint: the performance of SLOAD can be improved with mexSLOAD.mex which is part of biosig4c++.\n');
        FLAG_HINT_mexSLOAD = 1;     % turn off hint 
    end;
end;
 
if ~FlagLoaded,
    H = getfiletype(FILENAME);
 
    if isempty(H)
        fprintf(2,'Warning SLOAD: no file found\n');
        return;
    else    
        % FILENAME can be fn.name struct, or HDR struct. 
        FILENAME = H.FileName; 
    end;
    H.FLAG.UCAL = STATE.UCAL;
    H.FLAG.OVERFLOWDETECTION = STATE.OVERFLOWDETECTION;
    H.FLAG.EOG_CORRECTION = STATE.EOG_CORRECTION; 
    H.FLAG.OUTPUT = STATE.OUTPUT; 
 
 
if strncmp(H.TYPE,'IMAGE:',5)
    [H] = iopen(H);
    if H.FILE.OPEN,
        signal = iread(H);
        H.FILE.OPEN = 0; 
        fclose(H.FILE.FID);
    end;
    return;
end;
 
%%%%%%%%%%%%%%% --------- Load single file ------------%%%%%%%%%%%
H = sopen(H,'r',CHAN,MODE);
if 0, ~isnan(H.NS),
%------ ignore 'NaC'-channels
    NS = size(H.Calib,2);
    SelMx = speye(NS); 
    ch = 1:NS; %ch(strmatch('NaC',H.Label))=[];
    if length(ch)<NS,
        fprintf(2,'Warning SLOAD: Some NaC channels have been removed %s\n',H.FileName); 
    end; 
    SelMx = SelMx(:,ch); 
    H.Calib = H.Calib*SelMx; 
 
    % generate correction matrix 
    if H.FLAG.EOG_CORRECTION, 
        H.Calib = H.Calib*h.REGRESS.r0; 
        if ~isequal(H.Label(ch),h.Label)
            fprintf(2, 'Warning SLOAD (EOG correction): Channel labels in %s do not fit with arti* recording!\n',H.FileName); 
        end; 
    elseif STATE.EOG_CORRECTION, 
        fprintf(2, 'Warning: EOG correction not supported for this file (%s)!\n',H.FileName); 
    end; 
 
    %----- generate HDR.Calib -----
    if all(size(CHAN)>1) || any(floor(CHAN)~=CHAN) || any(CHAN<0) || (any(CHAN==0) && (numel(CHAN)>1));
            ReRefMx = CHAN; 
            CHAN = find(any(CHAN,2));
    elseif all(CHAN>0) && all(floor(CHAN)==CHAN), 
        if any(diff(CHAN)<=0),
        %   fprintf(HDR.FILE.FID,'Warning SOPEN: CHAN-argument not sorted - header information like Labels might not correspond to data.\n');
        end;    
            ReRefMx = sparse(CHAN,1:length(CHAN),1,H.NS,length(CHAN));
    else %if (CHAN==0)    
        ReRefMx = []; 
    end
    if ~isempty(ReRefMx),
        %[size(H.Calib),size(SelMx),size(h.REGRESS.r0),size(ReRefMx)] 
        ReRefMx = [ReRefMx(1:min(size(ReRefMx,1),size(H.Calib,2)),:); zeros(max(0,size(H.Calib,2)-size(ReRefMx,1)),size(ReRefMx,2))]; 
        H.Calib = H.Calib*ReRefMx; 
    end;    
 
    H.InChanSelect = find(any(H.Calib(2:end,:),2));
    H.Calib = H.Calib([1;1+H.InChanSelect(:)],:);
end
    
if 0,
        
elseif isfield(H,'data') && all(diag(H.Calib(2:end,1:end))==1) && ~H.FLAG.OVERFLOWDETECTION,
    % only a single copy of the data
    % important for large data sets, close to the available memory 
    signal = H.data; 
    H = rmfield(H,'data'); 
    if (CHAN>0)
        signal = signal(:,CHAN);
    end;     
 
 
elseif any(strmatch(H.TYPE,{'native','TFM_EXCEL_Beat_to_Beat','EEProbe-CNT','EEProbe-AVR'})); 
    [signal,H] = sread(H);  
        H = sclose(H);
 
elseif (H.FILE.OPEN > 0)
    [signal,H] = sread(H);  
        H = sclose(H);
 
elseif (H.FILE.OPEN > 0)
        signal = repmat(NaN,H.SPR*H.NRec,size(H.Calib,2));
    k1 = 0;
        while ~seof(H),
            [s,H] = sread(H,100);
            if 1,
            k2 = size(s,1);
            signal(k1+1:k1+k2,:)=s;
            k1 = k1+k2;
        else
                 %fprintf('%i/%i\n',stell(H),H.SPR*H.NRec);
                signal=[signal;s];
        end;
    end;      
        H = sclose(H);
 
        
elseif strncmp(H.TYPE,'EVENT',5)
        signal = H.EVENT;
        return; 
        
 
elseif strncmp(H.TYPE,'ELPOS',5)
        signal = H.ELEC.XYZ;
    return; 
 
elseif strcmp(H.TYPE,'DAQ')
    fprintf(1,'Loading a matlab DAQ data file - this can take a while.\n');
    tic;
        [signal, tmp, H.DAQ.T0, H.DAQ.events, DAQ.info] = daqread(H.FileName);
        fprintf(1,'Loading DAQ file finished after %.0f s.\n',toc);
        [H.SPR,H.NS] = size(signal);
        
        H.SampleRate = DAQ.info.ObjInfo.SampleRate;
        sz     = size(signal);
        if length(sz)==2, sz=[1,sz]; end;
        H.NRec = sz(1);
        H.Dur  = sz(2)/H.SampleRate;
        H.NS   = sz(3);
        H.FLAG.TRIGGERED = H.NRec>1;
        H.FLAG.UCAL = 1;
        
        H.PhysDim = {DAQ.info.ObjInfo.Channel.Units};
        H.DAQ   = DAQ.info.ObjInfo.Channel;
        
        H.Cal   = diff(cat(1,DAQ.info.ObjInfo.Channel.InputRange),[],2).*(2.^(-DAQ.info.HwInfo.Bits));
        H.Off   = cat(1,DAQ.info.ObjInfo.Channel.NativeOffset); 
        H.Calib = sparse([H.Off';eye(H.NS)]*diag(H.Cal));
        
        if CHAN<1,
                CHAN = 1:H.NS; 
        end;
        if ~H.FLAG.UCAL,
                Calib = H.Calib;    % Octave can not index sparse matrices within a struct
                signal = [ones(size(signal,1),1),signal]*Calib(:,CHAN);
        end;
        
        
elseif 0, strcmp(H.TYPE,'BIFF'),
    try, 
                [H.TFM.S,H.TFM.E] = xlsread(H.FileName,'Beat-To-Beat');
            if size(H.TFM.S,1)+1==size(H.TFM.E,1),
                        H.TFM.S = [repmat(NaN,1,size(H.TFM.S,2));H.TFM.S];
                end;
                H.TYPE = 'TFM_EXCEL_Beat_to_Beat'; 
    catch
        try
                    [H.TFM.S,H.TFM.E] = xlsread(H.FileName,'Beat-to-Beat');
                    H.TYPE = 'TFM_EXCEL_Beat_to_Beat'; 
            catch 
            0; 
            end;
    end;    
 
    if strcmp(H.TYPE, 'TFM_EXCEL_Beat_to_Beat');
                if ~isempty(strfind(H.TFM.E{3,1},'---'))
                        H.TFM.S(3,:) = [];    
                        H.TFM.E(3,:) = [];    
                end;
                
                H.Label   = strvcat(H.TFM.E{4,:});
                H.PhysDim = strvcat(H.TFM.E{5,:});
           
                H.TFM.S = H.TFM.S(6:end,:);
                H.TFM.E = H.TFM.E(6:end,:);
                
                ix = find(isnan(H.TFM.S(:,2)) & ~isnan(H.TFM.S(:,1)));
                H.EVENT.Desc = H.TFM.E(ix,2);
                H.EVENT.POS  = ix;
                
                if any(CHAN),
            H.TFM.S = H.TFM.S(:,CHAN);
            H.TFM.E = H.TFM.E(:,CHAN);
        end;
        [H.SPR,H.NS] = size(H.TFM.S);
        H.NRec = 1; 
        H.THRESHOLD  = repmat([0,NaN],H.NS,1);
 
        signal  = H.TFM.S;
        signal(signal==0) = NaN;
        end;
 
 
elseif strcmp(H.TYPE,'MatrixMarket'),
        H.FILE.FID = fopen(H.FileName,'rt','ieee-le');
 
        line = fgetl(H.FILE.FID);
    
    H.FLAG.Coordinate = ~isempty(strfind(line,'coordinate'));
    H.FLAG.Array      = ~isempty(strfind(line,'array'));
 
    H.FLAG.Complex = ~isempty(strfind(line,'complex'));
    H.FLAG.Real = ~isempty(strfind(line,'real'));
    H.FLAG.Integer = ~isempty(strfind(line,'integer'));
    H.FLAG.Pattern = ~isempty(strfind(line,'pattern'));
    
    H.FLAG.General = ~isempty(strfind(line,'general'));
    H.FLAG.Symmetric = ~isempty(strfind(line,' symmetric'));
    H.FLAG.SkewSymmetric = ~isempty(strfind(line,'skew-symmetric'));
    H.FLAG.Hermitian = ~isempty(strfind(lower(line),'hermitian'));
 
    while strncmp(line,'%',1)
            line = fgetl(H.FILE.FID);
    end;
 
    [tmp,status] = str2double(line);
    if any(status)
        fprintf(H.FILE.stderr,'SLOAD (MM): invalid size %s\n',line);
    else
        H.MATRIX.Size = tmp;
    end;    
 
    if length(H.MATRIX.Size)==3,
        H.Length = tmp(3);
        signal = sparse([],[],[],tmp(1),tmp(2),tmp(3));
        for k = 1:H.Length,
                line = fgetl(H.FILE.FID);
            [tmp,status] = str2double(line);
            if any(status)
                fprintf(H.FILE.stderr,'SLOAD (MM): invalid size %s\n',line);
            elseif length(tmp)==4,  
                    val = tmp(3) + i*tmp(4);
            elseif length(tmp)==3,  
                val = tmp(3);
            elseif length(tmp)==2,  
                val = 1;
            else
                fprintf(H.FILE.stderr,'SLOAD (MM): invalid size %s\n',line);
            end;
 
            if H.FLAG.General,
                signal(tmp(1),tmp(2)) = val;
            elseif H.FLAG.Symmetric,
                signal(tmp(1),tmp(2)) = val;
                signal(tmp(2),tmp(1)) = val;
            elseif H.FLAG.SkewSymmetric,
                signal(tmp(1),tmp(2)) = val;
                signal(tmp(2),tmp(1)) =-val;
            elseif H.FLAG.Hermitian,
                signal(tmp(1),tmp(2)) = val;
                signal(tmp(2),tmp(1)) = conj(val);
            else    
                fprintf(H.FILE.stderr,'SLOAD (MM): invalid size %s\n',line);
            end;    
        end;
                    
    elseif length(H.MATRIX.Size)==2
        H.Length = prod(tmp);
        signal = zeros(H.MATRIX.Size);
        if H.FLAG.General==1,
            [IX,IY]=find(ones(H.MATRIX.Size));
        else
            [IX,IY]=find(cumsum(eye(H.MATRIX.Size)));
        end;
                
        for k = 1:H.Length,
                line = fgetl(H.FILE.FID);
            [tmp,status] = str2double(line);
            if any(status)
                error('SLOAD (MM)');
            elseif length(tmp)==2,  
                val=tmp(1) + i*tmp(2);
            elseif length(tmp)==1,  
                val=tmp(1);
            else
                fprintf(H.FILE.stderr,'SLOAD (MM): invalid size %s\n',line);
            end;
 
            signal(IX(k),IY(k)) = val;
            if H.FLAG.Symmetric,
                signal(IY(k),IX(k)) = val;
            elseif H.FLAG.SkewSymmetric,
                signal(IY(k),IX(k)) =-val;
            elseif H.FLAG.Hermitian,
                signal(IY(k),IX(k)) = conj(val);
            else    
                fprintf(H.FILE.stderr,'SLOAD (MM): invalid size %s\n',line);
            end;    
        end;
        end;
        fclose(H.FILE.FID);
 
        
elseif strcmp(H.TYPE,'OFF'),
            H.FILE.FID = fopen(H.FileName,'rt','ieee-le');
        
                line1 = fgetl(H.FILE.FID);
                line2 = fgetl(H.FILE.FID);
        while ~feof(H.FILE.FID) && (line2(1)=='#'),
                    line2 = fgetl(H.FILE.FID);
        end;
        [tmp,status] = str2double(line2);
        if status || (size(tmp,2)~=3), 
            fclose(H.FILE.FID);
            error('SOPEN (OFF)');
        else
            H.VertexCount = tmp(1);
            H.FaceCount = tmp(2);
            H.EdgeCount = tmp(3);
        end 
        
        H.Vertex = repmat(NaN,H.VertexCount,3);
        for k = 1:H.VertexCount,
            line = '';
            while isempty(line) || strncmp(line,'#',1)
                line = fgetl(H.FILE.FID);
            end;
            len = min(length(line),min(find(line=='#')));
            tmp = str2double(line(1:len));
            H.Vertex(k,:) = tmp(1:H.ND);
        end;    
        
%       H.Face = repmat(NaN,H.FaceCount,3);
        for k = 1:H.FaceCount,
            line = '';
            while isempty(line) || strncmp(line,'#',1)
                line = fgetl(H.FILE.FID);
            end;
            len = min(length(line),min(find(line=='#')));
            tmp = str2double(line(1:len));
            H.Ngon(k) = tmp(1);
            H.Face{k} = tmp(2:tmp(1)+1) + 1;
        end;    
        if all(H.Ngon(1)==H.Ngon),
            H.Face = cat(1,H.Face{:});
        end;    
                fclose(H.FILE.FID);
        
        
elseif strcmp(H.TYPE,'POLY'),
        H.FILE.FID = fopen(H.FileName,'rt','ieee-le');
 
    K = 0;
    while ~feof(H.FILE.FID)
            line = fgetl(H.FILE.FID);
        if isempty(line),
        elseif line(1)=='#',
        else
            K = K + 1;
            
        end;
        end;
        fclose(H.FILE.FID);
        
        
elseif strcmp(H.TYPE,'SMF'),
            H.FILE.FID = fopen(H.FileName,'rt','ieee-le');
        
        VertexCount = 0;
        FaceCount = 0;
        PalLen = 0; 
        K = 1;
        while ~feof(H.FILE.FID)
                    line = fgetl(H.FILE.FID);
            if isempty(line)
            elseif line(1)=='#';
 
            elseif line(1)=='v';
                [tmp,status] = str2double(line(3:end));
                if ~any(status)
                    VertexCount = VertexCount + 1 ;
                    H.Vertex(VertexCount,:) = tmp;
                else
                    fprintf(H.FILE.stderr,'Warning SLOAD: could not read line %i in file %s\n',K,H.FileName);   
                end;    
 
            elseif line(1)=='f';
                [tmp,status] = str2double(line(3:end));
                if ~any(status)
                    FaceCount  = FaceCount + 1; 
                    H.Ngon(FaceCount) = length(tmp);
                    H.Face{FaceCount} = tmp;
                else
                    fprintf(H.FILE.stderr,'Warning SLOAD: could not read line %i in file %s\n',K,H.FileName);   
                end;    
 
            elseif line(1)=='n';
                [tmp,status] = str2double(line(3:end));
                if ~any(status)
                    H.NormalVector = tmp;
                else
                    fprintf(H.FILE.stderr,'Warning SLOAD: could not read line %i in file %s\n',K,H.FileName);   
                end;    
 
            elseif line(1)=='c';
                [tmp,status] = str2double(line(3:end));
                if ~any(status)
                    PalLen = PalLen +1; 
                    H.Palette(PalLen,:)= tmp;
                else
                    fprintf(H.FILE.stderr,'Warning SLOAD: could not read line %i in file %s\n',K,H.FileName);   
                end;    
            else
 
            end;
            K = K+1;
        end;
        fclose(H.FILE.FID);
        if all(H.Ngon(1)==H.Ngon),
            H.Face = cat(1,H.Face{:});
        end;    
    
        
elseif strcmp(H.TYPE,'TVF 1.1A'),
        H.FILE.FID = fopen(H.FileName,'rt');
 
    tmp = fgetl(H.FILE.FID);
    H.TVF.Name = fgetl(H.FILE.FID);
    H.TVF.Desc = fgetl(H.FILE.FID);
    [tmp,status] = str2double(fgetl(H.FILE.FID));
    H.TVF.NTR = tmp(1);
    H.TVF.NV  = tmp(2);
    H.FLAG.CULL = ~~tmp(3);
    [tmp,status] = str2double(fgetl(H.FILE.FID));
    H.TVF.NTRC  = tmp(1);
    H.TVF.NTRM  = tmp(2);
    [tmp,status] = str2double(fgetl(H.FILE.FID));
    H.TVF.NVC   = tmp(1);
    H.TVF.NVN   = tmp(2);
    [tmp,status] = str2double(fgetl(H.FILE.FID));
    H.TVF.GlobalColor = tmp;
    [tmp,status] = str2double(fgetl(H.FILE.FID));
    H.TVF.GlobalMtrlProps = tmp;
    
    H.TVF.Triangles = repmat(NaN,[H.TVF.NTR,3]);
    for k = 1:H.TVF.NTR,
        [tmp, status] = str2double(fgetl(H.FILE.FID));
        H.TVF.Triangles(k,:) = tmp;
    end;    
    H.TVF.TrColorSets = reshape(NaN,[H.TVF.NTRC,H.TVF.NTR]);
    for k = 1:H.TVF.NTRC,
        [tmp, status] = str2double(fgetl(H.FILE.FID));
        H.TVF.TrColorSets(k,:) = tmp;
    end;
    H.TVF.TrMtrlSets = repmat(NaN, [H.TVF.NTRM,H.TVF.NTR]);
    for k = 1:H.TVF.NTRM,
        [tmp, status] = str2double(fgetl(H.FILE.FID));
        H.TVF.TrMtrlSets(k,:) = tmp;
    end;
 
    H.TVF.Vertices   = repmat(NaN, [H.TVF.NV,3]);
    for k = 1:H.TVF.NV,
        [tmp, status] = str2double(fgetl(H.FILE.FID));
        H.TVF.Vertices(k,:) = tmp;
    end;
    H.TVF.VColorSets = repmat(NaN, [H.TVF.NVC,H.TVF.NV]);
    for k = 1:H.TVF.NVC,
        [tmp, status] = str2double(fgetl(H.FILE.FID));
        H.TVF.VColorSets(k,:) = tmp;
    end;
    H.TVF.VNrmlSets  = repmat(NaN, [H.TVF.NVN,3]);
    for k = 1:H.TVF.NVN,
        [tmp, status] = str2double(fgetl(H.FILE.FID));
        H.TVF.VNrmlSets(k,:) = tmp;
    end;
 
    fclose(H.FILE.FID); 
 
        
elseif strcmp(H.TYPE,'TVF 1.1B'),
        H.FILE.FID = fopen(H.FileName,'rb',H.Endianity);
 
    tmp = fread(H.FILE.FID,12,'uchar');
    H.TVF.Name = char(fread(H.FILE.FID,32,'uchar'));
    H.TVF.Desc = char(fread(H.FILE.FID,80,'uchar'));
    tmp = fread(H.FILE.FID,7,'uint32');
    H.TVF.NTR = tmp(1);
    H.TVF.NV  = tmp(2);
    H.FLAG.CULL = ~~tmp(3);
    H.TVF.NTRC  = tmp(4);
    H.TVF.NTRM  = tmp(5);
    H.TVF.NVC   = tmp(6);
    H.TVF.NVN   = tmp(7);
    tmp = fread(H.FILE.FID,[4,2],'float32');
    H.TVF.GlobalColor = tmp(:,1)';
    H.TVF.GlobalMtrlProps = tmp(:,2)';
 
    H.TVF.Triangles   = fread(H.FILE.FID, [3,H.TVF.NTR],'uint32')';
    H.TVF.TrColorSets = fread(H.FILE.FID, [H.TVF.NTR,H.TVF.NTRC],'uint32')';
    H.TVF.TrMtrlSets  = fread(H.FILE.FID, [H.TVF.NTR,H.TVF.NTRM],'uint32')';
 
    H.TVF.Vertices   = fread(H.FILE.FID, [3,H.TVF.NV],'float32')';
    H.TVF.VColorSets = fread(H.FILE.FID, [H.TVF.NV,H.TVF.NVC],'uint32')';
    H.TVF.VNrmlSets  = fread(H.FILE.FID, [3,H.TVF.NVN],'float32')';
 
    fclose(H.FILE.FID); 
 
        
elseif strcmp(H.TYPE,'VTK'),
                H.FILE.FID = fopen(H.FileName,'rt','ieee-le');
                
                H.VTK.version = fgetl(H.FILE.FID);
                H.VTK.Title   = fgetl(H.FILE.FID);
                H.VTK.type    = fgetl(H.FILE.FID);
        
        while ~feof(H.FILE.FID),
            tline = fgetl(H.FILE.FID);
            
            if 0, 
            elseif strncmp(tline,'CELLS',5);
            elseif strncmp(tline,'CELL_DATA',9);
            elseif strncmp(tline,'DATASET',7);
                        H.VTK.DATASET = tline(8:end);
            elseif strncmp(tline,'LINES',6);
            elseif strncmp(tline,'POINTS',6);
            elseif strncmp(tline,'POINT_DATA',10);
            elseif strncmp(tline,'POLYGON',7);
            elseif strncmp(tline,'SCALARS',7);
                [t1,r]=strtok(tline);
                [dataName,r]=strtok(r);
                [dataType,r]=strtok(r);
                [numComp ,r]=strtok(r);
                if isempty(numComp), numComp=1;
                else numComp = str2double(numComp); end;
                tline = fgetl(fid);
                if strcmp(tline,'LOOKUP_TABLE');
                        [t1,r]=strtok(tline);
                        [tableName,r]=strtok(tline);
                else    
                    tline = fgetl(fid);
                end;
            end;
 
        
        end;
                fclose(H.FILE.FID);
 
                fprintf(H.FILE.stderr,'Warning SOPEN: VTK-format not supported, yet.\n');
    
        
elseif strcmp(H.TYPE,'unknown')
    [p,f,e]=fileparts(H.FileName);
        TYPE = upper(e);
        if 0, 
        elseif strcmp(TYPE,'.RAW')
                loadraw;
        elseif strcmp(TYPE,'.RDT')
                [signal] = loadrdt(FILENAME,CHAN);
                fs = 128;
        elseif strcmp(TYPE,'.XLS')
                loadxls;
        elseif strcmp(TYPE,'.DA_')
                fprintf('Warning SLOAD: Format DA# in testing state and is not supported\n');
                loadda_;
        elseif strcmp(TYPE,'.RG64')
                [signal,H.SampleRate,H.Label,H.PhysDim,H.NS]=loadrg64(FILENAME,CHAN);
                %loadrg64;
        else
                fprintf('Error SLOAD: Unknown Data Format\n');
                signal = [];
                return;
        end;
end;
 
end; 
 
%%%%%%%%%% Post-Processing %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if strcmp(H.TYPE,'CNT');    
        f = fullfile(H.FILE.Path, [H.FILE.Name,'.txt']); 
        if exist(f,'file'),
                fid = fopen(f,'r');
        tmp = fread(fid,inf,'uint8');
        fclose(fid);
        [tmp,v] = str2double(char(tmp'));
        if ~any(v), 
                    H.Classlabel=tmp(:);                        
            end;
        end
        f = fullfile(H.FILE.Path, [H.FILE.Name,'.par']); 
        if exist(f,'file'),
                fid = fopen(f,'r');
        tmp = fread(fid,inf,'uint8');
        fclose(fid);
        [tmp,v] = str2double(char(tmp'));
        if ~any(v), 
                    H.Classlabel=tmp(:);                        
            end;
        end
        f = fullfile(H.FILE.Path, [H.FILE.Name,'.mat']);
        if exist(f,'file'),
                try,
                        tmp = load(f);
                catch
                        tmp = []; 
                end
                if isfield(tmp,'classlabel') && ~isfield(H,'Classlabel')
                        H.Classlabel=tmp.classlabel(:);                        
                elseif isfield(tmp,'classlabel') && isfield(tmp,'header') && isfield(tmp.header,'iniFile') && strcmp(tmp.header.iniFile,'oom.ini'), %%% for OOM study only. 
                        H.Classlabel=tmp.classlabel(:);                        
                end;
        end;
        f = fullfile(H.FILE.Path, [H.FILE.Name,'c.mat']);
        if exist(f,'file'),
                try,
                        tmp = load(f);
                catch
                        tmp = []; 
                end
                if isfield(tmp,'classlabel') && ~isfield(H,'Classlabel')
                        H.Classlabel=tmp.classlabel(:);                        
                end;
        end;
        f = fullfile(H.FILE.Path, [H.FILE.Name,'_classlabel.mat']);
        if exist(f,'file'),
                try,
                        tmp = load(f);
                catch
                        tmp = []; 
                end
                if isfield(tmp,'Classlabel') && (size(tmp.Classlabel,2)==4)
                        [x,H.Classlabel] = max(tmp.Classlabel,[],2);                        
                end;
                if isfield(tmp,'classlabel') && (size(tmp.classlabel,2)==4)
                        [x,H.Classlabel] = max(tmp.classlabel,[],2);                        
                end;
        end;
        
        f=fullfile(H.FILE.Path,[H.FILE.Name,'.sel']);
        if ~exist(f,'file'),
                f=fullfile(H.FILE.Path,[H.FILE.Name,'.SEL']);
        end
        if exist(f,'file'),
                fid = fopen(f,'r');
        tmp = fread(fid,inf,'uint8');
        fclose(fid);
        [tmp,v] = str2double(char(tmp'));
        if any(isnan(tmp)) || any(tmp~=ceil(tmp)) || any(tmp<0) || (any(tmp==0) && any(tmp>1))
                        fprintf(2,'Warning SLOAD(CNT): corrupted SEL-file %s\n',f);
                else
                        if ~isfield(H,'Classlabel'), 
                                H.Classlabel = H.EVENT.TYP;
                        end;        
                        n = length(H.Classlabel);
                        
                        H.ArtifactSelection = zeros(n,1);
                        if all((tmp==0) | (tmp==1)) && (length(tmp)>1) && (sum(diff(sort(tmp))~=0) ~= length(tmp)-1)
                                H.ArtifactSelection = logical(tmp);         
                        else
                                H.ArtifactSelection(tmp) = 1;         
                        end;
                end;
        end;
end;
 
if ~isempty(strfind(upper(MODE),'TSD'));
        f = fullfile(H.FILE.Path, [H.FILE.Name,'.tsd']);
        if ~exist(f,'file'),
                        fprintf(2,'Warning SLOAD-TSD: file %s.tsd found\n',H.FILE(1).Name,H.FILE(1).Name);
        else
                fid = fopen(f,'rb');
                tsd = fread(fid,inf,'float');
                fclose(fid);
                nc = size(signal,1)\size(tsd,1);
                if (nc == round(nc)),
                        signal = [signal, reshape(tsd,nc,size(tsd,1)/nc)'];
                else
                        fprintf(2,'Warning SLOAD: size of %s.tsd does not fit to size of %s.bkr\n',H.FILE(1).Name,H.FILE(1).Name);
                end;
        end;
end;
 
 
if strcmp(H.TYPE,'GDF')
        if isfield(H.EVENT,'TYP')
                if isempty(H.EVENT.TYP)
                        %%%%% if possible, load Reinhold's configuration files
                        f = fullfile(H.FILE.Path, [H.FILE.Name,'.mat']);
                        if exist(f,'file'),
                                try
                                        x = load(f,'header');
                                catch;
                                        x=[];
                                end;
                                if isfield(x,'header'),
                                        if isfield(x.header,'Paradigm')
                                                fprintf(2,'Warning SLOAD (GDF): Obsolete feature is used (header information loaded from MAT-file %s).\n',H.FileName);
                                                fprintf(2,'This feature will be removed in future. Contact <a.schloegl@ieee.org> if you need this.\n ');
                                                H.BCI.Paradigm = x.header.Paradigm;
                                                if isfield(H.BCI.Paradigm,'TriggerTiming');
                                                        H.TriggerOffset = H.BCI.Paradigm.TriggerTiming;
                                                elseif isfield(H.BCI.Paradigm,'TriggerOnset');
                                                        H.TriggerOffset = H.BCI.Paradigm.TriggerOnset;
                                                end;
 
                                                if isfield(H,'Classlabel') && isempty(H.Classlabel),
                                                        H.Classlabel = x.header.Paradigm.Classlabel;
                                                end;
                                        end;
                                end;
                        end;
                end;
        end;
 
        fid=fopen(fullfile(H.FILE.Path,[H.FILE.Name,'.sel']),'r');
        if fid<0,
            fid=fopen(fullfile(H.FILE.Path,[H.FILE.Name,'.SEL']),'r');
        end
        if fid>0,
                [tmp,c] = fread(fid,[1,inf],'uint8');
                fclose(fid);
                [tmp,v,sa] = str2double(tmp);
                if isempty(sa) || isempty(sa{1})
                        H.ArtifactSelection = repmat(0,length(H.TRIG),1);
                elseif any(isnan(tmp)) || any(tmp~=ceil(tmp)) || any(tmp<0) || (any(tmp==0) && any(tmp>1))
                        fprintf(2,'Warning SLOAD(GDF): corrupted SEL-file %s\n',fullfile(H.FILE.Path,[H.FILE.Name,'.sel']));
        else                        
                        H.ArtifactSelection = zeros(length(H.TRIG),1);
                        if all((tmp==0) | (tmp==1)) && (length(tmp)>1) && (sum(diff(sort(tmp))~=0) ~= length(tmp)-1)
                                H.ArtifactSelection = logical(tmp);
                        else
                                H.ArtifactSelection(tmp) = 1;
                        end;
                end;
        end;
end;
 
 
%%%%% introduce NaNs between segments %%%%%%%
if STATE.FLAG_NUM_NAN,
    %% NaN's between segments/sweeps are introduced, but only when argument NUMBER_OF_NAN_IN_BREAK is explicitely defined in the list of input arguments
    winlen = STATE.NUMBER_OF_NAN_IN_BREAK;
    ix = [1; H.EVENT.POS(H.EVENT.TYP==hex2dec('7ffe')); size(signal,1)+1];
    n  = length(ix) - 1;
    signal  = [signal; repmat(NaN, (n-1)*winlen, size(signal, 2))]; 
    for k = n:-1:1,
        H.EVENT.POS(H.EVENT.POS >= ix(k+1)) = H.EVENT.POS(H.EVENT.POS >= ix(k+1) ) + winlen;
        signal([ix(k) : ix(k+1) - 1] + (k-1) * winlen, :) = signal(ix(k) : ix(k+1) - 1,:);
        if k>1,
            signal(ix(k) + (k-1) * winlen + [-winlen:-1], :) = NaN;
        end; 
    end;
    if H.SPR * H.NRec ~= size(signal,1),
        H.SPR = H.SPR + winlen * (n - 1);
    end;
 
    if ~STATE.FLAG_NUM_NAN && (H.SPR * H.NRec ~= size(signal,1))
        warning('Undefined state: HDR.SPR*HDR.NRec does not match size of data - if you see this message contact <alois.schloegl@gmail.com> and provide the details %s', HDR.FileName); 
    end; 
    % [HDR.SPR, HDR.NRec, size(signal)]
end;
 
 
%%%%% Resampling %%%%%
if ~isnan(Fs) && (H.SampleRate~=Fs);
        tmp = ~mod(H.SampleRate,Fs) || ~mod(Fs,H.SampleRate);
        tmp2= ~mod(H.SampleRate,Fs*2.56);
        if tmp,
                signal = rs(signal,H.SampleRate,Fs);
                H.EVENT.POS = H.EVENT.POS/H.SampleRate*Fs;
                if isfield(H.EVENT,'DUR');
                        H.EVENT.DUR = H.EVENT.DUR/H.SampleRate*Fs;
                end;
                if isfield(H,'TRIG');
                        H.TRIG = H.TRIG/H.SampleRate*Fs;
                end;
                H.SampleRate = Fs;
        elseif tmp2,
                x = load('resample_matrix.mat');
                signal = rs(signal,x.T256100);
                if H.SampleRate*100~=Fs*256,
                        signal = rs(signal,H.SampleRate/(Fs*2.56),1);
                end;
                H.EVENT.POS = H.EVENT.POS/H.SampleRate*Fs;
                if isfield(H.EVENT,'DUR');
                        H.EVENT.DUR = H.EVENT.DUR/H.SampleRate*Fs;
                end;
                if isfield(H,'TRIG');
                        H.TRIG = H.TRIG/H.SampleRate*Fs;
                end;
                H.SampleRate = Fs;
        else 
                fprintf(2,'Warning SLOAD: resampling %f Hz to %f Hz not implemented.\n',H.SampleRate,Fs);
        end;                
end;
 
return; 
ratiomissing = mean(isnan(signal)); 
if any(ratiomissing>.1)
    fprintf(2,'Warning SLOAD: ratio of missing samples exceeds 10%% in file %s.\n',H.FileName);
    ix = find(ratiomissing); 
    fprintf(1,'#%3i:  %4.1f%%\n',[ix;ratiomissing(ix)*100])
end;