save2bkr.m

matlab数字信号处理工具箱

        %fwrite(HDR.FILE.FID,data','int16');  	% WRITE BKR FILE
        HDR = sclose(HDR);            % CLOSE BKR FILE
        
        % save Classlabels
        if isfield(HDR,'Classlabel'),
                fid = fopen([HDR.FileName(1:length(HDR.FileName)-4) '.par'],'wt');
                fprintf(fid, '%i\n', HDR.Classlabel);
                fclose(fid);
        end;
        if isfield(HDR,'ArtifactSelection'),
                fid = fopen([HDR.FileName(1:length(HDR.FileName)-4) '.sel'],'w');
                fprintf(fid, '%i\r\n', HDR.ArtifactSelection);
                fclose(fid);
        end;
        
        % final test 
        try
                HDR = sopen(HDR.FileName,'r');
                HDR = sclose(HDR);
        catch
                fprintf(2,'Error SAVE2BKR: saving file %s failed\n',HDR.FileName);
        end;
        return;
end;

for k=1:length(infile);
        filename = fullfile(inpath,infile(k).name);
        [pf,fn,ext] = fileparts(filename);
        
        % load eeg data 
        [y,HDR] = sload(filename);
        
        % load classlabels if the exist
        tmp = fullfile(HDR.FILE.Path,[HDR.FILE.Name,'.mat']);
        if exist(tmp),
                tmp=load(tmp);
                if isfield(tmp,'classlabel');
                        HDR.Classlabel = tmp.classlabel;
                end;
        end;
        
        if isempty(y), 
                fprintf(2,'Error SAVE2BKR: file %s not found\n',filename);
                return; 
        end; 
        
        if ~isfield(HDR,'NS'),
                warning(['number of channels undefined in ',filename]);
                HDR.NS = size(y,2);
        end;
        if ~HDR.FLAG.TRIGGERED,
                HDR.NRec = 1; 
                HDR.SPR = size(y,1);
        end;
        if ~isfield(HDR,'NRec'),
                HDR.NRec = 1;
        end;
        if ~isfield(HDR,'SPR'),
                HDR.SPR = size(y,1)/HDR.NRec;
        elseif length(HDR.SPR)>1,       % use just one sampling rate 
                HDR.SPR = HDR.AS.MAXSPR;
                HDR.SampleRate = HDR.AS.MAXSPR/HDR.Dur;
                FLAG_PHYSMAX = 1; 
                PHYSMAX = max(abs(y(:)));
                HDR.DigMax  = 2^15-1;
        end;

        if FLAG_REGRESS_EOG,
                fprintf(1,'\tREGRESS_EOG \n');
                [R,y] = regress_eog(y,chansel_dt3,chansel_dt4);
        end;

        if FLAG_REMOVE_DC,
                fprintf(1,'\tREMOVE_DC \n');
                y = y - repmat(mean(y,1),size(y,1),1);
        end;
        if FLAG_DETREND,
                B = -ones(1,HDR.SampleRate)/HDR.SampleRate;
                B(HDR.SampleRate/2) = B(HDR.SampleRate/2)+1;
                HDR.Filter.B = B;
                HDR.Filter.A = 1;
                %HDR.Filter.B=B;%conv(-B, HDR.Filter.B);
                Delay = length(B)/2;        
                HDR.FLAG.FILT = 1;
                HDR.Filter.HighPass = .5;
                
                for k = chansel_dt,
                        tmp = filter(B,1,[y(:,k);zeros(length(B),1)]);
                        y(:,k) = tmp(Delay+1:size(y,1)+Delay);
                end;                
        end;
        
        if FLAG_removeDrift,
                B = .5*(1 - cos(2*pi*(1:4*HDR.SampleRate+1)'/(4*HDR.SampleRate+2))); 
                B = -B/sum(B);
                B(2*HDR.SampleRate) = B(HDR.SampleRate)+1;
                
                B = -ones(1,HDR.SampleRate)/HDR.SampleRate;
                B(HDR.SampleRate/2) = B(HDR.SampleRate/2)+1;
                %B(1) = B(1)+1;
                
                HDR.Filter.B = B;
                HDR.Filter.A = 1;
                %HDR.Filter.B=B;%conv(-B, HDR.Filter.B);
                Delay = (length(B)-1)/2;        
                HDR.FLAG.FILT = 1;
                HDR.Filter.HighPass = .5;
                
                for k = chansel_dt2,
                        y(:,k) = filtfilt(B,1,y(:,k));
                        %y(:,k) = tmp(Delay+1:size(y,1)+Delay);
                end;                
        end;
        
        if chansel == 0;
                chansel=1:HDR.NS;
        end;
        
        % add event channel 
        if isfield(HDR,'EVENT')
                if HDR.EVENT.N <= 0,
                elseif 0,
                        % TypeList = unique(HDR.EVENT.TYP); but ignores NaN's
                        [sY ,idx] = sort(HDR.EVENT.TYP(:));
                        TypeList  = sY([1;find(diff(sY,1)>0)+1]);
                        
                        event = zeros(size(y,1),length(TypeList));
                        for k2 = 1:length(TypeList),
                                tmp = (HDR.EVENT.TYP==TypeList(k2));
                                event(HDR.EVENT.POS(tmp),k2) = HDR.EVENT.TYP(tmp);        
                        end;
                        HDR.NS = HDR.NS + size(event,2);
                        y = [y, event];
                elseif all(HDR.EVENT.TYP < 256),  % only NeuroScan Events are converted into separate channels
                        K = 0; event = [];
                        for k1 = 0:7,
                                tmp = bitand(HDR.EVENT.TYP,2^k1);
                                if any(tmp),
                                        K = K+1;
                                        event(size(y,1),K) = 0;        
                                        event(HDR.EVENT.POS(tmp>0),K) = 1;        
                                end;				
                        end;    
                        if any(sum(event,2)>1),
                                fprintf(2,'Warning SAVE2BKR: simulateneous events occur. \n');
                        end;	
                        HDR.NS = HDR.NS + size(event,2);
                        y = [y, event];
                end;
        end;
        
        % re-scale data to account for the scaling factor in the header
        HDR.DigMax = 2^15-1;
        if FLAG_PHYSMAX,
                HDR.PhysMax = PHYSMAX;
        else
                tmp = y(:,chansel);
                HDR.PhysMax = max(abs(tmp(:))); %gives max of the whole matrix
        end;
        for k = 1:HDR.NS,
                if any(k==chansel),
                        y(:,k) = y(:,k)*HDR.DigMax/HDR.PhysMax; % keep correct scaling factor 
                else
                        mm = max(abs(y(:,k)));
                        y(:,k) = y(:,k)*HDR.DigMax/mm;          % scale to maximum resolution
                end;
        end;
        HDR.FLAG.UCAL = 1;              % data is de-calibrated, no rescaling within SWRITE 
        
        tmp = round(HDR.PhysMax);
        fprintf(1,'Rounding of PhysMax yields %f%% error.\n',abs((HDR.PhysMax-tmp)/tmp)*100);
        HDR.PhysMax = tmp;
        HDR.TYPE = 'BKR';
        HDR.FLAG.REFERENCE = ' ';
        HDR.FLAG.TRIGGERED = (HDR.NRec>1);
        
        if isempty(outfile), 	% default destination directory  
                ix = max(find(filename=='.'));
                %HDR.FileName = [filename(1:ix-1),'.bkr'];  % destination directory is same as source directory 
                HDR.FileName  = [HDR.FILE.Name,'.bkr'];     % destination directory is current working directory 
        elseif isdir(outfile),	% output file
                HDR.FILE.Path = outfile;            
                HDR.FileName  = fullfile(outfile,[HDR.FILE.Name,'.bkr']);
        else
                [HDR.FILE.Path,HDR.FILE.Name,Ext] = fileparts(outfile);
                HDR.FileName = fullfile(HDR.FILE.Path,[HDR.FILE.Name,Ext]);
        end;
        %HDR = eegchkhdr(HDR);
        
        HDR = sopen(HDR,'w');
        if HDR.FILE.FID < 0,
                fprintf(1,'Error SAVE2BKR: couldnot open file %s.\n',HDR.FileName);
                return;
        end;
        % writes data
        HDR = swrite(HDR,y(:,1:HDR.NS));  	% WRITE BKR FILE
        %count = fwrite(HDR.FILE.FID,y','short');
        HDR = sclose(HDR);
        
        % save classlabels
        if isfield(HDR,'Classlabel'),
                if ~isempty(HDR.Classlabel),
                        fid = fopen([HDR.FileName(1:length(HDR.FileName)-4) '.par'],'w');
                        fprintf(fid, '%i\r\n', HDR.Classlabel);
                        fclose(fid);
                end;
        end;
        if isfield(HDR,'ArtifactSelection'),
                fid = fopen([HDR.FileName(1:length(HDR.FileName)-4) '.sel'],'w');
                fprintf(fid, '%i\r\n', HDR.ArtifactSelection);
                fclose(fid);
        end;
        
        % final test 
        try
                HDR = sopen(HDR.FileName,'r');
                HDR = sclose(HDR);
        catch
                fprintf(2,'Error SAVE2BKR: saving file %s failed\n',HDR.FileName);
        end;
end;