matlab对常见特殊数据文件读取：edf、cnf、bdf、tdms、e |
您所在的位置：网站首页 › patran导入bdf文件 › matlab对常见特殊数据文件读取：edf、cnf、bdf、tdms、e |
matlab对常见特殊数据文件读取：edf、cnf、bdf、tdms、e
|
/公众号：好玩的Matlab
小编是个从事脑电、肌电、心电方面的，在处理数据时候遇到edf文件格式的处理，经过查询资料终于找到了数据处理方法。 网址：European Data Format (EDF) EEGLAB工具箱：EEGLAB 目录 Matlab读取edf文件方法 Matlab读取bdf数据格式的方法 Matlab读取cnf文件方法 tdms文件 e文件 Matlab读取edf文件方法 function [hdr, record] = edfread(fname, varargin) % Read European Data Format file into MATLAB % % [hdr, record] = edfread(fname) % Reads data from ALL RECORDS of file fname ('*.edf'). Header % information is returned in structure hdr, and the signals % (waveforms) are returned in structure record, with waveforms % associated with the records returned as fields titled 'data' of % structure record. % % [...] = edfread(fname, 'assignToVariables', assignToVariables) % If assignToVariables is true, triggers writing of individual % output variables, as defined by field 'labels', into the caller % workspace. % % [...] = edfread(...,'targetSignals',targetSignals) % Allows user to specify the names (or position numbers) of the % subset of signals to be read. |targetSignals| may be either a % string, a cell array of comma-separated strings, or a vector of % numbers. (Default behavior is to read all signals.) % E.g.: % data = edfread(mydata.edf,'targetSignals','Thoracic'); % data = edfread(mydata.edf,'targetSignals',{'Thoracic1','Abdominal'}); % or % data = edfread(mydata.edf,'targetSignals',[2,4,6:13]); % % FORMAT SPEC: Source: http://www.edfplus.info/specs/edf.html SEE ALSO: % http://www.dpmi.tu-graz.ac.at/~schloegl/matlab/eeg/edf_spec.htm % % The first 256 bytes of the header record specify the version number of % this format, local patient and recording identification, time information % about the recording, the number of data records and finally the number of % signals (ns) in each data record. Then for each signal another 256 bytes % follow in the header record, each specifying the type of signal (e.g. % EEG, body temperature, etc.), amplitude calibration and the number of % samples in each data record (from which the sampling frequency can be % derived since the duration of a data record is also known). In this way, % the format allows for different gains and sampling frequencies for each % signal. The header record contains 256 + (ns * 256) bytes. % % Following the header record, each of the subsequent data records contains % 'duration' seconds of 'ns' signals, with each signal being represented by % the specified (in the header) number of samples. In order to reduce data % size and adapt to commonly used software for acquisition, processing and % graphical display of polygraphic signals, each sample value is % represented as a 2-byte integer in 2's complement format. Figure 1 shows % the detailed format of each data record. % % DATA SOURCE: Signals of various types (including the sample signal used % below) are available from PHYSIONET: http://www.physionet.org/ % % % % EXAMPLE 1: % % Read all waveforms/data associated with file 'ecgca998.edf': % % [header, recorddata] = edfread('ecgca998.edf'); % % % EXAMPLE 2: % % Read records 3 and 5, associated with file 'ecgca998.edf': % % % E % header = edfread('ecgca998.edf','AssignToVariables',true); % % Header file specifies data labels 'label_1'...'label_n'; these are % % created as variables in the caller workspace. % % Coded 8/27/09 by Brett Shoelson, PhD % [email protected] % Copyright 2009 - 2012 MathWorks, Inc. % % Modifications: % 5/6/13 Fixed a problem with a poorly subscripted variable. (Under certain % conditions, data were being improperly written to the 'records' variable. % Thanks to Hisham El Moaqet for reporting the problem and for sharing a % file that helped me track it down.) % % 5/22/13 Enabled import of a user-selected subset of signals. Thanks to % Farid and Cindy for pointing out the deficiency. Also fixed the import of % signals that had "bad" characters (spaces, etc) in their names. % % 10/30/14 Now outputs frequency field directly, and (thanks to Olga Imas % for the suggestion.) % % 3/23/2015 Fixed a doc problem with an improperly named variable % (desiredSignals->targetSignals). % HEADER RECORD % 8 ascii : version of this data format (0) % 80 ascii : local patient identification % 80 ascii : local recording identification % 8 ascii : startdate of recording (dd.mm.yy) % 8 ascii : starttime of recording (hh.mm.ss) % 8 ascii : number of bytes in header record % 44 ascii : reserved % 8 ascii : number of data records (-1 if unknown) % 8 ascii : duration of a data record, in seconds % 4 ascii : number of signals (ns) in data record % ns * 16 ascii : ns * label (e.g. EEG FpzCz or Body temp) % ns * 80 ascii : ns * transducer type (e.g. AgAgCl electrode) % ns * 8 ascii : ns * physical dimension (e.g. uV or degreeC) % ns * 8 ascii : ns * physical minimum (e.g. -500 or 34) % ns * 8 ascii : ns * physical maximum (e.g. 500 or 40) % ns * 8 ascii : ns * digital minimum (e.g. -2048) % ns * 8 ascii : ns * digital maximum (e.g. 2047) % ns * 80 ascii : ns * prefiltering (e.g. HP:0.1Hz LP:75Hz) % ns * 8 ascii : ns * nr of samples in each data record % ns * 32 ascii : ns * reserved % DATA RECORD % nr of samples[1] * integer : first signal in the data record % nr of samples[2] * integer : second signal % .. % .. % nr of samples[ns] * integer : last signal if nargin > 5 error('EDFREAD: Too many input arguments.'); end if ~nargin error('EDFREAD: Requires at least one input argument (filename to read).'); end [fid,msg] = fopen(fname,'r'); if fid == -1 error(msg) end assignToVariables = false; %Default targetSignals = []; %Default for ii = 1:2:numel(varargin) switch lower(varargin{ii}) case 'assigntovariables' assignToVariables = varargin{ii+1}; case 'targetsignals' targetSignals = varargin{ii+1}; otherwise error('EDFREAD: Unrecognized parameter-value pair specified. Valid values are ''assignToVariables'' and ''targetSignals''.') end end % HEADER hdr.ver = str2double(char(fread(fid,8)')); hdr.patientID = fread(fid,80,'*char')'; hdr.recordID = fread(fid,80,'*char')'; hdr.startdate = fread(fid,8,'*char')';% (dd.mm.yy) % hdr.startdate = datestr(datenum(fread(fid,8,'*char')','dd.mm.yy'), 29); %'yyyy-mm-dd' (ISO 8601) hdr.starttime = fread(fid,8,'*char')';% (hh.mm.ss) % hdr.starttime = datestr(datenum(fread(fid,8,'*char')','hh.mm.ss'), 13); %'HH:MM:SS' (ISO 8601) hdr.bytes = str2double(fread(fid,8,'*char')'); reserved = fread(fid,44);%#ok hdr.records = str2double(fread(fid,8,'*char')'); if hdr.records == -1 beep disp('There appears to be a problem with this file; it returns an out-of-spec value of -1 for ''numberOfRecords.''') disp('Attempting to read the file with ''edfReadUntilDone'' instead....'); [hdr, record] = edfreadUntilDone(fname, varargin); return end hdr.duration = str2double(fread(fid,8,'*char')'); % Number of signals hdr.ns = str2double(fread(fid,4,'*char')'); for ii = 1:hdr.ns hdr.label{ii} = regexprep(fread(fid,16,'*char')','\W',''); end if isempty(targetSignals) targetSignals = 1:numel(hdr.label); elseif iscell(targetSignals)||ischar(targetSignals) targetSignals = find(ismember(hdr.label,regexprep(targetSignals,'\W',''))); end if isempty(targetSignals) error('EDFREAD: The signal(s) you requested were not detected.') end for ii = 1:hdr.ns hdr.transducer{ii} = fread(fid,80,'*char')'; end % Physical dimension for ii = 1:hdr.ns hdr.units{ii} = fread(fid,8,'*char')'; end % Physical minimum for ii = 1:hdr.ns hdr.physicalMin(ii) = str2double(fread(fid,8,'*char')'); end % Physical maximum for ii = 1:hdr.ns hdr.physicalMax(ii) = str2double(fread(fid,8,'*char')'); end % Digital minimum for ii = 1:hdr.ns hdr.digitalMin(ii) = str2double(fread(fid,8,'*char')'); end % Digital maximum for ii = 1:hdr.ns hdr.digitalMax(ii) = str2double(fread(fid,8,'*char')'); end for ii = 1:hdr.ns hdr.prefilter{ii} = fread(fid,80,'*char')'; end for ii = 1:hdr.ns hdr.samples(ii) = str2double(fread(fid,8,'*char')'); end for ii = 1:hdr.ns reserved = fread(fid,32,'*char')';%#ok end hdr.label = hdr.label(targetSignals); hdr.label = regexprep(hdr.label,'\W',''); hdr.units = regexprep(hdr.units,'\W',''); hdr.frequency = hdr.samples./hdr.duration; disp('Step 1 of 2: Reading requested records. (This may take a few minutes.)...'); if nargout > 1 || assignToVariables % Scale data (linear scaling) scalefac = (hdr.physicalMax - hdr.physicalMin)./(hdr.digitalMax - hdr.digitalMin); dc = hdr.physicalMax - scalefac .* hdr.digitalMax; % RECORD DATA REQUESTED tmpdata = struct; for recnum = 1:hdr.records for ii = 1:hdr.ns % Read or skip the appropriate number of data points if ismember(ii,targetSignals) % Use a cell array for DATA because number of samples may vary % from sample to sample tmpdata(recnum).data{ii} = fread(fid,hdr.samples(ii),'int16') * scalefac(ii) + dc(ii); else fseek(fid,hdr.samples(ii)*2,0); end end end hdr.units = hdr.units(targetSignals); hdr.physicalMin = hdr.physicalMin(targetSignals); hdr.physicalMax = hdr.physicalMax(targetSignals); hdr.digitalMin = hdr.digitalMin(targetSignals); hdr.digitalMax = hdr.digitalMax(targetSignals); hdr.prefilter = hdr.prefilter(targetSignals); hdr.transducer = hdr.transducer(targetSignals); record = zeros(numel(hdr.label), hdr.samples(1)*hdr.records); % NOTE: 5/6/13 Modified for loop below to change instances of hdr.samples to % hdr.samples(ii). I think this underscored a problem with the reader. % Sort by requested targetSignals order: disp('Step 2 of 2: Parsing data...'); recnum = 1; for ii = 1:hdr.ns if ismember(ii,targetSignals) ctr = 1; for jj = 1:hdr.records try %#ok record(recnum, ctr : ctr + hdr.samples(ii) - 1) = tmpdata(jj).data{ii}; end ctr = ctr + hdr.samples(ii); end recnum = recnum + 1; end end hdr.ns = numel(hdr.label); hdr.samples = hdr.samples(targetSignals); if assignToVariables for ii = 1:numel(hdr.label) try %#ok eval(['assignin(''caller'',''',hdr.label{ii},''',record(ii,:))']) end end % Uncomment this line to duplicate output in a single matrix % ''record''. (Could be memory intensive if dataset is large.) record = [];% No need to output record data as a matrix? end end fclose(fid); Matlab读取bdf数据格式的方法 function EEG = readbdfdata(filename, pathname) % % Syntax: EEG = readbdfdata(pathname) % % % Inputs: % filename: % pathname: path to data files % Outputs: % EEG data structure % % Example: % % % % Subfunctions: read_bdf, readLowLevel, readEvents % MAT-files required: none % % Author: Xiaoshan Huang, [email protected] % % Versions: % v0.1: 2017-09-27, orignal % v0.2: 2018-08-07, fix update events bug by FANG Junying % Copyright (c) 2017 Neuracle, Inc. All Rights Reserved. http://neuracle.cn %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %find the index of data file index = []; if ~iscell(filename) filename = {filename}; end for i = 1:length(filename) if length(filename{i}) >= 4 %length('data')==4 if isequal(filename{i}(1:4), 'data') index = [index, i]; end end end datafilelength = length(index); %read data files if datafilelength >0 %add information for eeglab structure EEG.setname = 'BDF file'; EEG.filename = char(filename); EEG.filepath = pathname; EEG.subject = ''; EEG.group = ''; EEG.condition = ''; EEG.session = []; EEG.comments = ['Original files:' pathname]; EEG.icaact = []; EEG.icawinv = []; EEG.icasphere = []; EEG.icaweights = []; EEG.icachansind = []; chaninfo.plotrad = []; chaninfo.shrink = []; chaninfo.nosedir = '+X'; chaninfo.icachansind = []; EEG.chaninfo = chaninfo; EEG.ref = 'common'; EEG.specdata = []; EEG.specicaact = []; EEG.splinefile = ''; EEG.icasplinefile = ''; EEG.dipfit = []; EEG.xmin = 0; EEG.xmax = 0; EEG.pnts = 0; T0 = []; EEG.urevent = []; EEG.eventdescription = {}; EEG.epoch = []; EEG.epochdescription = {}; EEG.data = []; Data = {}; datapnts = []; for i = 1:datafilelength datafilename = [pathname filename{index(i)}]; %read header hdr = read_bdf(datafilename); EEG.srate = hdr.Fs; EEG.nbchan = hdr.nChans; datapnts = [datapnts,hdr.nSamples]; EEG.pnts = EEG.pnts+hdr.nSamples; EEG.chanlocs= hdr.chanlocs; EEG.trials = hdr.nTrials; EEG.xmax = EEG.xmax+EEG.pnts/EEG.srate; %in seconds T0 = [T0; hdr.T0]; %compatible with older version data file which has event channel if hdr.Annotation == 1 evt = hdr.event; evt = cell2mat(evt); EEG.event = struct([]); if isstruct(evt) for kk = 1:size(evt,2) event = struct('type',evt(kk).eventvalue,'latency',round(evt(kk).offset_in_sec*EEG.srate)); EEG.event = [EEG.event;event]; end end end if datafilelength == 1 EEG.data = read_bdf(datafilename,hdr,1,hdr.nSamples); else Data{i} = read_bdf(datafilename,hdr,1,hdr.nSamples); end end EEG.times = [1:EEG.pnts]*1000/EEG.srate; for j=1:size(T0,1) startsecs(j) = T0(j,1)*946080000+T0(j,2)*2592000+T0(j,3)*86400+T0(j,4)*3600+T0(j,5)*60+T0(j,6); end [startsecs,indexsorted] = sort(startsecs); T0 = T0(indexsorted,:); if datafilelength > 1 for k = 1:datafilelength EEG.data = [EEG.data,Data{indexsorted(k)}]; end end etc.T0 = T0(1,:); EEG.etc = etc; %read event ind = find(ismember(filename,'evt.bdf')); if numel(ind) eventfilename = [pathname filename{ind}]; hdr = read_bdf(eventfilename); evt = hdr.event; evt = cell2mat(evt); EEG.event = struct([]); if isstruct(evt) for i = 1:size(evt,2) event = struct('type',evt(i).eventvalue,'latency',round(evt(i).offset_in_sec*EEG.srate)); EEG.event = [EEG.event;event]; end end %shift event latency for multiple data files only if datafilelength > 1 startpnts = (startsecs-startsecs(1))*EEG.srate; %relative time endpnts = zeros(1,datafilelength); pausepnts = 0; Event = struct([]); datapnts = datapnts(indexsorted); for j = 1:datafilelength endpnts(j) = startpnts(j)+datapnts(j); if j > 1 pausepnts = pausepnts+(startpnts(j)-endpnts(j-1)); end for k = 1:length(EEG.event) if startpnts(j) = EEG.event(k).latency EEG.event(k).latency = EEG.event(k).latency-pausepnts; Event = [Event;EEG.event(k)]; end end if j < datafilelength eventBoundary = struct('type','boundary','latency',endpnts(j)+1); Event = [Event;eventBoundary]; end end EEG.event = Event; end elseif exist([pathname 'evt.bdf'],'file') == 0 && hdr.Annotation == 0 EEG.event = struct([]); end %read mems data and combine with EEG data memstype = {'acc.edf','gyro.edf','mag.edf'}; for gg = 1:length(memstype) %initialize indexmems = []; T0mems = []; memsData = {}; mems.data = []; %find the index of the mems file for ii = 1:length(filename) lennametype = length(memstype{gg})-4; %ignore '.edf' if length(filename{ii}) >= lennametype if isequal(filename{ii}(1:lennametype), memstype{gg}(1:lennametype)) indexmems = [indexmems, ii]; end end end filelength = length(indexmems); if filelength >0 for ii = 1:filelength memsfilename = [pathname filename{indexmems(ii)}]; %read header hdrmems = read_bdf(memsfilename); chanlocs= hdrmems.chanlocs; nbchan = hdrmems.nChans; T0mems = [T0mems; hdrmems.T0]; %read data data = read_bdf(memsfilename,hdrmems,1,hdrmems.nSamples); %resample datarsp = []; for jj = 1:nbchan datamems = resample(data(jj,:), EEG.srate, hdrmems.Fs); datarsp = [datarsp;datamems]; end if filelength == 1 mems.data = datarsp; else memsData{ii} = datarsp; end end if filelength > 1 %sort the mems data files based on start time (T0) for aa=1:size(T0mems,1) startsecs(aa) = T0mems(aa,1)*946080000+T0mems(aa,2)*2592000+T0mems(aa,3)*86400+T0mems(aa,4)*3600+T0mems(aa,5)*60+T0mems(aa,6); end [startsecs,indexsorted] = sort(startsecs); for bb = 1:filelength mems.data = [mems.data,memsData{indexsorted(bb)}]; end end %make sure EEG and mems data have same length nbpnts = size(mems.data,2) ; if nbpnts == EEG.pnts EEG.data = [EEG.data; mems.data]; elseif nbpnts > EEG.pnts EEG.data = [EEG.data; mems.data(:,1:EEG.pnts)]; elseif nbpnts < EEG.pnts EEG.data = [EEG.data(:,1:hdr.nSamples); mems.data]; end EEG.chanlocs = [EEG.chanlocs ; chanlocs ]; EEG.nbchan = EEG.nbchan + nbchan; end end else error('Please select EEG datasets with the filenames starting by "data" ') end function dat = read_bdf(filename, hdr, begsample, endsample, chanindx) if nargin==1 % read the header, this code is from EEGLAB's openbdf FILENAME = filename; % defines Seperator for Subdirectories SLASH='/'; BSLASH=char(92); cname=computer; if cname(1:2)=='PC' SLASH=BSLASH; end; fid=fopen(FILENAME,'r','ieee-le'); if fid= EDF.DigMax)) fprintf(2,'Warning OPENEDF: Digital Minimum larger than Maximum\n'); end % check validity of PhysMin and PhysMax if (length(EDF.PhysMin) ~= EDF.NS) fprintf(2,'Warning OPENEDF: Failing Physical Minimum\n'); EDF.PhysMin = EDF.DigMin; end if (length(EDF.PhysMax) ~= EDF.NS) fprintf(2,'Warning OPENEDF: Failing Physical Maximum\n'); EDF.PhysMax = EDF.DigMax; end if (any(EDF.PhysMin >= EDF.PhysMax)) fprintf(2,'Warning OPENEDF: Physical Minimum larger than Maximum\n'); EDF.PhysMin = EDF.DigMin; EDF.PhysMax = EDF.DigMax; end EDF.PreFilt= char(fread(EDF.FILE.FID,[80,EDF.NS],'char')'); % tmp = fread(EDF.FILE.FID,[8,EDF.NS],'char')'; % samples per data record EDF.SPR = str2num(char(tmp)); % samples per data record fseek(EDF.FILE.FID,32*EDF.NS,0); EDF.Cal = (EDF.PhysMax-EDF.PhysMin)./(EDF.DigMax-EDF.DigMin); EDF.Off = EDF.PhysMin - EDF.Cal .* EDF.DigMin; tmp = find(EDF.Cal < 0); EDF.Cal(tmp) = ones(size(tmp)); EDF.Off(tmp) = zeros(size(tmp)); EDF.Calib=[EDF.Off';(diag(EDF.Cal))]; %EDF.Calib=sparse(diag([1; EDF.Cal])); %EDF.Calib(1,2:EDF.NS+1)=EDF.Off'; EDF.SampleRate = EDF.SPR / EDF.Dur; EDF.FILE.POS = ftell(EDF.FILE.FID); if EDF.NRec == -1 % unknown record size, determine correct NRec fseek(EDF.FILE.FID, 0, 'eof'); endpos = ftell(EDF.FILE.FID); EDF.NRec = floor((endpos - EDF.FILE.POS) / (sum(EDF.SPR) * 2)); fseek(EDF.FILE.FID, EDF.FILE.POS, 'bof'); H1(237:244)=sprintf('%-8i',EDF.NRec); % write number of records end; EDF.Chan_Select=(EDF.SPR==max(EDF.SPR)); for k=1:EDF.NS if EDF.Chan_Select(k) EDF.ChanTyp(k)='N'; else EDF.ChanTyp(k)=' '; end; if findstr(upper(EDF.Label(k,:)),'ECG') EDF.ChanTyp(k)='C'; elseif findstr(upper(EDF.Label(k,:)),'EKG') EDF.ChanTyp(k)='C'; elseif findstr(upper(EDF.Label(k,:)),'EEG') EDF.ChanTyp(k)='E'; elseif findstr(upper(EDF.Label(k,:)),'EOG') EDF.ChanTyp(k)='O'; elseif findstr(upper(EDF.Label(k,:)),'EMG') EDF.ChanTyp(k)='M'; end; end; EDF.AS.spb = sum(EDF.SPR); % Samples per Block hdr.Annotation = 0; hdr.AnnotationChn = -1; if any(EDF.SampleRate~=EDF.SampleRate(1)) chn_indx = find(EDF.SampleRate~=EDF.SampleRate(1)); if(length(chn_indx) > 0 && (strcmp(EDF.Label(chn_indx,:),repmat('BDF Annotations',length(chn_indx),1))) || strcmp(EDF.Label(chn_indx,:),repmat('BDF Annotations ',length(chn_indx),1))) disp('BDF+ file format detected, the BDF Annotation channel will be skipped when reading data'); hdr.Annotation = 1; hdr.AnnotationChn = chn_indx; EDF.SampleRateOrg = EDF.SampleRate; EDF.SampleRate(chn_indx) = EDF.SampleRate(1);%force it to be the same as others, in order to pass the FieldTrip test EDF.Label(chn_indx,:) = [];%delete the info about the Annotation channel EDF.NS = EDF.NS - length(chn_indx); else error('Channels with different sampling rate but not BDF Annotation found'); end end %read all the events as well if(hdr.Annotation) epoch_total = EDF.Dur * sum(EDF.SampleRateOrg); %the number of total data per record, usually one second % allocate memory to hold the data n_annotation = sum(EDF.Dur*EDF.SampleRateOrg(hdr.AnnotationChn)); %the number of annotation per record dat_event = zeros(EDF.NRec,n_annotation*3); % read and concatenate all required data epochs for i=1:EDF.NRec %offset: the bit index right before the current first event offset = EDF.HeadLen + (i-1)*epoch_total*3 + EDF.Dur*sum(EDF.SampleRateOrg(1:EDF.NS))*3; %extract the bits for events buf = readEvents(filename, offset, n_annotation); % see below in subfunction % buf = readLowLevel(filename, offset, n_annotation); dat_event(i,:) = buf; end %decoding the events event_cnt = 1; event = []; for i = 1:EDF.NRec %last index char20_index = find(dat_event(i,:)==20); TAL_start = char20_index(2) + 1;%the first event right after the second char 20 char20_index = char20_index(3:end);%remove the first two 20 (belonging to the TAL time index) num_event = length(char20_index)/2;%two char 20 per event if(num_event == 0) continue; end for j = 1:num_event %structure of a single event: [offset in sec] [char21][Duration in sec] [char20][Trigger Code] [char20][char0] %the field [char21][Duration in sec] can be skipped %multiple events can be stored within one Annotation block, the last ends with [char0][char0] singleTAL = dat_event(i,TAL_start:char20_index(2*j)-1); char21_one = find(singleTAL == 21); char20_one = find(singleTAL == 20); event{event_cnt}.eventtype = 'trigger';%fixed info event{event_cnt}.eventvalue = char(singleTAL(char20_one+1:end));%trigger code, char if(~isempty(char21_one))%if Duration field exist event{event_cnt}.offset_in_sec = str2num(char(singleTAL(1:char21_one-1)));%in sec event{event_cnt}.duration = str2num(char(singleTAL(char21_one+1:char20_one-1)));%in sec; else%if no Duration field event{event_cnt}.offset_in_sec = str2num(char(singleTAL(1:char20_one-1)));%in sec event{event_cnt}.duration = 0; end event{event_cnt}.offset = round(event{event_cnt}.offset_in_sec*EDF.SampleRateOrg(1));%in sampling points event_cnt = event_cnt + 1; TAL_start = char20_index(2*j) + 1; if(dat_event(i,char20_index(2*j)+1) ~= 0)%check if the TAL is complete error('BDF+ event error'); end end %check if all TALs are read if ~(dat_event(i,char20_index(2*num_event)+1) == 0 && dat_event(i,char20_index(2*num_event)+2) == 0) error('BDF+ final event error'); end end hdr.event = event; end % close the file fclose(EDF.FILE.FID); % convert the header to Fieldtrip-style hdr.Fs = EDF.SampleRate(1); hdr.nChans = EDF.NS; hdr.label = cellstr(EDF.Label); hdr.chanlocs = struct([]); for i = 1:hdr.nChans chanlocs = struct('labels',hdr.label(i),'ref',[],'theta',[],'radius',[],'X',[],'Y',[],'Z',[],'sph_theta',[],'sph_phi',[],'sph_radius',[],'type',[],'urchan',[] ); hdr.chanlocs = [hdr.chanlocs;chanlocs]; end % it is continuous data, therefore append all records in one trial hdr.nTrials = 1; hdr.nSamples = EDF.NRec * EDF.Dur * EDF.SampleRate(1); hdr.nSamplesPre = 0; hdr.orig = EDF; % return the header dat = hdr; else % read the data % retrieve the original header EDF = hdr.orig; % determine the trial containing the begin and end sample epochlength = EDF.Dur * EDF.SampleRate(1); begepoch = floor((begsample-1)/epochlength) + 1; endepoch = floor((endsample-1)/epochlength) + 1; nepochs = endepoch - begepoch + 1; nchans = EDF.NS; if(hdr.Annotation) epoch_total = sum(EDF.Dur * EDF.SampleRateOrg); else epoch_total = EDF.Dur * EDF.SampleRate(1)*nchans; end epoch_data = EDF.Dur * EDF.SampleRate(1)*nchans; %TODO HERE if nargin= hdr.AnnotationChn) chanindx(tmp_indx) = chanindx(tmp_indx) + 1; end end end % allocate memory to hold the data dat = zeros(length(chanindx),nepochs*epochlength); % read and concatenate all required data epochs for i=begepoch:endepoch%number of records if hdr.filetype == 0 offset = EDF.HeadLen + (i-1)*epoch_total*2; else offset = EDF.HeadLen + (i-1)*epoch_total*3; end % read the data from all channels and then select the desired channels buf = readLowLevel(filename, offset, epoch_data,hdr.filetype); % see below in subfunction buf = reshape(buf, epochlength, nchans); dat(:,((i-begepoch)*epochlength+1):((i-begepoch+1)*epochlength)) = buf(:,chanindx)'; end % select the desired samples begsample = begsample - (begepoch-1)*epochlength; % correct for the number of bytes that were skipped endsample = endsample - (begepoch-1)*epochlength; % correct for the number of bytes that were skipped dat = dat(:, begsample:endsample); % Calibrate the data calib = diag(EDF.Cal(chanindx)); if length(chanindx)>1 % using a sparse matrix speeds up the multiplication dat = sparse(calib) * dat; else % in case of one channel the sparse multiplication would result in a sparse array dat = calib * dat; end end % SUBFUNCTION for reading the 24 bit values function buf = readLowLevel(filename, offset, numwords,filetype) % if offset < 2*1024^3 % % use the external mex file, only works for cnt = loadcnt(file, varargin) % % Inputs: % filename - name of the file with extension % % Optional inputs: % 't1' - start at time t1, default 0. Warning, events latency % might be innacurate (this is an open issue). % 'sample1' - start at sample1, default 0, overrides t1. Warning, % events latency might be innacurate. % 'lddur' - duration of segment to load, default = whole file % 'ldnsamples' - number of samples to load, default = whole file, % overrides lddur % 'scale' - ['on'|'off'] scale data to microvolt (default:'on') % 'dataformat' - ['int16'|'int32'] default is 'int16' for 16-bit data. % Use 'int32' for 32-bit data. % 'blockread' - [integer] by default it is automatically determined % from the file header, though sometimes it finds an % incorect value, so you may want to enter a value manually % here (1 is the most standard value). % 'memmapfile' - ['memmapfile_name'] use this option if the .cnt file % is too large to read in conventially. The suffix of % the memmapfile_name must be .fdt. The memmapfile % functions process files based on their suffix, and an % error will occur if you use a different suffix. % % Outputs: % cnt - structure with the continuous data and other informations % cnt.header % cnt.electloc % cnt.data % cnt.tag % % Authors: Sean Fitzgibbon, Arnaud Delorme, 2000- % % Note: function original name was load_scan41.m % % Known limitations: % For more see http://www.cnl.salk.edu/~arno/cntload/index.html % Copyright (C) 2000 Sean Fitzgibbon, % Copyright (C) 2003 Arnaud Delorme, Salk Institute, [email protected] % % This program 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 2 of the License, or % (at your option) any later version. % % This program 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 this program; if not, write to the Free Software % Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA function [f,lab,ev2p] = loadcnt(filename,varargin) if ~isempty(varargin) r=struct(varargin{:}); else r = []; end; try, r.t1; catch, r.t1=0; end try, r.sample1; catch, r.sample1=[]; end try, r.lddur; catch, r.lddur=[]; end try, r.ldnsamples; catch, r.ldnsamples=[]; end try, r.scale; catch, r.scale='on'; end try, r.blockread; catch, r.blockread = []; end try, r.dataformat; catch, r.dataformat = 'int32'; end try, r.memmapfile; catch, r.memmapfile = ''; end sizeEvent1 = 8 ; %%% 8 bytes for Event1 sizeEvent2 = 19 ; %%% 19 bytes for Event2 sizeEvent3 = 19 ; %%% 19 bytes for Event3 type='cnt'; if nargin ==1 scan=0; end fid = fopen(filename,'r','l'); disp(['Loading file ' filename ' ...']) h.rev = fread(fid,12,'char'); h.nextfile = fread(fid,1,'long'); h.prevfile = fread(fid,1,'ulong'); h.type = fread(fid,1,'char'); h.id = fread(fid,20,'char'); h.oper = fread(fid,20,'char'); h.doctor = fread(fid,20,'char'); h.referral = fread(fid,20,'char'); h.hospital = fread(fid,20,'char'); h.patient = fread(fid,20,'char'); h.age = fread(fid,1,'short'); h.sex = fread(fid,1,'char'); h.hand = fread(fid,1,'char'); h.med = fread(fid,20, 'char'); h.category = fread(fid,20, 'char'); h.state = fread(fid,20, 'char'); h.label = fread(fid,20, 'char'); h.date = fread(fid,10, 'char'); h.time = fread(fid,12, 'char'); h.mean_age = fread(fid,1,'float'); h.stdev = fread(fid,1,'float'); h.n = fread(fid,1,'short'); h.compfile = fread(fid,38,'char'); h.spectwincomp = fread(fid,1,'float'); h.meanaccuracy = fread(fid,1,'float'); h.meanlatency = fread(fid,1,'float'); h.sortfile = fread(fid,46,'char'); h.numevents = fread(fid,1,'int'); h.compoper = fread(fid,1,'char'); h.avgmode = fread(fid,1,'char'); h.review = fread(fid,1,'char'); h.nsweeps = fread(fid,1,'ushort'); h.compsweeps = fread(fid,1,'ushort'); h.acceptcnt = fread(fid,1,'ushort'); h.rejectcnt = fread(fid,1,'ushort'); h.pnts = fread(fid,1,'ushort'); h.nchannels = fread(fid,1,'ushort'); h.avgupdate = fread(fid,1,'ushort'); h.domain = fread(fid,1,'char'); h.variance = fread(fid,1,'char'); h.rate = fread(fid,1,'ushort'); % A USER CLAIMS THAT SAMPLING RATE CAN BE h.scale = fread(fid,1,'double'); % FRACTIONAL IN NEUROSCAN WHICH IS h.veogcorrect = fread(fid,1,'char'); % OBVIOUSLY NOT POSSIBLE HERE (BUG 606) h.heogcorrect = fread(fid,1,'char'); h.aux1correct = fread(fid,1,'char'); h.aux2correct = fread(fid,1,'char'); h.veogtrig = fread(fid,1,'float'); h.heogtrig = fread(fid,1,'float'); h.aux1trig = fread(fid,1,'float'); h.aux2trig = fread(fid,1,'float'); h.heogchnl = fread(fid,1,'short'); h.veogchnl = fread(fid,1,'short'); h.aux1chnl = fread(fid,1,'short'); h.aux2chnl = fread(fid,1,'short'); h.veogdir = fread(fid,1,'char'); h.heogdir = fread(fid,1,'char'); h.aux1dir = fread(fid,1,'char'); h.aux2dir = fread(fid,1,'char'); h.veog_n = fread(fid,1,'short'); h.heog_n = fread(fid,1,'short'); h.aux1_n = fread(fid,1,'short'); h.aux2_n = fread(fid,1,'short'); h.veogmaxcnt = fread(fid,1,'short'); h.heogmaxcnt = fread(fid,1,'short'); h.aux1maxcnt = fread(fid,1,'short'); h.aux2maxcnt = fread(fid,1,'short'); h.veogmethod = fread(fid,1,'char'); h.heogmethod = fread(fid,1,'char'); h.aux1method = fread(fid,1,'char'); h.aux2method = fread(fid,1,'char'); h.ampsensitivity = fread(fid,1,'float'); h.lowpass = fread(fid,1,'char'); h.highpass = fread(fid,1,'char'); h.notch = fread(fid,1,'char'); h.autoclipadd = fread(fid,1,'char'); h.baseline = fread(fid,1,'char'); h.offstart = fread(fid,1,'float'); h.offstop = fread(fid,1,'float'); h.reject = fread(fid,1,'char'); h.rejstart = fread(fid,1,'float'); h.rejstop = fread(fid,1,'float'); h.rejmin = fread(fid,1,'float'); h.rejmax = fread(fid,1,'float'); h.trigtype = fread(fid,1,'char'); h.trigval = fread(fid,1,'float'); h.trigchnl = fread(fid,1,'char'); h.trigmask = fread(fid,1,'short'); h.trigisi = fread(fid,1,'float'); h.trigmin = fread(fid,1,'float'); h.trigmax = fread(fid,1,'float'); h.trigdir = fread(fid,1,'char'); h.autoscale = fread(fid,1,'char'); h.n2 = fread(fid,1,'short'); h.dir = fread(fid,1,'char'); h.dispmin = fread(fid,1,'float'); h.dispmax = fread(fid,1,'float'); h.xmin = fread(fid,1,'float'); h.xmax = fread(fid,1,'float'); h.automin = fread(fid,1,'float'); h.automax = fread(fid,1,'float'); h.zmin = fread(fid,1,'float'); h.zmax = fread(fid,1,'float'); h.lowcut = fread(fid,1,'float'); h.highcut = fread(fid,1,'float'); h.common = fread(fid,1,'char'); h.savemode = fread(fid,1,'char'); h.manmode = fread(fid,1,'char'); h.ref = fread(fid,10,'char'); h.rectify = fread(fid,1,'char'); h.displayxmin = fread(fid,1,'float'); h.displayxmax = fread(fid,1,'float'); h.phase = fread(fid,1,'char'); h.screen = fread(fid,16,'char'); h.calmode = fread(fid,1,'short'); h.calmethod = fread(fid,1,'short'); h.calupdate = fread(fid,1,'short'); h.calbaseline = fread(fid,1,'short'); h.calsweeps = fread(fid,1,'short'); h.calattenuator = fread(fid,1,'float'); h.calpulsevolt = fread(fid,1,'float'); h.calpulsestart = fread(fid,1,'float'); h.calpulsestop = fread(fid,1,'float'); h.calfreq = fread(fid,1,'float'); h.taskfile = fread(fid,34,'char'); h.seqfile = fread(fid,34,'char'); h.spectmethod = fread(fid,1,'char'); h.spectscaling = fread(fid,1,'char'); h.spectwindow = fread(fid,1,'char'); h.spectwinlength = fread(fid,1,'float'); h.spectorder = fread(fid,1,'char'); h.notchfilter = fread(fid,1,'char'); h.headgain = fread(fid,1,'short'); h.additionalfiles = fread(fid,1,'int'); h.unused = fread(fid,5,'char'); h.fspstopmethod = fread(fid,1,'short'); h.fspstopmode = fread(fid,1,'short'); h.fspfvalue = fread(fid,1,'float'); h.fsppoint = fread(fid,1,'short'); h.fspblocksize = fread(fid,1,'short'); h.fspp1 = fread(fid,1,'ushort'); h.fspp2 = fread(fid,1,'ushort'); h.fspalpha = fread(fid,1,'float'); h.fspnoise = fread(fid,1,'float'); h.fspv1 = fread(fid,1,'short'); h.montage = fread(fid,40,'char'); h.eventfile = fread(fid,40,'char'); h.fratio = fread(fid,1,'float'); h.minor_rev = fread(fid,1,'char'); h.eegupdate = fread(fid,1,'short'); h.compressed = fread(fid,1,'char'); h.xscale = fread(fid,1,'float'); h.yscale = fread(fid,1,'float'); h.xsize = fread(fid,1,'float'); h.ysize = fread(fid,1,'float'); h.acmode = fread(fid,1,'char'); h.commonchnl = fread(fid,1,'uchar'); h.xtics = fread(fid,1,'char'); h.xrange = fread(fid,1,'char'); h.ytics = fread(fid,1,'char'); h.yrange = fread(fid,1,'char'); h.xscalevalue = fread(fid,1,'float'); h.xscaleinterval = fread(fid,1,'float'); h.yscalevalue = fread(fid,1,'float'); h.yscaleinterval = fread(fid,1,'float'); h.scaletoolx1 = fread(fid,1,'float'); h.scaletooly1 = fread(fid,1,'float'); h.scaletoolx2 = fread(fid,1,'float'); h.scaletooly2 = fread(fid,1,'float'); h.port = fread(fid,1,'short'); h.numsamples = fread(fid,1,'ulong'); h.filterflag = fread(fid,1,'char'); h.lowcutoff = fread(fid,1,'float'); h.lowpoles = fread(fid,1,'short'); h.highcutoff = fread(fid,1,'float'); h.highpoles = fread(fid,1,'short'); h.filtertype = fread(fid,1,'char'); h.filterdomain = fread(fid,1,'char'); h.snrflag = fread(fid,1,'char'); h.coherenceflag = fread(fid,1,'char'); h.continuoustype = fread(fid,1,'char'); h.eventtablepos = fread(fid,1,'ulong'); h.continuousseconds = fread(fid,1,'float'); h.channeloffset = fread(fid,1,'long'); h.autocorrectflag = fread(fid,1,'char'); h.dcthreshold = fread(fid,1,'uchar'); for n = 1:h.nchannels e(n).lab = deblank(char(fread(fid,10,'char')')); e(n).reference = fread(fid,1,'char'); e(n).skip = fread(fid,1,'char'); e(n).reject = fread(fid,1,'char'); e(n).display = fread(fid,1,'char'); e(n).bad = fread(fid,1,'char'); e(n).n = fread(fid,1,'ushort'); e(n).avg_reference = fread(fid,1,'char'); e(n).clipadd = fread(fid,1,'char'); e(n).x_coord = fread(fid,1,'float'); e(n).y_coord = fread(fid,1,'float'); e(n).veog_wt = fread(fid,1,'float'); e(n).veog_std = fread(fid,1,'float'); e(n).snr = fread(fid,1,'float'); e(n).heog_wt = fread(fid,1,'float'); e(n).heog_std = fread(fid,1,'float'); e(n).baseline = fread(fid,1,'short'); e(n).filtered = fread(fid,1,'char'); e(n).fsp = fread(fid,1,'char'); e(n).aux1_wt = fread(fid,1,'float'); e(n).aux1_std = fread(fid,1,'float'); e(n).senstivity = fread(fid,1,'float'); e(n).gain = fread(fid,1,'char'); e(n).hipass = fread(fid,1,'char'); e(n).lopass = fread(fid,1,'char'); e(n).page = fread(fid,1,'uchar'); e(n).size = fread(fid,1,'uchar'); e(n).impedance = fread(fid,1,'uchar'); e(n).physicalchnl = fread(fid,1,'uchar'); e(n).rectify = fread(fid,1,'char'); e(n).calib = fread(fid,1,'float'); end % finding if 32-bits of 16-bits file % ---------------------------------- begdata = ftell(fid); if strcmpi(r.dataformat, 'auto') r.dataformat = 'int16'; if (h.nextfile > 0) fseek(fid,h.nextfile+52,'bof'); is32bit = fread(fid,1,'char'); if (is32bit == 1) r.dataformat = 'int32' end; fseek(fid,begdata,'bof'); end; end; enddata = h.eventtablepos; % after data if strcmpi(r.dataformat, 'int16') nums = (enddata-begdata)/h.nchannels/2; else nums = (enddata-begdata)/h.nchannels/4; end; % number of sample to read % ------------------------ if ~isempty(r.sample1) r.t1 = r.sample1/h.rate; else r.sample1 = r.t1*h.rate; end; if strcmpi(r.dataformat, 'int16') startpos = r.t1*h.rate*2*h.nchannels; else startpos = r.t1*h.rate*4*h.nchannels; end; if isempty(r.ldnsamples) if ~isempty(r.lddur) r.ldnsamples = round(r.lddur*h.rate); else r.ldnsamples = nums; end; end; % channel offset % -------------- if ~isempty(r.blockread) h.channeloffset = r.blockread; end; if h.channeloffset > 1 fprintf('WARNING: reading data in blocks of %d, if this fails, try using option "''blockread'', 1"\n', ... h.channeloffset); end; disp('Reading data .....') if type == 'cnt' % while (ftell(fid) +1 < h.eventtablepos) %d(:,i)=fread(fid,h.nchannels,'int16'); %end fseek(fid, startpos, 0); % **** This marks the beginning of the code modified for reading % large .cnt files % Switched to r.memmapfile for continuity. Check to see if the % variable exists. If it does, then the user has indicated the % file is too large to be processed in memory. If the variable % is blank, the file is processed in memory. if (~isempty(r.memmapfile)) % open a file for writing foutid = fopen(r.memmapfile, 'w') ; % This portion of the routine reads in a section of the EEG file % and then writes it out to the harddisk. samples_left = h.nchannels * r.ldnsamples ; % the size of the data block to be read is limited to 4M % samples. This equates to 16MB and 32MB of memory for % 16 and 32 bit files, respectively. data_block = 4000000 ; max_rows = data_block / h.nchannels ; %warning off ; max_written = h.nchannels * uint32(max_rows) ; %warning on ; % This while look tracks the remaining samples. The % data is processed in chunks rather than put into % memory whole. while (samples_left > 0) % Check to see if the remaining data is smaller than % the general processing block by looking at the % remaining number of rows. to_read = max_rows ; if (data_block > samples_left) to_read = samples_left / h.nchannels ; end ; % Read data in a relatively small chunk temp_dat = fread(fid, [h.nchannels to_read], r.dataformat) ; % The data is then scaled using the original routine. % In the original routine, the entire data set was scaled % after being read in. For this version, scaling occurs % after every chunk is read. if strcmpi(r.scale, 'on') disp('Scaling data .....') %%% scaling to microvolts for i=1:h.nchannels bas=e(i).baseline;sen=e(i).senstivity;cal=e(i).calib; mf=sen*(cal/204.8); temp_dat(i,:)=(temp_dat(i,:)-bas).*mf; end end % Write out data in float32 form to the file name % supplied by the user. written = fwrite (foutid, temp_dat, 'float32') ; if (written ~= max_written) samples_left = 0 ; else samples_left = samples_left - written ; end ; end ; fclose (foutid) ; % Set the dat variable. This gets used later by other % EEGLAB functions. dat = r.memmapfile ; % This variable tracks how the data should be read. bReadIntoMemory = false ; else % The memmapfile variable is empty, read into memory. bReadIntoMemory = true ; end % This ends the modifications made to read large files. % Everything contained within the following if statement is the % original code. if (bReadIntoMemory == true) if h.channeloffset 0 ev2(nevents).stimtype = []; for i=1:nevents ev2(i).stimtype = fread(fid,1,'ushort'); ev2(i).keyboard = fread(fid,1,'char'); temp = fread(fid,1,'uint8'); ev2(i).keypad_accept = bitand(15,temp); ev2(i).accept_ev1 = bitshift(temp,-4); ev2(i).offset = fread(fid,1,'long'); ev2(i).type = fread(fid,1,'short'); ev2(i).code = fread(fid,1,'short'); ev2(i).latency = fread(fid,1,'float'); ev2(i).epochevent = fread(fid,1,'char'); ev2(i).accept = fread(fid,1,'char'); ev2(i).accuracy = fread(fid,1,'char'); end else ev2 = []; end; elseif eT.teeg==3 % type 3 is similar to type 2 except the offset field encodes the global sample frame nevents=eT.size/sizeEvent3; if nevents > 0 ev2(nevents).stimtype = []; if r.dataformat == 'int32' bytes_per_samp = 4; % I only have 32 bit data, unable to check whether this is necessary, else % perhaps there is no type 3 file with 16 bit data bytes_per_samp = 2; end for i=1:nevents ev2(i).stimtype = fread(fid,1,'ushort'); ev2(i).keyboard = fread(fid,1,'char'); temp = fread(fid,1,'uint8'); ev2(i).keypad_accept = bitand(15,temp); ev2(i).accept_ev1 = bitshift(temp,-4); os = fread(fid,1,'ulong'); ev2(i).offset = os * bytes_per_samp * h.nchannels; ev2(i).type = fread(fid,1,'short'); ev2(i).code = fread(fid,1,'short'); ev2(i).latency = fread(fid,1,'float'); ev2(i).epochevent = fread(fid,1,'char'); ev2(i).accept = fread(fid,1,'char'); ev2(i).accuracy = fread(fid,1,'char'); end else ev2 = []; end; elseif eT.teeg==1 nevents=eT.size/sizeEvent1; if nevents > 0 ev2(nevents).stimtype = []; for i=1:nevents ev2(i).stimtype = fread(fid,1,'ushort'); ev2(i).keyboard = fread(fid,1,'char'); % modified by Andreas Widmann 2005/05/12 14:15:00 %ev2(i).keypad_accept = fread(fid,1,'char'); temp = fread(fid,1,'uint8'); ev2(i).keypad_accept = bitand(15,temp); ev2(i).accept_ev1 = bitshift(temp,-4); % end modification ev2(i).offset = fread(fid,1,'long'); end; else ev2 = []; end; else disp('Skipping event table (tag != 1,2,3 ; theoritically impossible)'); ev2 = []; end fseek(fid, -1, 'eof'); t = fread(fid,'char'); f.header = h; f.electloc = e; f.data = dat; f.Teeg = eT; f.event = ev2; f.tag=t; % Surgical addition of number of samples f.ldnsamples = r.ldnsamples ; %%%% channels labels for i=1:h.nchannels plab=sprintf('%c',f.electloc(i).lab); if i>1 lab=str2mat(lab,plab); else lab=plab; end end %%%% to change offest in bytes to points if ~isempty(ev2) if r.sample1 ~= 0 fprintf(2,'Warning: events imported with a time shift might be innacurate\n'); end; ev2p=ev2; ioff=900+(h.nchannels*75); %% initial offset : header + electordes desc if strcmpi(r.dataformat, 'int16') for i=1:nevents ev2p(i).offset=(ev2p(i).offset-ioff)/(2*h.nchannels) - r.sample1; %% 2 short int end end else % 32 bits for i=1:nevents ev2p(i).offset=(ev2p(i).offset-ioff)/(4*h.nchannels) - r.sample1; %% 4 short int end end end; f.event = ev2p; end; frewind(fid); fclose(fid); end tdms文件 function [ConvertedData,ConvertVer,ChanNames,GroupNames,ci]=convertTDMS(varargin) %Function to load LabView TDMS data file(s) into variables in the MATLAB workspace. %An *.MAT file can also be created. If called with one input, the user selects %a data file. % % TDMS format is based on information provided by National Instruments at: % http://zone.ni.com/devzone/cda/tut/p/id/5696 % % [ConvertedData,ConvertVer,ChanNames]=convertTDMS(SaveConvertedFile,filename) % % Inputs: % SaveConvertedFile (required) - Logical flag (true/false) that % determines whether a MAT file is created. The MAT file's name % is the same as 'filename' except that the 'TDMS' file extension is % replaced with 'MAT'. The MAT file is saved in the same folder % and will overwrite an existing file without warning. The % MAT file contains all the output variables. % % filename (optional) - Filename (fully defined) to be converted. % If not supplied, the user is provided a 'File Open' dialog box % to navigate to a file. Can be a cell array of files for bulk % conversion. % % Outputs: % ConvertedData (required) - Structure of all of the data objects. % ConvertVer (optional) - Version number of this function. % ChanNames (optional) - Cell array of channel names % GroupNames (optional) - Cell array of group names % ci (optional) - Structure of the channel index (an index to % where all of the information for a channel resides in a % file. % % %'ConvertedData' is a structure with 'FileName', 'FileFolder', 'SegTDMSVerNum', %'NumOfSegments' and 'Data' fields'. The 'Data' field is a structure. % %'ConvertedData.SegTDMSVerNum' is a vector of the TDMS version number for each %segment. % %'ConvertedData.Data' is a structure with 'Root' and 'MeasuredData' fields. % %'ConvertedData.Data.Root' is a structure with 'Name' and 'Property' fields. %The 'Property' field is also a structure; it contains all the specified properties %(1 entry for each 'Property) for the 'Root' group. For each 'Property' there are %'Name' and 'Value' fields. To display a list of all the property names, input %'{ConvertedData.Data.Root.Property.Name}'' in the Command Window. % %'ConvertedData.Data.MeasuredData' is a structure containing all the channel/group %information. For each index (for example, 'ConvertedData.Data.MeasuredData(1)'), %there are 'Name', 'Data' and 'Property' fields. The list of channel names can %be displayed by typing 'ChanNames' in the Command Window. Similarly, the list %of group names can be displayed by typing 'GroupNames' in the Command Window. %The 'Property' field is also a structure; it contains all the specified properties %for that index (1 entry in the structure for each 'Property'). Any LabView waveform %attributes ('wf_start_time', 'wf_start_offset', 'wf_increment' and 'wf_samples') that %may exist are also included in the properties. For each 'Property' there are 'Name' %and 'Value' fields. To display a list of all the property names, input %'{ConvertedData.Data.MeasuredData(#).Property.Name}'' in the Command Window %where '#' is the index of interest. % %If you recieve an error that DAQmxRaw data cannot be converted, this is due %to the details on parsing this data type not being published by NI. the %work around is to use a VI that reads in the data file and then writes it %to a new file. See: https://decibel.ni.com/content/docs/DOC-32817 % % See Also: simpleconvertTDMS %------------------------------------------------------------------------- %Brad Humphreys - v1.0 2008-04-23 %ZIN Technologies %------------------------------------------------------------------------- %------------------------------------------------------------------------- %Brad Humphreys - v1.1 2008-07-03 %ZIN Technologies %-Added abilty for timestamp to be a raw data type, not just meta data. %-Addressed an issue with having a default nsmaples entry for new objects. %-Added Error trap if file name not found. %-Corrected significant problem where it was assumed that once an object % existsed, it would in in every subsequent segement. This is not true. %------------------------------------------------------------------------- %------------------------------------------------------------------------- %Grant Lohsen - v1.2 2009-11-15 %Georgia Tech Research Institute %-Converts TDMS v2 files %Folks, it's not pretty but I don't have time to make it pretty. Enjoy. %------------------------------------------------------------------------- %------------------------------------------------------------------------- %Jeff Sitterle - v1.3 2010-01-10 %Georgia Tech Research Institute %Modified to return all information stored in the TDMS file to inlcude %name, start time, start time offset, samples per read, total samples, unit %description, and unit string. Also provides event time and event %description in text form %Vast speed improvement as save was the previous longest task %------------------------------------------------------------------------- %------------------------------------------------------------------------- %Grant Lohsen - v1.4 2009-04-15 %Georgia Tech Research Institute %Reads file header info and stores in the Root Structure. %------------------------------------------------------------------------- %------------------------------------------------------------------------- %Robert Seltzer - v1.5 2010-07-14 %BorgWarner Morse TEC %-Tested in MATLAB 2007b and 2010a. %-APPEARS to now be compatible with TDMS version 1.1 (a.k.a 4712) files; % although, this has not been extensively tested. For some unknown % reason, the version 1.2 (4713) files process noticeably faster. I think % that it may be related to the 'TDSm' tag. %-"Time Stamp" data type was not tested. %-"Waveform" fields was not tested. %-Fixed an error in the 'LV2MatlabDataType' function where LabView data type % 'tdsTypeSingleFloat' was defined as MATLAB data type 'float64' . Changed % to 'float32'. %-Added error trapping. %-Added feature to count the number of segments for pre-allocation as % opposed to estimating the number of segments. %-Added option to save the data in a MAT file. %-Fixed "invalid field name" error caused by excessive string lengths. %------------------------------------------------------------------------- %------------------------------------------------------------------------- %Robert Seltzer - v1.6 2010-09-01 %BorgWarner Morse TEC %-Tested in MATLAB 2010a. %-Fixed the "Coversion to cell from char is not possible" error found % by Francisco Botero in version 1.5. %-Added capability to process both fragmented or defragmented data. %-Fixed the "field" error found by Lawrence. %------------------------------------------------------------------------- %------------------------------------------------------------------------- %Christian Buxel - V1.7 2010-09-17 %RWTH Aachen %-Tested in Matlab2007b. %-Added support for german umlauts (�,�,�,�,�,�,�) in 'propsName' %------------------------------------------------------------------------- %------------------------------------------------------------------------- %Andr� R�egg - V1.7 2010-09-29 %Supercomputing Systems AG %-Tested in MATLAB 2006a & 2010b %-Make sure that data can be loaded correctly independently of character % encoding set in matlab. %-Fixed error if object consists of several segments with identical segment % information (if rawdataindex==0, not all segments were loaded) %------------------------------------------------------------------------- %------------------------------------------------------------------------- %Robert Seltzer - v1.7 2010-09-30 %BorgWarner Morse TEC %-Tested in MATLAB 2010b. %-Added 'error trapping' to the 'fixcharformatlab' function for group and % channel names that contain characters that are not 'A' through 'Z', % 'a' through 'z', 0 through 9 or underscore. The characters are replaced % with an underscore and a message is written to the Command Window % explaining to the user what has happened and how to fix it. Only tested % with a very limited number of "special" characters. %------------------------------------------------------------------------- %------------------------------------------------------------------------- %Robert Seltzer - v1.8 2010-10-12 %BorgWarner Morse TEC %-As a result of an error found by Peter Sulcs when loading data with very % long channel names, I have re-written the sections of the function that % creates the channel and property names that are used within the body of % the function to make them robust against long strings and strings % containing non-UTF8 characters. The original channel and property % names (no truncation or character replacement) are now retained and % included in the output structure. In order to implement this improvement, % I added a 'Property' field as a structure to the 'ConvertedData' output % structure. %-Added a more detailed 'help' description ('doc convertTDMS') of the % returned structure. %-List of channel names added as an output parameter of the function. %-Corrected an error in the time stamp converion. It was off by exactly % 1 hour. %-Tested in MATLAB 2010b with a limited number of diverse TDMS files. %------------------------------------------------------------------------- %------------------------------------------------------------------------- %Robert Seltzer - v1.8 2010-10-19 %BorgWarner Morse TEC %-Fixed an error found by Terenzio Girotto with the 'save' routine. %------------------------------------------------------------------------- %------------------------------------------------------------------------- %Robert Seltzer - v1.8 2010-10-25 %BorgWarner Morse TEC %-Fixed an error with channels that contain no data. Previously, if a % channel contained no data, then it was not passed to the output structure % even if it did contain properties. %-Added 'GroupNames' as an optional output variable. %-Fixed an error with capturing the properties of the Root object %------------------------------------------------------------------------- %------------------------------------------------------------------------- %Philip Top - v1.9 2010-11-09 %John Breneman %-restructured code as function calls %-seperated metadata reads from data reads %-preallocated space for SegInfo with two pass file read %-preallocated index information and defined segdataoffset for each segment %-preallocate space for data for speedup in case of fragmented files %-used matlab input parser instead of nargin switch %-vectorized timestamp reads for substantial speedup %------------------------------------------------------------------------- %------------------------------------------------------------------------- %Robert Seltzer - v1.9 2010-11-10 %BorgWarner Morse TEC %-Fixed an error error in the 'offset' calculation for strings %------------------------------------------------------------------------- %------------------------------------------------------------------------- %Philip Top - v1.95 2011-5-10 %Fix Bug with out of order file segments %Fix some issues with string array reads for newer version files, %------------------------------------------------------------------------- %------------------------------------------------------------------------- %Brad Humphreys - v1.96 2013-10-9 %Fixed problem with error catch messages themselves (interleaved, version, %big endian) using the TDMSFileName variable which was not available %(passed into) the getSegInfo function. Added ability to work with %interleaved and big endian files. %So function now covers: % -v1.0-v2.0 % -Interleaved and Decimated Data Formats % -Big and Little Endian storage % Does not work with DAQmxRawData %------------------------------------------------------------------------- %------------------------------------------------------------------------- %Brad Humphreys - v1.97 2013-11-13 %Added information to help documentation on how to deal with DAQmxRaw Data. %------------------------------------------------------------------------- %------------------------------------------------------------------------- %Brad Humphreys - v1.98 2014-5-27 %Per G. Lohsen's suggestion, added check to verify that first caharters are % TDMs. If not, errors out and lets user know that the selected file is % not a TDMS file. %------------------------------------------------------------------------- %------------------------------------------------------------------------- %Sebastian Schwarzendahl (alias Haench) - v1.99 2014-10-23 % Added support for complex data types % CSG - tdsTypeComplexSingleFloat=0x08000c % CDB - tdsTypeComplexDoubleFloat=0x10000d) % This feature was added in LV2013 (I believ) and produced an error in % the previous version of this code. %------------------------------------------------------------------------- %Initialize outputs ConvertVer='1.98'; %Version number of this conversion function ConvertedData=[]; p=inputParser(); p.addRequired('SaveConvertedFile',@(x) islogical(x)||(ismember(x,[0,1]))); p.addOptional('filename','',@(x) iscell(x)||exist(x,'file')); p.parse(varargin{:}); filename=p.Results.filename; SaveConvertedFile=p.Results.SaveConvertedFile; if isempty(filename) %Prompt the user for the file [filename,pathname]=uigetfile({'*.tdms','All Files (*.tdms)'},'Choose a TDMS File'); if filename==0 return end filename=fullfile(pathname,filename); end if iscell(filename) %For a list of files infilename=filename; else infilename=cellstr(filename); end for fnum=1:numel(infilename) if ~exist(infilename{fnum},'file') e=errordlg(sprintf('File ''%s'' not found.',infilename{fnum}),'File Not Found'); uiwait(e) return end FileNameLong=infilename{fnum}; [pathstr,name,ext]=fileparts(FileNameLong); FileNameShort=sprintf('%s%s',name,ext); FileNameNoExt=name; FileFolder=pathstr; if fnum==1 fprintf('\n\n') end fprintf('Converting ''%s''...',FileNameShort) fid=fopen(FileNameLong); if fid==-1 e=errordlg(sprintf('Could not open ''%s''.',FileNameLong),'File Cannot Be Opened'); uiwait(e) fprintf('\n\n') return end % Build a database with segment info [SegInfo,NumOfSeg]=getSegInfo(fid); % Build a database with channel info [channelinfo SegInfo]=getChannelInfo(fid,SegInfo,NumOfSeg); % Add channel count to SegInfo %SegInfo=addChanCount(SegInfo,channelinfo); % Get the raw data ob=getData(fid,channelinfo,SegInfo); %Returns the objects which have data. See postProcess function (appends to all of the objects) fclose(fid); %Assign the outputs ConvertedData(fnum).FileName=FileNameShort; ConvertedData(fnum).FileFolder=FileFolder; ConvertedData(fnum).SegTDMSVerNum=SegInfo.vernum; ConvertedData(fnum).NumOfSegments=NumOfSeg; [ConvertedData(fnum).Data,CurrGroupNames]=postProcess(ob,channelinfo); GroupNames(fnum)={CurrGroupNames}; TempChanNames={ConvertedData(fnum).Data.MeasuredData.Name}; TempChanNames(strcmpi(TempChanNames,'Root'))=[]; ChanNames(fnum)={sort(setdiff(TempChanNames',CurrGroupNames))}; if SaveConvertedFile MATFileNameShort=sprintf('%s.mat',FileNameNoExt); MATFileNameLong=fullfile(FileFolder,MATFileNameShort); try save(MATFileNameLong,'ConvertedData','ConvertVer','ChanNames') fprintf('\n\nConversion complete (saved in ''%s'').\n\n',MATFileNameShort) catch exception fprintf('\n\nConversion complete (could not save ''%s'').\n\t%s: %s\n\n',MATFileNameShort,exception.identifier,... exception.message) end else fprintf('\n\nConversion complete.\n\n') end end ci=channelinfo; end %% function [SegInfo,NumOfSeg]=getSegInfo(fid) % Go through the whole file and build a database of the segment lead-in % information: %1) Count the total number of segments in the file. (NumOfSeg) %2) Get all of the "Lead In" Headers from all of the segements and store % that in the SegInfo variable: % SegInfo.SegStartPosn: Absolute Starting Position of Segment in file % SegInfo.MetaStartPosn: Absolute Starting Position of Meta Data in file % SegInfo.DataStartPosn: Absolute Starting Position of Raw Data in file % SegInfo.DataLength: Number of bytes of Data in segment % SegInfo.vernum: LV version number (4712 is v1.0, 4713 is v2.0) % SegInfo.NumChan: number of channels in this segement. This is % only instatated in this function to 0. The addChanCount function % updates this to the actual value later. % % SegInfo.SegHasMetaData: There is Meta Data in the Segement % SegInfo.SegHasRawData: There is Raw Data in the Segment % SegInfo.SegHasDaqMxRaw: There is DAQmxRaw Data in the Segment % SegInfo.SegInterleaved: 0: Contigous Data, 1: Interleaved Data % SegInfo.SegBigEndian: 0: Little, 1:Big (numeric, leadin, raw, and meta) % SegInfo.SegHasNewObjList: New object list in Segment % %3) While doing the above, also include error trapping for incompatibity %Find the end of the file fseek(fid,0,'eof'); eoff=ftell(fid); frewind(fid); segCnt=0; CurrPosn=0; LeadInByteCount=28; %From the National Instruments web page (http://zone.ni.com/devzone/cda/tut/p/id/5696) under %the 'Lead In' description on page 2: Counted the bytes shown in the table. while (ftell(fid) ~= eoff) Ttag=fread(fid,1,'uint8','l'); Dtag=fread(fid,1,'uint8','l'); Stag=fread(fid,1,'uint8','l'); mtag=fread(fid,1,'uint8','l'); if Ttag==84 && Dtag==68 && Stag==83 && mtag==109 %Apparently, this sequence of numbers identifies the start of a new segment. segCnt=segCnt+1; %ToC Field ToC=fread(fid,1,'uint32','l'); kTocBigEndian=bitget(ToC,7); if kTocBigEndian kTocEndian='b'; else kTocEndian='l'; end %TDMS format version number vernum=fread(fid,1,'uint32',kTocEndian); %From the National Instruments web page (http://zone.ni.com/devzone/cda/tut/p/id/5696) under the 'Lead In' %description on page 2: %The next eight bytes (64-bit unsigned integer) describe the length of the remaining segment (overall length of the %segment minus length of the lead in). If further segments are appended to the file, this number can be used to %locate the starting point of the following segment. If an application encountered a severe problem while writing %to a TDMS file (crash, power outage), all bytes of this integer can be 0xFF. This can only happen to the last %segment in a file. nlen=fread(fid,1,'uint64',kTocEndian); if (nlen>2^63) break; else segLength=nlen; end TotalLength=segLength+LeadInByteCount; CurrPosn=CurrPosn+TotalLength; status=fseek(fid,CurrPosn,'bof'); %Move to the beginning position of the next segment if (status2^63) fseek(fid,0,'eof'); flen=ftell(fid); segLength=flen-LeadInByteCount-TotalLength; TotalLength=segLength+LeadInByteCount; else TotalLength=segLength+LeadInByteCount; CurrPosn=CurrPosn+TotalLength; fseek(fid,CurrPosn,'bof'); %Move to the beginning position of the next segment end SegInfo.SegStartPosn(segCnt)=StartPosn; SegInfo.MetaStartPosn(segCnt)=StartPosn+LeadInByteCount; SegInfo.DataStartPosn(segCnt)=SegInfo.MetaStartPosn(segCnt)+metaLength; SegInfo.DataLength(segCnt)=segLength-metaLength; SegInfo.vernum(segCnt)=vernum; SegInfo.SegHasMetaData(segCnt)=kTocMetaData; SegInfo.SegHasRawData(segCnt)=kTocRawData; SegInfo.SegHasDaqMxRaw(segCnt)=kTocDaqMxRawData; SegInfo.SegInterleaved(segCnt)=kTocInterleavedData; SegInfo.SegBigEndian(segCnt)=kTocBigEndian; SegInfo.SegHasNewObjList(segCnt)=kTocNewObjectList; end end NumOfSeg=segCnt; end %% function [index SegInfo]=getChannelInfo(fid,SegInfo,NumOfSeg) %Loop through the segements and get all of the object information. % name: Short name such as Object3 % long_name: The path directory form of the object name % rawdatacount: number of segments that have raw data for this object % datastartindex: Absolute position in file of the begining of all raw % data for the segment. Add the rawdata offset to get to this object's % start point. % arrayDim: array dimension of raw data (for now should always be one as % NI does not support multi-dimension % nValues: number of values of raw data in segment % byteSize: number of bytes for string/char data % index: segment numbers that contain raw data for this object % rawdataoffset: offset from datastartindex for the begining of raw data % for this object % multiplier: As part of LV's optimization, it will appenend raw data writes to the % TDMS file when meta data does not change, but forget to update the % nValues. multiplier is the number of times the base pattern % repeats. % skip: When multiple raw data writes are appended, this is the number of % bytes between the end of a channel in one block to the begining of that channel in % the next block. % dataType: The data type for this object. See LV2MatlabDataType function. % datasize: Number of bytes for each raw data entry. % Property Structures: Structure containg properties. %Initialize variables for the file conversion index=struct(); objOrderList={}; for segCnt=1:NumOfSeg %Go to the segment starting position of the segment fseek(fid,SegInfo.SegStartPosn(segCnt)+28,'bof'); % +28 bytes: TDSm (4) + Toc (4) + segVer# (4) + segLength (8) + metaLength (8) %segVersionNum=SegInfo.vernum(segCnt); kTocMetaData=SegInfo.SegHasMetaData(segCnt); kTocNewObjectList=SegInfo.SegHasNewObjList(segCnt); kTocRawData=SegInfo.SegHasRawData(segCnt); kTocInterleavedData=SegInfo.SegInterleaved(segCnt); kTocBigEndian=SegInfo.SegBigEndian(segCnt); if kTocBigEndian kTocEndian='b'; else kTocEndian='l'; end %% Process Meta Data %If the object list from the last segment should be used.... if (kTocNewObjectList==0) fnm=fieldnames(index); %Get a list of the objects/channels to loop through for kk=1:length(fnm) ccnt=index.(fnm{kk}).rawdatacount; if (ccnt>0) %If there is raw data info in the previous segement, copy it into the new segement if (index.(fnm{kk}).index(ccnt)==segCnt-1) ccnt=ccnt+1; index.(fnm{kk}).rawdatacount=ccnt; index.(fnm{kk}).datastartindex(ccnt)=SegInfo.DataStartPosn(segCnt); index.(fnm{kk}).arrayDim(ccnt)=index.(fnm{kk}).arrayDim(ccnt-1); index.(fnm{kk}).nValues(ccnt)=index.(fnm{kk}).nValues(ccnt-1); index.(fnm{kk}).byteSize(ccnt)=index.(fnm{kk}).byteSize(ccnt-1); index.(fnm{kk}).index(ccnt)=segCnt; index.(fnm{kk}).rawdataoffset(ccnt)=index.(fnm{kk}).rawdataoffset(ccnt-1); SegInfo.NumChan(segCnt)=SegInfo.NumChan(segCnt)+1; end end end end %If there is Meta data in the segement if kTocMetaData numObjInSeg=fread(fid,1,'uint32',kTocEndian); if (kTocNewObjectList) objOrderList=cell(numObjInSeg,1); end for q=1:numObjInSeg obLength=fread(fid,1,'uint32',kTocEndian); %Get the length of the objects name ObjName=convertToText(fread(fid,obLength,'uint8','l'))'; %Get the objects name if strcmp(ObjName,'/') long_obname='Root'; else long_obname=ObjName; %Delete any apostrophes. If the first character is a slash (forward or backward), delete it too. long_obname(strfind(long_obname,''''))=[]; if strcmpi(long_obname(1),'/') || strcmpi(long_obname(1),'\') long_obname(1)=[]; end end newob=0; %Create object's name. Use a generic field name to avoid issues with strings that are too long and/or %characters that cannot be used in MATLAB variable names. The actual channel name is retained for the final %output structure. if exist('ObjNameList','var') %Check to see if the object already exists NameIndex=find(strcmpi({ObjNameList.LongName},long_obname)==1,1,'first'); if isempty(NameIndex) newob=1; %It does not exist, so create the generic name field name ObjNameList(end+1).FieldName=sprintf('Object%.0f',numel(ObjNameList)+1); ObjNameList(end).LongName=long_obname; NameIndex=numel(ObjNameList); end else %No objects exist, so create the first one using a generic name field name. ObjNameList.FieldName='Object1'; ObjNameList.LongName=long_obname; NameIndex=1; newob=1; end %Assign the generic field name obname=ObjNameList(NameIndex).FieldName; %Create the 'index' structure if (~isfield(index,obname)) index.(obname).name=obname; index.(obname).long_name=long_obname; index.(obname).rawdatacount=0; index.(obname).datastartindex=zeros(NumOfSeg,1); index.(obname).arrayDim=zeros(NumOfSeg,1); index.(obname).nValues=zeros(NumOfSeg,1); index.(obname).byteSize=zeros(NumOfSeg,1); index.(obname).index=zeros(NumOfSeg,1); index.(obname).rawdataoffset=zeros(NumOfSeg,1); index.(obname).multiplier=ones(NumOfSeg,1); index.(obname).skip=zeros(NumOfSeg,1); end if (kTocNewObjectList) objOrderList{q}=obname; else if ~ismember(obname,objOrderList) objOrderList{end+1}=obname; end end %Get the raw data Index rawdataindex=fread(fid,1,'uint32',kTocEndian); if rawdataindex==0 if segCnt==0 e=errordlg(sprintf('Seqment %.0f within ''%s'' has ''rawdataindex'' value of 0 (%s.m).',segCnt,... TDMSFileNameShort,mfilename),'Incorrect ''rawdataindex'''); uiwait(e) end if kTocRawData if (kTocNewObjectList) ccnt=index.(obname).rawdatacount+1; else ccnt=index.(obname).rawdatacount; end index.(obname).rawdatacount=ccnt; index.(obname).datastartindex(ccnt)=SegInfo.DataStartPosn(segCnt); index.(obname).arrayDim(ccnt)=index.(obname).arrayDim(ccnt-1); index.(obname).nValues(ccnt)=index.(obname).nValues(ccnt-1); index.(obname).byteSize(ccnt)=index.(obname).byteSize(ccnt-1); index.(obname).index(ccnt)=segCnt; SegInfo.NumChan(segCnt)=SegInfo.NumChan(segCnt)+1; end elseif rawdataindex+1==2^32 %Objects raw data index matches previous index - no changes. The root object will always have an %'FFFFFFFF' entry if strcmpi(index.(obname).long_name,'Root') index.(obname).rawdataindex=0; else %Need to account for the case where an object (besides the 'root') is added that has no data but reports %using previous. if newob index.(obname).rawdataindex=0; else if kTocRawData if (kTocNewObjectList) ccnt=index.(obname).rawdatacount+1; else ccnt=index.(obname).rawdatacount; end index.(obname).rawdatacount=ccnt; index.(obname).datastartindex(ccnt)=SegInfo.DataStartPosn(segCnt); index.(obname).arrayDim(ccnt)=index.(obname).arrayDim(ccnt-1); index.(obname).nValues(ccnt)=index.(obname).nValues(ccnt-1); index.(obname).byteSize(ccnt)=index.(obname).byteSize(ccnt-1); index.(obname).index(ccnt)=segCnt; SegInfo.NumChan(segCnt)=SegInfo.NumChan(segCnt)+1; end end end else %Get new object information if (kTocNewObjectList) ccnt=index.(obname).rawdatacount+1; else ccnt=index.(obname).rawdatacount; if (ccnt==0) ccnt=1; end end index.(obname).rawdatacount=ccnt; index.(obname).datastartindex(ccnt)=SegInfo.DataStartPosn(segCnt); %index(end).lenOfIndexInfo=fread(fid,1,'uint32'); index.(obname).dataType=fread(fid,1,'uint32',kTocEndian); if (index.(obname).dataType~=32) index.(obname).datasize=getDataSize(index.(obname).dataType); end index.(obname).arrayDim(ccnt)=fread(fid,1,'uint32',kTocEndian); index.(obname).nValues(ccnt)=fread(fid,1,'uint64',kTocEndian); index.(obname).index(ccnt)=segCnt; SegInfo.NumChan(segCnt)=SegInfo.NumChan(segCnt)+1; if index.(obname).dataType==32 %Datatype is a string index.(obname).byteSize(ccnt)=fread(fid,1,'uint64',kTocEndian); else index.(obname).byteSize(ccnt)=0; end end %Get the properties numProps=fread(fid,1,'uint32',kTocEndian); if numProps>0 if isfield(index.(obname),'PropertyInfo') PropertyInfo=index.(obname).PropertyInfo; else clear PropertyInfo end for p=1:numProps propNameLength=fread(fid,1,'uint32',kTocEndian); switch 1 case 1 PropName=fread(fid,propNameLength,'*uint8','l')'; PropName=native2unicode(PropName,'UTF-8'); case 2 PropName=fread(fid,propNameLength,'uint8=>char','l')'; otherwise end propsDataType=fread(fid,1,'uint32',kTocEndian); %Create property's name. Use a generic field name to avoid issues with strings that are too long and/or %characters that cannot be used in MATLAB variable names. The actual property name is retained for the %final output structure. if exist('PropertyInfo','var') %Check to see if the property already exists for this object. Need to get the existing 'PropertyInfo' %structure for this object. The 'PropertyInfo' structure is not necessarily the same for every %object in the data file. PropIndex=find(strcmpi({PropertyInfo.Name},PropName)); if isempty(PropIndex) %Is does not exist, so create the generic name field name propExists=false; PropIndex=numel(PropertyInfo)+1; propsName=sprintf('Property%.0f',PropIndex); PropertyInfo(PropIndex).Name=PropName; PropertyInfo(PropIndex).FieldName=propsName; else %Assign the generic field name propExists=true; propsName=PropertyInfo(PropIndex).FieldName; end else %No properties exist for this object, so create the first one using a generic name field name. propExists=false; PropIndex=p; propsName=sprintf('Property%.0f',PropIndex); PropertyInfo(PropIndex).Name=PropName; PropertyInfo(PropIndex).FieldName=propsName; end dataExists=isfield(index.(obname),'data'); %Get the number of samples already found and in the object if dataExists nsamps=index.(obname).nsamples+1; else nsamps=0; end if propsDataType==32 %String data type PropertyInfo(PropIndex).DataType='String'; propsValueLength=fread(fid,1,'uint32',kTocEndian); propsValue=convertToText(fread(fid,propsValueLength,'uint8=>char',kTocEndian))'; if propExists if isfield(index.(obname).(propsName),'cnt') cnt=index.(obname).(propsName).cnt+1; else cnt=1; end index.(obname).(propsName).cnt=cnt; index.(obname).(propsName).value{cnt}=propsValue; index.(obname).(propsName).samples(cnt)=nsamps; else if strcmp(index.(obname).long_name,'Root') %Header data index.(obname).(propsName).name=index.(obname).long_name; index.(obname).(propsName).value={propsValue}; index.(obname).(propsName).cnt=1; else index.(obname).(propsName).name=PropertyInfo(PropIndex).Name; index.(obname).(propsName).datatype=PropertyInfo(PropIndex).DataType; index.(obname).(propsName).cnt=1; index.(obname).(propsName).value=cell(nsamps,1); %Pre-allocation index.(obname).(propsName).samples=zeros(nsamps,1); %Pre-allocation if iscell(propsValue) index.(obname).(propsName).value(1)=propsValue; else index.(obname).(propsName).value(1)={propsValue}; end index.(obname).(propsName).samples(1)=nsamps; end end else %Numeric data type if propsDataType==68 %Timestamp PropertyInfo(PropIndex).DataType='Time'; %Timestamp data type if kTocBigEndian tsec=fread(fid,1,'uint64','b')+fread(fid,1,'uint64','b')/2^64; %time since Jan-1-1904 in seconds else tsec=fread(fid,1,'uint64','l')/2^64+fread(fid,1,'uint64','l'); %time since Jan-1-1904 in seconds end %tsec=fread(fid,1,'uint64',kTocEndian)/2^64+fread(fid,1,'uint64',kTocEndian); %time since Jan-1-1904 in seconds %R. Seltzer: Not sure why '5/24' (5 hours) is subtracted from the time value. That's how it was %coded in the original function I downloaded from MATLAB Central. But I found it to be 1 hour too %much. So, I changed it to '4/24'. %propsValue=tsec/86400+695422-5/24; %/864000 convert to days; +695422 days from Jan-0-0000 to Jan-1-1904 propsValue=tsec/86400+695422-4/24; %/864000 convert to days; +695422 days from Jan-0-0000 to Jan-1-1904 else %Numeric PropertyInfo(PropIndex).DataType='Numeric'; matType=LV2MatlabDataType(propsDataType); if strcmp(matType,'Undefined') e=errordlg(sprintf('No MATLAB data type defined for a ''Property Data Type'' value of ''%.0f''.',... propsDataType),'Undefined Property Data Type'); uiwait(e) fclose(fid); return end if strcmp(matType,'uint8=>char') propsValue=convertToText(fread(fid,1,'uint8',kTocEndian)); else propsValue=fread(fid,1,matType,kTocEndian); end end if propExists cnt=index.(obname).(propsName).cnt+1; index.(obname).(propsName).cnt=cnt; index.(obname).(propsName).value(cnt)=propsValue; index.(obname).(propsName).samples(cnt)=nsamps; else index.(obname).(propsName).name=PropertyInfo(PropIndex).Name; index.(obname).(propsName).datatype=PropertyInfo(PropIndex).DataType; index.(obname).(propsName).cnt=1; index.(obname).(propsName).value=NaN(nsamps,1); %Pre-allocation index.(obname).(propsName).samples=zeros(nsamps,1); %Pre-allocation index.(obname).(propsName).value(1)=propsValue; index.(obname).(propsName).samples(1)=nsamps; end end end %'end' for the 'Property' loop index.(obname).PropertyInfo=PropertyInfo; end end %'end' for the 'Objects' loop end %Address Decimation and Interleaving if (kTocRawData) singleSegDataSize=0; rawDataBytes=0; chanBytes=0; rawDataOffset=0; for kk=1:numel(objOrderList) obname=objOrderList{kk}; ccnt=hasRawDataInSeg(segCnt,index.(obname)); if ccnt>0 %If segement has raw data index.(obname).rawdataoffset(ccnt)=rawDataOffset; if index.(obname).dataType==32 %Datatype is a string rawDataBytes=index.(obname).byteSize(ccnt); else rawDataBytes=index.(obname).nValues(ccnt)*index.(obname).datasize; if kTocInterleavedData chanBytes=index.(obname).datasize; else chanBytes=rawDataBytes; end end rawDataOffset=rawDataOffset+chanBytes; singleSegDataSize=singleSegDataSize+rawDataBytes; end end %Calculate the offset for each segment. The offset is the amount %of bytes for one non appened (non optimized) segement % for kk=1:numel(objOrderList) % obname=objOrderList{kk}; % if hasRawDataInSeg(segCnt,index.(obname)) %If segement has raw data % index.(obname).rawdataoffset(ccnt)=rawdataoffset; %This will set the offset correctly for dec data. % % %Update the singleSegDataSize value for the next object % if index.(obname).dataType==32 %Datatype is a string % singleSegDataSize=singleSegDataSize+index.(obname).byteSize(ccnt); % else % singleSegDataSize=singleSegDataSize+index.(obname).nValues(ccnt)*index.(obname).datasize; % end % % if kTocInterleavedData % rawdataoffset= % else % rawdataoffset=singleSegDataSize; % end % end % end % %As part of LV's optimization, it will append back to back %segements into one segement (when nothing changes in meta data). if (singleSegDataSize~=SegInfo.DataLength(segCnt)) %Multiple appended raw segement (offset did not grow to be greater than the DataLength) numAppSegs=floor(SegInfo.DataLength(segCnt)/singleSegDataSize); %Number of appened (optimized segements) for kk=1:numel(objOrderList) %Loop through all the objects obname=objOrderList{kk}; ccnt=hasRawDataInSeg(segCnt,index.(obname)); if ccnt>0 %If segement has raw data if kTocInterleavedData if index.(obname).dataType==32 %Datatype is a string error('Interleaved string channels are not supported.') end index.(obname).multiplier(ccnt)=numAppSegs*index.(obname).nValues(ccnt); index.(obname).skip(ccnt)=singleSegDataSize/index.(obname).nValues(ccnt)-index.(obname).datasize; index.(obname).nValues(ccnt)=1; else %Decimated Data index.(obname).multiplier(ccnt)=numAppSegs; if index.(obname).dataType==32 %Datatype is a string index.(obname).skip(ccnt)=singleSegDataSize-index.(obname).byteSize(ccnt); else index.(obname).skip(ccnt)=singleSegDataSize-index.(obname).nValues(ccnt)*index.(obname).datasize; end end end end else %If single segement if kTocInterleavedData for kk=1:numel(objOrderList) %Loop through all the objects obname=objOrderList{kk}; ccnt=hasRawDataInSeg(segCnt,index.(obname)); if ccnt>0 %If segement has raw data if (index.(obname).index(ccnt)==segCnt) %If the object has rawdata in the current segment if index.(obname).dataType==32 %Datatype is a string error('Interleaved string channels are not supported.') end index.(obname).multiplier(ccnt)=index.(obname).nValues(ccnt); index.(obname).skip(ccnt)=singleSegDataSize/index.(obname).nValues(ccnt)-index.(obname).datasize; index.(obname).nValues(ccnt)=1; end end end end end % if (offset~=SegInfo.DataLength(segCnt)) %If the offset grew to be greater than the DataLength (from Meta Data) % multiplier=floor(SegInfo.DataLength(segCnt)/offset); % for kk=1:numel(objOrderList) %Loop through the objects % obname=objOrderList{kk}; % ccnt=index.(obname).rawdatacount; % if ccnt>0 % index.(obname).multiplier(segCnt)=multiplier; % if index.(obname).dataType==32 %Datatype is a string % index.(obname).skip(ccnt)=offset-index.(obname).byteSize(ccnt); % else % index.(obname).skip(ccnt)=offset-index.(obname).nValues(ccnt)*index.(obname).datasize; % end % end % end % % end aa=1; % %Now adjust multiplier and skip if the data is interleaved % if kTocInterleavedData % if muliplier>0 %Muliple raw data segements appened % % else %Single raw segement % % end % end end end %clean up the index if it has to much data fnm=fieldnames(index); for kk=1:numel(fnm) ccnt=index.(fnm{kk}).rawdatacount+1; index.(fnm{kk}).datastartindex(ccnt:end)=[]; index.(fnm{kk}).arrayDim(ccnt:end)=[]; index.(fnm{kk}).nValues(ccnt:end)=[]; index.(fnm{kk}).byteSize(ccnt:end)=[]; index.(fnm{kk}).index(ccnt:end)=[]; index.(fnm{kk}).rawdataoffset(ccnt:end)=[]; index.(fnm{kk}).multiplier(ccnt:end)=[]; index.(fnm{kk}).skip(ccnt:end)=[]; end end %% function ob=getData(fid,index,SegInfo) %Using the file id (fid) and the index database, get the raw data. %Returns the data in the ob structure. The fields in the structure are the %generic object names: % ob.Object3.data - raw data % ob.Object3.nsamples - number of samples % % Note that not all of the objects in the index may not be in ob as it % contains only those objects which have raw data. ob=[]; fnm=fieldnames(index); %Get the object names for kk=1:length(fnm) %Loop through objects id=index.(fnm{kk}); nsamples=sum(id.nValues.*id.multiplier); if id.rawdatacount>0 %Only work with channels with raw data cname=id.name; ob.(cname).nsamples=0; %Initialize the data matrix if id.dataType==32 ob.(cname).data=cell(nsamples,1); else ob.(cname).data=zeros(nsamples,1); end for rr=1:id.rawdatacount %Loop through each of the segments and read the raw data %Move to the raw data start position fseek(fid,id.datastartindex(rr)+id.rawdataoffset(rr),'bof'); nvals=id.nValues(rr); segmentNum=index.(cname).index(rr); segInterleaved=SegInfo.SegInterleaved(segmentNum); if SegInfo.SegBigEndian(segmentNum) kTocEndian='b'; else kTocEndian='l'; end numChan=SegInfo.NumChan(segmentNum); if nvals>0 %If there is data in this segement switch id.dataType case 32 %String %From the National Instruments web page (http://zone.ni.com/devzone/cda/tut/p/id/5696) under the %'Raw Data' description on page 4: %String type channels are preprocessed for fast random access. All strings are concatenated to a %contiguous piece of memory. The offset of the first character of each string in this contiguous piece %of memory is stored to an array of unsigned 32-bit integers. This array of offset values is stored %first, followed by the concatenated string values. This layout allows client applications to access %any string value from anywhere in the file by repositioning the file pointer a maximum of three times %and without reading any data that is not needed by the client. data=cell(1,nvals*id.multiplier(rr)); %Pre-allocation for mm=1:id.multiplier(rr) StrOffsetArray=fread(fid,nvals,'uint32','l'); for dcnt=1:nvals if dcnt==1 StrLength=StrOffsetArray(dcnt); else StrLength=StrOffsetArray(dcnt)-StrOffsetArray(dcnt-1); end data{1,dcnt+(mm-1)*nvals}=char(convertToText(fread(fid,StrLength,'uint8=>char','l'))'); end if (id.multiplier(rr)>1)&&(id.skip(rr)>0) fseek(fid,id.skip(rr),'cof'); end end cnt=nvals*id.multiplier(rr); case 68 %Timestamp %data=NaN(1,nvals); %Pre-allocation data=NaN(1,nvals*id.multiplier(rr)); for mm=1:id.multiplier(rr) dn=fread(fid,2*nvals,'uint64',kTocEndian); tsec=dn(1:2:end)/2^64+dn(2:2:end); data((mm-1)*nvals+1:(mm)*nvals)=tsec/86400+695422-4/24; fseek(fid,id.skip(rr),'cof'); end %{ for dcnt=1:nvals tsec=fread(fid,1,'uint64')/2^64+fread(fid,1,'uint64'); %time since Jan-1-1904 in seconds %R. Seltzer: Not sure why '5/24' (5 hours) is subtracted from the time value. That's how it was %coded in the original function I downloaded from MATLAB Central. But I found it to be 1 hour too %much. So, I changed it to '4/24'. data(1,dcnt)=tsec/86400+695422-5/24; %/864000 convert to days; +695422 days from Jan-0-0000 to Jan-1-1904 data(1,dcnt)=tsec/86400+695422-4/24; %/864000 convert to days; +695422 days from Jan-0-0000 to Jan-1-1904 end %} cnt=nvals*id.multiplier(rr); otherwise %Numeric matType=LV2MatlabDataType(id.dataType); if strcmp(matType,'Undefined') %Bad Data types catch e=errordlg(sprintf('No MATLAB data type defined for a ''Raw Data Type'' value of ''%.0f''.',... id.dataType),'Undefined Raw Data Type'); uiwait(e) fclose(fid); return end if (id.skip(rr)>0) ntype=sprintf('%d*%s',nvals,matType); [data,cnt]=fread(fid,nvals*id.multiplier(rr),ntype,id.skip(rr),kTocEndian); if strcmp(matType,'uint8=>char') data=convertToText(data); end else % Added by Haench start if (id.dataType == 524300) || (id.dataType == 1048589) % complex CDB data [data,cnt]=fread(fid,2*nvals*id.multiplier(rr),matType,kTocEndian); data= data(1:2:end)+1i*data(2:2:end); cnt = cnt/2; else [data,cnt]=fread(fid,nvals*id.multiplier(rr),matType,kTocEndian); end % Haench end % Original: [data,cnt]=fread(fid,nvals*id.multiplier(rr),matType,kTocEndian); end end %Update the sample counter if isfield(ob.(cname),'nsamples') ssamples=ob.(cname).nsamples; else ssamples=0; end if (cnt>0) ob.(cname).data(ssamples+1:ssamples+cnt,1)=data; ob.(cname).nsamples=ssamples+cnt; end end end end end end function [DataStructure,GroupNames]=postProcess(ob,index) %Re-organize the 'ob' structure into a more user friendly format for output. DataStructure.Root=[]; DataStructure.MeasuredData.Name=[]; DataStructure.MeasuredData.Data=[]; obFieldNames=fieldnames(index); cntData=1; for i=1:numel(obFieldNames) cname=obFieldNames{i}; if strcmp(index.(cname).long_name,'Root') DataStructure.Root.Name=index.(cname).long_name; %Assign all the 'Property' values if isfield(index.(cname),'PropertyInfo') for p=1:numel(index.(cname).PropertyInfo) cfield=index.(cname).PropertyInfo(p).FieldName; if isfield(index.(cname).(cfield),'datatype') DataType=index.(cname).(cfield).datatype; else %ASSUME a 'string' data type DataType='String'; end DataStructure.Root.Property(p).Name=index.(cname).PropertyInfo(p).Name; switch DataType case 'String' if iscell(index.(cname).(cfield).value) Value=index.(cname).(cfield).value'; else Value=cellstr(index.(cname).(cfield).value); end case 'Time' clear Value if index.(cname).(cfield).cnt==1 if iscell(index.(cname).(cfield).value) Value=datestr(cell2mat(index.(cname).(cfield).value),'dd-mmm-yyyy HH:MM:SS'); else Value=datestr(index.(cname).(cfield).value,'dd-mmm-yyyy HH:MM:SS'); end else Value=cell(index.(cname).(cfield).cnt,1); for c=1:index.(cname).(cfield).cnt if iscell(index.(cname).(cfield).value) Value(c)={datestr(cell2mat(index.(cname).(cfield).value),'dd-mmm-yyyy HH:MM:SS')}; else Value(c)={datestr(index.(cname).(cfield).value,'dd-mmm-yyyy HH:MM:SS')}; end end end case 'Numeric' if isfield(index.(cname).(cfield),'cnt') Value=NaN(index.(cname).(cfield).cnt,1); else if iscell(index.(cname).(cfield).value) Value=NaN(numel(cell2mat(index.(cname).(cfield).value)),1); else Value=NaN(numel(index.(cname).(cfield).value),1); end end for c=1:numel(Value) if iscell(index.(cname).(cfield).value) Value(c)=index.(cname).(cfield).value{c}; else Value(c)=index.(cname).(cfield).value(c); end end otherwise e=errordlg(sprintf(['No format defined for Data Type ''%s'' in the private function ''postProcess'' '... 'within %s.m.'],index.(cname).(cfield).datatype,mfilename),'Undefined Property Data Type'); uiwait(e) return end if isempty(Value) DataStructure.Root.Property(p).Value=[]; else DataStructure.Root.Property(p).Value=Value; end end end end DataStructure.MeasuredData(cntData).Name=index.(cname).long_name; %Should only need the 'ShortName' for debugging the function %DataStructure.MeasuredData(cntData).ShortName=cname; if (isfield(ob,cname)) if isfield(ob.(cname),'data') DataStructure.MeasuredData(cntData).Data=ob.(cname).data; %The following field is redundant because the information can be obtained from the size of the 'Data' field. DataStructure.MeasuredData(cntData).Total_Samples=ob.(cname).nsamples; else DataStructure.MeasuredData(cntData).Data=[]; DataStructure.MeasuredData(cntData).Total_Samples=0; end else DataStructure.MeasuredData(cntData).Data=[]; DataStructure.MeasuredData(cntData).Total_Samples=0; end %Assign all the 'Property' values if isfield(index.(cname),'PropertyInfo') for p=1:numel(index.(cname).PropertyInfo) cfield=index.(cname).PropertyInfo(p).FieldName; DataStructure.MeasuredData(cntData).Property(p).Name=index.(cname).(cfield).name; if strcmpi(DataStructure.MeasuredData(cntData).Property(p).Name,'Root') Value=index.(cname).(cfield).value; else switch index.(cname).(cfield).datatype case 'String' clear Value if index.(cname).(cfield).cnt==1 if iscell(index.(cname).(cfield).value) Value=char(index.(cname).(cfield).value); else Value=index.(cname).(cfield).value; end else Value=cell(index.(cname).(cfield).cnt,1); for c=1:index.(cname).(cfield).cnt if iscell(index.(cname).(cfield).value) Value(c)=index.(cname).(cfield).value; else Value(c)={index.(cname).(cfield).value}; end end end case 'Time' clear Value if index.(cname).(cfield).cnt==1 if iscell(index.(cname).(cfield).value) Value=datestr(cell2mat(index.(cname).(cfield).value),'dd-mmm-yyyy HH:MM:SS'); else Value=datestr(index.(cname).(cfield).value,'dd-mmm-yyyy HH:MM:SS'); end else Value=cell(index.(cname).(cfield).cnt,1); for c=1:index.(cname).(cfield).cnt if iscell(index.(cname).(cfield).value) Value(c)={datestr(cell2mat(index.(cname).(cfield).value),'dd-mmm-yyyy HH:MM:SS')}; else Value(c)={datestr(index.(cname).(cfield).value,'dd-mmm-yyyy HH:MM:SS')}; end end end case 'Numeric' if isfield(index.(cname).(cfield),'cnt') Value=NaN(index.(cname).(cfield).cnt,1); else if iscell(index.(cname).(cfield).value) Value=NaN(numel(cell2mat(index.(cname).(cfield).value)),1); else Value=NaN(numel(index.(cname).(cfield).value),1); end end for c=1:numel(Value) if iscell(index.(cname).(cfield).value) Value(c)=index.(cname).(cfield).value{c}; else Value(c)=index.(cname).(cfield).value(c); end end otherwise e=errordlg(sprintf(['No format defined for Data Type ''%s'' in the private function ''postProcess'' '... 'within %s.m.'],index.(cname).(cfield).datatype,mfilename),'Undefined Property Data Type'); uiwait(e) return end end if isempty(Value) DataStructure.MeasuredData(cntData).Property(p).Value=[]; else DataStructure.MeasuredData(cntData).Property(p).Value=Value; end end else DataStructure.MeasuredData(cntData).Property=[]; end cntData = cntData + 1; end %'end' for the 'groups/channels' loop %Extract the Group names GroupIndices=false(numel(DataStructure.MeasuredData),1); for d=1:numel(DataStructure.MeasuredData) if ~strcmpi(DataStructure.MeasuredData(d).Name,'Root') if (DataStructure.MeasuredData(d).Total_Samples==0) fs=strfind(DataStructure.MeasuredData(d).Name,'/'); if (isempty(fs)) GroupIndices(d)=true; end end end end if any(GroupIndices) GroupNames=sort({DataStructure.MeasuredData(GroupIndices).Name})'; else GroupNames=[]; end end function SegInfo=addChanCount(SegInfo,channelinfo) %This function determines the number of channels in each segement by %looping through the channels in channelinfo. It looks at their index, %which conatins the segments numbers in which that channel has data. %It then increments the NumChan field in the segment info database %(SegInfo). fnm=fieldnames(channelinfo); %Get the channels for ccnt=1:numel(fnm) %Loop through the channels channelName=fnm{ccnt}; ind=channelinfo.(channelName).index; %Get the index for icnt=1:numel(ind) segNum=ind(icnt); SegInfo.NumChan(segNum)=SegInfo.NumChan(segNum)+1; end end end function sz=getDataSize(LVType) %Get the number of bytes for each LV data type. See LV2MatlabDataType. switch(LVType) case 0 sz=0; case {1,5,33} sz=1; case 68 sz=16; case {8,10} sz=8; case {3,7,9} sz=4; case {2,6} sz=2; case 32 e=errordlg('Do not call the getDataSize function for strings. Their size is written in the data file','Error'); uiwait(e) sz=NaN; case 11 sz=10; % Added by Haench for tdsTypeComplexSingleFloat=0x08000c,tdsTypeComplexDoubleFloat=0x10000d, case 524300 sz=8; case 1048589 sz=16; % end add haench otherwise error('LVData type %d is not defined',LVType) end end function matType=LV2MatlabDataType(LVType) %Cross Refernce Labview TDMS Data type to MATLAB switch LVType case 0 %tdsTypeVoid matType=''; case 1 %tdsTypeI8 matType='int8'; case 2 %tdsTypeI16 matType='int16'; case 3 %tdsTypeI32 matType='int32'; case 4 %tdsTypeI64 matType='int64'; case 5 %tdsTypeU8 matType='uint8'; case 6 %tdsTypeU16 matType='uint16'; case 7 %tdsTypeU32 matType='uint32'; case 8 %tdsTypeU64 matType='uint64'; case 9 %tdsTypeSingleFloat matType='single'; case 10 %tdsTypeDoubleFloat matType='double'; case 11 %tdsTypeExtendedFloat matType='10*char'; case 25 %tdsTypeSingleFloat with units matType='Undefined'; case 26 %tdsTypeDoubleFloat with units matType='Undefined'; case 27 %tdsTypeextendedFloat with units matType='Undefined'; case 32 %tdsTypeString matType='uint8=>char'; case 33 %tdsTypeBoolean matType='bit1'; case 68 %tdsTypeTimeStamp matType='2*int64'; % Added by Haench for tdsTypeComplexSingleFloat=0x08000c,tdsTypeComplexDoubleFloat=0x10000d, case 524300 matType='single'; case 1048589 matType='double'; % end add haench otherwise matType='Undefined'; end end function ccnt=hasRawDataInSeg(segCnt,ob) %Function to check and see if a channel has raw data in a particulr %segment. If it does, return the index of index (ccnt). If not present, %return ccn=0 indexList=ob.index; ccnt=find(indexList==segCnt); if isempty(ccnt) ccnt=0; else ccnt=ccnt(1); end end function text=convertToText(bytes) %Convert numeric bytes to the character encoding localy set in MATLAB (TDMS uses UTF-8) text=native2unicode(bytes,'UTF-8'); end e文件 classdef NicoletFile < handle % NICOLETFILE Reading Nicolet .e files. % % GETTING DATA % You can load data from the .e file using the GETDATA method. The % inputs to the method are the object, the segment of the data file % that you want to load data from, the min, and max index you want to % retrieve, and a vector of channels that you want to retrieve. % % Example: % OUT = GETDATA(OBJ, 1, [1 1000], 1:10) will return the first 1000 % values on the first 10 channels of the first segment of the file. % % GETTING Nr OF SAMPLES % Use the GETNRSAMPLES method to find the number of samples per % channel in each data segment. % % WARNING! % The .e format allows for changes in the TimeSeries map during the % recording. This results in multiple TSINFO structures. Depending on % where these structures are located in the .e file, the appropriate % TSINFO structure should be used. However, there seems to be a bug in % the Nicolet .e file writer which sometimes renders the TSINFO structures % unreadable on disk (verify with hex-edit). %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Copyright 2013 Trustees of the University of Pennsylvania % % Licensed under the Apache License, Version 2.0 (the "License"); % you may not use this file except in compliance with the License. % You may obtain a copy of the License at % % http://www.apache.org/licenses/LICENSE-2.0 % % Unless required by applicable law or agreed to in writing, software % distributed under the License is distributed on an "AS IS" BASIS, % WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. % See the License for the specific language governing permissions and % limitations under the License. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Joost Wagenaar, Jan 2015 % Cristian Donos, Dec 2015 % Jan Brogger, Jun 2016 % Callum Stewart, Apr 2018 properties fileName patientInfo segments eventMarkers end properties (Hidden) sections index sigInfo tsInfos chInfo notchFreq montage Qi Qii allIndexIDs useTSinfoIdx = 1 end properties (Constant, Hidden) LABELSIZE = 32 TSLABELSIZE = 64 UNITSIZE = 16 ITEMNAMESIZE = 64 end methods function obj = NicoletFile(filename) h = fopen(filename,'r','ieee-le'); [folder, name, ext] = fileparts(filename); assert(strcmp(ext,'.e'), 'File extention must be .e'); if isempty(folder) filename = fullfile(pwd,filename); end obj.fileName = filename; % Get init misc1 = fread(h,5, 'uint32'); %#ok unknown = fread(h,1,'uint32'); %#ok indexIdx = fread(h,1,'uint32'); % Get TAGS structure and Channel IDS fseek(h, 172,'bof'); nrTags = fread(h,1, 'uint32'); Tags = struct(); for i = 1:nrTags Tags(i).tag = deblank(cast(fread(h, 40, 'uint16'),'char')'); Tags(i).index = fread(h,1,'uint32'); switch Tags(i).tag case 'ExtraDataTags' Tags(i).IDStr = 'ExtraDataTags'; case 'SegmentStream' Tags(i).IDStr = 'SegmentStream'; case 'DataStream' Tags(i).IDStr = 'DataStream'; case 'InfoChangeStream' Tags(i).IDStr = 'InfoChangeStream'; case 'InfoGuids' Tags(i).IDStr = 'InfoGuids'; case '{A271CCCB-515D-4590-B6A1-DC170C8D6EE2}' Tags(i).IDStr = 'TSGUID'; case '{8A19AA48-BEA0-40D5-B89F-667FC578D635}' Tags(i).IDStr = 'DERIVATIONGUID'; case '{F824D60C-995E-4D94-9578-893C755ECB99}' Tags(i).IDStr = 'FILTERGUID'; case '{02950361-35BB-4A22-9F0B-C78AAA5DB094}' Tags(i).IDStr = 'DISPLAYGUID'; case '{8E94EF21-70F5-11D3-8F72-00105A9AFD56}' Tags(i).IDStr = 'FILEINFOGUID'; case '{E4138BC0-7733-11D3-8685-0050044DAAB1}' Tags(i).IDStr = 'SRINFOGUID'; case '{C728E565-E5A0-4419-93D2-F6CFC69F3B8F}' Tags(i).IDStr = 'EVENTTYPEINFOGUID'; case '{D01B34A0-9DBD-11D3-93D3-00500400C148}' Tags(i).IDStr = 'AUDIOINFOGUID'; case '{BF7C95EF-6C3B-4E70-9E11-779BFFF58EA7}' Tags(i).IDStr = 'CHANNELGUID'; case '{2DEB82A1-D15F-4770-A4A4-CF03815F52DE}' Tags(i).IDStr = 'INPUTGUID'; case '{5B036022-2EDC-465F-86EC-C0A4AB1A7A91}' Tags(i).IDStr = 'INPUTSETTINGSGUID'; case '{99A636F2-51F7-4B9D-9569-C7D45058431A}' Tags(i).IDStr = 'PHOTICGUID'; case '{55C5E044-5541-4594-9E35-5B3004EF7647}' Tags(i).IDStr = 'ERRORGUID'; case '{223A3CA0-B5AC-43FB-B0A8-74CF8752BDBE}' Tags(i).IDStr = 'VIDEOGUID'; case '{0623B545-38BE-4939-B9D0-55F5E241278D}' Tags(i).IDStr = 'DETECTIONPARAMSGUID'; case '{CE06297D-D9D6-4E4B-8EAC-305EA1243EAB}' Tags(i).IDStr = 'PAGEGUID'; case '{782B34E8-8E51-4BB9-9701-3227BB882A23}' Tags(i).IDStr = 'ACCINFOGUID'; case '{3A6E8546-D144-4B55-A2C7-40DF579ED11E}' Tags(i).IDStr = 'RECCTRLGUID'; case '{D046F2B0-5130-41B1-ABD7-38C12B32FAC3}' Tags(i).IDStr = 'GUID TRENDINFOGUID'; case '{CBEBA8E6-1CDA-4509-B6C2-6AC2EA7DB8F8}' Tags(i).IDStr = 'HWINFOGUID'; case '{E11C4CBA-0753-4655-A1E9-2B2309D1545B}' Tags(i).IDStr = 'VIDEOSYNCGUID'; case '{B9344241-7AC1-42B5-BE9B-B7AFA16CBFA5}' Tags(i).IDStr = 'SLEEPSCOREINFOGUID'; case '{15B41C32-0294-440E-ADFF-DD8B61C8B5AE}' Tags(i).IDStr = 'FOURIERSETTINGSGUID'; case '{024FA81F-6A83-43C8-8C82-241A5501F0A1}' Tags(i).IDStr = 'SPECTRUMGUID'; case '{8032E68A-EA3E-42E8-893E-6E93C59ED515}' Tags(i).IDStr = 'SIGNALINFOGUID'; case '{30950D98-C39C-4352-AF3E-CB17D5B93DED}' Tags(i).IDStr = 'SENSORINFOGUID'; case '{F5D39CD3-A340-4172-A1A3-78B2CDBCCB9F}' Tags(i).IDStr = 'DERIVEDSIGNALINFOGUID'; case '{969FBB89-EE8E-4501-AD40-FB5A448BC4F9}' Tags(i).IDStr = 'ARTIFACTINFOGUID'; case '{02948284-17EC-4538-A7FA-8E18BD65E167}' Tags(i).IDStr = 'STUDYINFOGUID'; case '{D0B3FD0B-49D9-4BF0-8929-296DE5A55910}' Tags(i).IDStr = 'PATIENTINFOGUID'; case '{7842FEF5-A686-459D-8196-769FC0AD99B3}' Tags(i).IDStr = 'DOCUMENTINFOGUID'; case '{BCDAEE87-2496-4DF4-B07C-8B4E31E3C495}' Tags(i).IDStr = 'USERSINFOGUID'; case '{B799F680-72A4-11D3-93D3-00500400C148}' Tags(i).IDStr = 'EVENTGUID'; case '{AF2B3281-7FCE-11D2-B2DE-00104B6FC652}' Tags(i).IDStr = 'SHORTSAMPLESGUID'; case '{89A091B3-972E-4DA2-9266-261B186302A9}' Tags(i).IDStr = 'DELAYLINESAMPLESGUID'; case '{291E2381-B3B4-44D1-BB77-8CF5C24420D7}' Tags(i).IDStr = 'GENERALSAMPLESGUID'; case '{5F11C628-FCCC-4FDD-B429-5EC94CB3AFEB}' Tags(i).IDStr = 'FILTERSAMPLESGUID'; case '{728087F8-73E1-44D1-8882-C770976478A2}' Tags(i).IDStr = 'DATEXDATAGUID'; case '{35F356D9-0F1C-4DFE-8286-D3DB3346FD75}' Tags(i).IDStr = 'TESTINFOGUID'; otherwise if isstrprop(Tags(i).tag, 'digit') Tags(i).IDStr = num2str(Tags(i).tag); else Tags(i).IDStr = 'UNKNOWN'; end end end obj.sections = Tags; %% QI index fseek(h, 172208,'bof'); obj.Qi=struct(); obj.Qi.nrEntries = fread(h,1,'uint32'); obj.Qi.misc1 = fread(h,1,'uint32'); obj.Qi.indexIdx = fread(h,1,'uint32'); obj.Qi.misc3 = fread(h,1,'uint32'); obj.Qi.LQi = fread(h,1,'uint64')'; obj.Qi.firstIdx = fread(h,nrTags,'uint64'); % Don't know what this index is for... Not required to get data and % can be huge... % fseek(h, 188664,'bof'); % Qindex = struct(); % for i = 1:obj.Qi.LQi % % Qindex(i).ftel = ftell(h); % Qindex(i).index = fread(h,2,'uint16')'; %4 % Qindex(i).misc1 = fread(h,1,'uint32'); %8 % Qindex(i).indexIdx = fread(h,1,'uint32'); %12 % Qindex(i).misc2 = fread(h,3,'uint32')'; %24 % Qindex(i).sectionIdx = fread(h,1,'uint32');%28 % Qindex(i).misc3 = fread(h,1,'uint32'); %32 % Qindex(i).offset = fread(h,1,'uint64'); % 40 % Qindex(i).blockL = fread(h,1,'uint32');%44 % Qindex(i).dataL = fread(h,1,'uint32')';%48 % end % obj.Qi.index = Qindex; %% Get Main Index: % Index consists of multiple blocks, after each block is the pointer % to the next block. Total number of entries is in obj.Qi.nrEntries Index = struct(); curIdx = 0; nextIndexPointer = indexIdx; fprintf('Parsing index '); curIdx2 = 1; while curIdx < obj.Qi.nrEntries if mod(curIdx2,20) fprintf('.'); else fprintf('\n.'); end fseek(h, nextIndexPointer, 'bof'); nrIdx = fread(h,1, 'uint64'); Index(curIdx + nrIdx).sectionIdx = 0; % Preallocate next set of indices var = fread(h,3*nrIdx, 'uint64'); for i = 1: nrIdx % Index(curIdx + i).sectionIdx = fread(h,1, 'uint64'); % Index(curIdx + i).offset = fread(h,1, 'uint64'); % Index(curIdx + i).blockL = fread(h,1, 'uint32'); % Index(curIdx + i).sectionL = fread(h,1, 'uint32'); Index(curIdx + i).sectionIdx = var(3*(i-1)+1); Index(curIdx + i).offset = var(3*(i-1)+2); Index(curIdx + i).blockL = mod(var(3*(i-1)+3),2^32); Index(curIdx + i).sectionL = round(var(3*(i-1)+3)/2^32); end nextIndexPointer = fread(h,1, 'uint64'); curIdx = curIdx + i; curIdx2=curIdx2+1; end fprintf('done\n'); obj.index = Index; obj.allIndexIDs = [obj.index.sectionIdx]; %---READ DYNAMIC PACKETS---% dynamicPackets = struct(); indexIdx = Tags(find(strcmp({Tags.IDStr},'InfoChangeStream'),1)).index; offset = Index(indexIdx).offset; nrDynamicPackets = Index(indexIdx).sectionL / 48; fseek(h, offset, 'bof'); %Read first only the dynamic packets structure without actual data for i = 1: nrDynamicPackets dynamicPackets(i).offset = offset+i*48; guidmixed = fread(h,16, 'uint8')'; guidnonmixed = [guidmixed(04), guidmixed(03), guidmixed(02), guidmixed(01), ... guidmixed(06), guidmixed(05), guidmixed(08), guidmixed(07), ... guidmixed(09), guidmixed(10), guidmixed(11), guidmixed(12), ... guidmixed(13), guidmixed(14), guidmixed(15), guidmixed(16)]; dynamicPackets(i).guid = num2str(guidnonmixed, '%02X'); dynamicPackets(i).guidAsStr = sprintf('{%02X%02X%02X%02X-%02X%02X-%02X%02X-%02X%02X-%02X%02X%02X%02X%02X%02X}', guidnonmixed); dynamicPackets(i).date = datenum(1899,12,31) + fread(h,1,'double'); dynamicPackets(i).datefrac = fread(h,1,'double'); dynamicPackets(i).internalOffsetStart = fread(h,1, 'uint64')'; dynamicPackets(i).packetSize = fread(h,1, 'uint64')'; dynamicPackets(i).data = zeros(0, 1,'uint8'); switch dynamicPackets(i).guid case 'BF7C95EF6C3B4E709E11779BFFF58EA7' dynamicPackets(i).IDStr = 'CHANNELGUID'; case '8A19AA48BEA040D5B89F667FC578D635' dynamicPackets(i).IDStr = 'DERIVATIONGUID'; case 'F824D60C995E4D949578893C755ECB99' dynamicPackets(i).IDStr = 'FILTERGUID'; case '0295036135BB4A229F0BC78AAA5DB094' dynamicPackets(i).IDStr = 'DISPLAYGUID'; case '782B34E88E514BB997013227BB882A23' dynamicPackets(i).IDStr = 'ACCINFOGUID'; case 'A271CCCB515D4590B6A1DC170C8D6EE2' dynamicPackets(i).IDStr = 'TSGUID'; case 'D01B34A09DBD11D393D300500400C148' dynamicPackets(i).IDStr = 'AUDIOINFOGUID'; otherwise dynamicPackets(i).IDStr = 'UNKNOWN'; end end %Then read the actual data from the pointers above for i = 1: nrDynamicPackets %Look up the GUID of this dynamic packet in the Tags % to find the section index infoIdx = Tags(find(strcmp({Tags.tag},dynamicPackets(i).guidAsStr),1)).index; %Matching index segments indexInstances = Index([Index.sectionIdx] == infoIdx); %Then, treat all these sections as one contiguous memory block % and grab this packet across these instances internalOffset = 0; remainingDataToRead = dynamicPackets(i).packetSize; %disp(['Target packet ' dynamicPackets(i).IDStr ' : ' num2str(dynamicPackets(i).internalOffsetStart) ' to ' num2str(dynamicPackets(i).internalOffsetStart+dynamicPackets(i).packetSize) ' target read length ' num2str(remainingDataToRead)]); currentTargetStart = dynamicPackets(i).internalOffsetStart; for j = 1: size(indexInstances,2) currentInstance = indexInstances(j); %hitInThisSegment = ''; if (internalOffset = currentTargetStart startAt = currentTargetStart; stopAt = min(startAt+remainingDataToRead, internalOffset+currentInstance.sectionL); readLength = stopAt-startAt; filePosStart = currentInstance.offset+startAt-internalOffset; fseek(h,filePosStart, 'bof'); dataPart = fread(h,readLength,'uint8=>uint8'); dynamicPackets(i).data = cat(1, dynamicPackets(i).data, dataPart); %hitInThisSegment = ['HIT at ' num2str(startAt) ' to ' num2str(stopAt)]; %if (readLength < remainingDataToRead) % hitInThisSegment = [hitInThisSegment ' (partial ' num2str(readLength) ' )']; %else % hitInThisSegment = [hitInThisSegment ' (finished - this segment contributed ' num2str(readLength) ' )']; %end %hitInThisSegment = [hitInThisSegment ' abs file pos ' num2str(filePosStart) ' - ' num2str(filePosStart+readLength)]; remainingDataToRead = remainingDataToRead-readLength; currentTargetStart = currentTargetStart + readLength; end %disp([' Index ' num2str(j) ' Offset: ' num2str(internalOffset) ' to ' num2str(internalOffset+currentInstance.sectionL) ' ' num2str(hitInThisSegment)]); internalOffset = internalOffset + currentInstance.sectionL; end end %% Get PatientGUID info = struct(); infoProps = { 'patientID', 'firstName','middleName','lastName',... 'altID','mothersMaidenName','DOB','DOD','street','sexID','phone',... 'notes','dominance','siteID','suffix','prefix','degree','apartment',... 'city','state','country','language','height','weight','race','religion',... 'maritalStatus'}; infoIdxStruct = Tags(find(strcmp({Tags.IDStr},'PATIENTINFOGUID'),1)); if ~isempty(infoIdxStruct) infoIdx = infoIdxStruct.index; indexInstance = Index(find([Index.sectionIdx]==infoIdx,1)); fseek(h, indexInstance.offset,'bof'); guid = fread(h, 16, 'uint8'); %#ok lSection = fread(h, 1, 'uint64'); %#ok % reserved = fread(h, 3, 'uint16'); %#ok nrValues = fread(h,1,'uint64'); nrBstr = fread(h,1,'uint64'); for i = 1:nrValues id = fread(h,1,'uint64'); switch id case {7,8} unix_time = (fread(h,1, 'double')*(3600*24)) - 2209161600;% 2208988800; %8 obj.segments(i).dateStr = datestr(unix_time/86400 + datenum(1970,1,1)); value = datevec( obj.segments(i).dateStr ); value = value([3 2 1]); case {23,24} value = fread(h,1,'double'); otherwise value = 0; end info.(infoProps{id}) = value; end strSetup = fread(h,nrBstr*2,'uint64'); for i=1:2:(nrBstr*2) id = strSetup(i); value = deblank(cast(fread(h, strSetup(i+1) + 1, 'uint16'),'char')'); info.(infoProps{id}) = value; end end obj.patientInfo = info; %% Get INFOGUID infoIdx = Tags(find(strcmp({Tags.IDStr},'InfoGuids'),1)).index; indexInstance = Index(find([Index.sectionIdx]==infoIdx,1)); fseek(h, indexInstance.offset,'bof'); % Ignoring, is list of GUIDS in file. %% Get SignalInfo (SIGNALINFOGUID): One per file SIG_struct = struct(); sensorIdx = Tags(find(strcmp({Tags.IDStr},'SIGNALINFOGUID'),1)).index; indexInstance = Index(find([Index.sectionIdx]==sensorIdx,1)); fseek(h, indexInstance.offset,'bof'); SIG_struct.guid = fread(h, 16, 'uint8'); SIG_struct.name = fread(h, obj.ITEMNAMESIZE, '*char'); unkown = fread(h, 152, '*char'); %#ok fseek(h, 512, 'cof'); nrIdx = fread(h,1, 'uint16'); %783 misc1 = fread(h,3, 'uint16'); %#ok obj.sigInfo = struct(); for i = 1: nrIdx obj.sigInfo(i).sensorName = deblank(cast(fread(h, obj.LABELSIZE, 'uint16'),'char')'); obj.sigInfo(i).transducer = deblank(cast(fread(h, obj.UNITSIZE, 'uint16'),'char')'); obj.sigInfo(i).guid = fread(h, 16, '*uint8'); obj.sigInfo(i).bBiPolar = logical(fread(h, 1 ,'uint32')); obj.sigInfo(i).bAC = logical(fread(h, 1 ,'uint32')); obj.sigInfo(i).bHighFilter = logical(fread(h, 1 ,'uint32')); obj.sigInfo(i).color = fread(h, 1 ,'uint32'); reserved = fread(h, 256, '*char'); %#ok end %% Get CHANNELINFO (CHANNELGUID) CH_struct = struct(); sensorIdx = Tags(find(strcmp({Tags.IDStr},'CHANNELGUID'),1)).index; indexInstance = Index(find([Index.sectionIdx]==sensorIdx,1)); fseek(h, indexInstance.offset,'bof'); CH_struct.guid = fread(h, 16, 'uint8'); CH_struct.name = fread(h, obj.ITEMNAMESIZE, '*char'); fseek(h, 152, 'cof'); CH_struct.reserved = fread(h, 16, 'uint8'); CH_struct.deviceID = fread(h, 16, 'uint8'); fseek(h, 488, 'cof'); nrIdx = fread(h,2, 'uint32'); %783 obj.chInfo = struct(); for i = 1: nrIdx(2) obj.chInfo(i).sensor = deblank(cast(fread(h, obj.LABELSIZE, 'uint16'),'char')'); obj.chInfo(i).samplingRate = fread(h,1,'double'); obj.chInfo(i).bOn = logical(fread(h, 1 ,'uint32')); obj.chInfo(i).lInputID = fread(h, 1 ,'uint32'); obj.chInfo(i).lInputSettingID = fread(h,1,'uint32'); obj.chInfo(i).reserved = fread(h,4,'char'); fseek(h, 128, 'cof'); end curIdx = 0; for i = 1: length(obj.chInfo) if obj.chInfo(i).bOn obj.chInfo(i).indexID = curIdx; curIdx = curIdx+1; else obj.chInfo(i).indexID = -1; end end %% Get TS info (TSGUID):(One per segment, last used if no new for segment) %% To simplify things, we only read the first TSINFO. tsPackets = dynamicPackets(strcmp({dynamicPackets.IDStr},'TSGUID')); if isempty(tsPackets) warning(['No TSINFO found']); else obj.tsInfos = {}; for j = 1:length(tsPackets) tsPacket = tsPackets(j); elems = typecast(tsPacket.data(753:756),'uint32'); alloc = typecast(tsPacket.data(757:760),'uint32'); offset = 761; tsInfo = struct(); for i = 1:elems internalOffset = 0; tsInfo(i).label = deblank(char(typecast(tsPacket.data(offset:(offset+obj.TSLABELSIZE-1))','uint16'))); internalOffset = internalOffset + obj.TSLABELSIZE*2; tsInfo(i).activeSensor = deblank(char(typecast(tsPacket.data(offset+internalOffset:(offset+internalOffset-1+obj.LABELSIZE))','uint16'))); internalOffset = internalOffset + obj.TSLABELSIZE; tsInfo(i).refSensor = deblank(char(typecast(tsPacket.data(offset+internalOffset:(offset+internalOffset-1+8))','uint16'))); internalOffset = internalOffset + 8; internalOffset = internalOffset + 56; tsInfo(i).dLowCut = typecast(tsPacket.data(offset+internalOffset:(offset+internalOffset-1+8))','double'); internalOffset = internalOffset + 8; tsInfo(i).dHighCut = typecast(tsPacket.data(offset+internalOffset:(offset+internalOffset-1+8))','double'); internalOffset = internalOffset + 8; tsInfo(i).dSamplingRate = typecast(tsPacket.data(offset+internalOffset:(offset+internalOffset-1+8))','double'); internalOffset = internalOffset + 8; tsInfo(i).dResolution = typecast(tsPacket.data(offset+internalOffset:(offset+internalOffset-1+8))','double'); internalOffset = internalOffset + 8; tsInfo(i).bMark = typecast(tsPacket.data(offset+internalOffset:(offset+internalOffset-1+2))','uint16'); internalOffset = internalOffset + 2; tsInfo(i).bNotch = typecast(tsPacket.data(offset+internalOffset:(offset+internalOffset-1+2))','uint16'); internalOffset = internalOffset + 2; tsInfo(i).dEegOffset = typecast(tsPacket.data(offset+internalOffset:(offset+internalOffset-1+8))','double'); offset = offset + 552; %disp([num2str(i) ' : ' TSInfo(i).label ' : ' TSInfo(i).activeSensor ' : ' TSInfo(i).refSensor ' : ' num2str(TSInfo(i).samplingRate)]); end obj.tsInfos{j} = tsInfo; end end % -- -- -- %% Get Segment Start Times segmentIdx = Tags(find(strcmp({Tags.IDStr}, 'SegmentStream'),1)).index; indexIdx = find([Index.sectionIdx] == segmentIdx, 1); segmentInstance = Index(indexIdx); nrSegments = segmentInstance.sectionL/152; fseek(h, segmentInstance.offset,'bof'); obj.segments = struct(); for i = 1: nrSegments dateOLE = fread(h,1, 'double'); obj.segments(i).dateOLE = dateOLE; unix_time = (dateOLE*(3600*24)) - 2209161600;% 2208988800; %8 obj.segments(i).dateStr = datestr(unix_time/86400 + datenum(1970,1,1)); datev = datevec( obj.segments(i).dateStr ); obj.segments(i).startDate = datev(1:3); obj.segments(i).startTime = datev(4:6); fseek(h, 8 , 'cof'); %16 obj.segments(i).duration = fread(h,1, 'double');%24 fseek(h, 128 , 'cof'); %152 end % Get nrValues per segment and channel for iSeg = 1:length(obj.segments) if iSeg > length(obj.tsInfos) obj.tsInfos{iSeg} = obj.tsInfos{iSeg-1} end % Add Channel Names to segments obj.segments(iSeg).chName = {obj.tsInfos{iSeg}.label}; obj.segments(iSeg).refName = {obj.tsInfos{iSeg}.refSensor}; obj.segments(iSeg).samplingRate = [obj.tsInfos{iSeg}.dSamplingRate]; obj.segments(iSeg).scale = [obj.tsInfos{iSeg}.dResolution]; end %% Get events - Andrei Barborica, Dec 2015 % Find sequence of events, that are stored in the section tagged 'Events' idxSection = find(strcmp('Events',{Tags.tag})); indexIdx = find([obj.index.sectionIdx] == obj.sections(idxSection).index); offset = obj.index(indexIdx).offset; ePktLen = 272; % Event packet length, see EVENTPACKET definition eMrkLen = 240; % Event marker length, see EVENTMARKER definition evtPktGUID = hex2dec({'80', 'F6', '99', 'B7', 'A4', '72', 'D3', '11', '93', 'D3', '00', '50', '04', '00', 'C1', '48'}); % GUID for event packet header HCEVENT_ANNOTATION = '{A5A95612-A7F8-11CF-831A-0800091B5BDA}'; HCEVENT_SEIZURE = '{A5A95646-A7F8-11CF-831A-0800091B5BDA}'; HCEVENT_FORMATCHANGE = '{08784382-C765-11D3-90CE-00104B6F4F70}'; HCEVENT_PHOTIC = '{6FF394DA-D1B8-46DA-B78F-866C67CF02AF}'; HCEVENT_POSTHYPERVENT = '{481DFC97-013C-4BC5-A203-871B0375A519}'; HCEVENT_REVIEWPROGRESS = '{725798BF-CD1C-4909-B793-6C7864C27AB7}'; HCEVENT_EXAMSTART = '{96315D79-5C24-4A65-B334-E31A95088D55}'; HCEVENT_HYPERVENTILATION = '{A5A95608-A7F8-11CF-831A-0800091B5BDA}'; HCEVENT_IMPEDANCE = '{A5A95617-A7F8-11CF-831A-0800091B5BDA}'; HCEVENT_AMPLIFIERDISCONNECT = '{A71A6DB5-4150-48BF-B462-1C40521EBD6F}'; HCEVENT_AMPLIFIERRECONNECT = '{6387C7C8-6F98-4886-9AF4-FA750ED300DE}'; HCEVENT_PAUSED = '{71EECE80-EBC4-41C7-BF26-E56911426FB4}'; DAYSECS = 86400.0; % From nrvdate.h fseek(h,offset,'bof'); pktGUID = fread(h,16,'uint8'); pktLen = fread(h,1,'uint64'); obj.eventMarkers = struct(); i = 0; % Event counter while (pktGUID == evtPktGUID) i = i + 1; % Please refer to EVENTMARKER structure in the Nervus file documentation fseek(h,8,'cof'); % Skip eventID, not used evtDate = fread(h,1,'double'); evtDateFraction = fread(h,1,'double'); obj.eventMarkers(i).dateOLE = evtDate; obj.eventMarkers(i).dateFraction = evtDateFraction; evtPOSIXTime = evtDate*DAYSECS + evtDateFraction - 2209161600;% 2208988800; %8 obj.eventMarkers(i).dateStr = datestr(evtPOSIXTime/DAYSECS + datenum(1970,1,1),'dd-mmmm-yyyy HH:MM:SS.FFF'); % Save fractions of seconds, as well obj.eventMarkers(i).duration = fread(h,1,'double'); fseek(h,48,'cof'); evtUser = fread(h,12,'uint16'); obj.eventMarkers(i).user = deblank(char(evtUser).'); evtTextLen = fread(h,1,'uint64'); evtGUID = fread(h,16,'uint8'); obj.eventMarkers(i).GUID = sprintf('{%.2X%.2X%.2X%.2X-%.2X%.2X-%.2X%.2X-%.2X%.2X-%.2X%.2X%.2X%.2X%.2X%.2X}',evtGUID([4 3 2 1 6 5 8 7 9:16])); fseek(h,16,'cof'); % Skip Reserved4 array evtLabel = fread(h,32,'uint16'); % LABELSIZE = 32; evtLabel = deblank(char(evtLabel).'); % Not used eventMarkers(i).label = evtLabel; %disp(sprintf('Offset: %s, TypeGUID:%s, User:%s, Label:%s',dec2hex(offset),evtGUID,evtUser,evtLabel)); % Only a subset of all event types are dealt with switch obj.eventMarkers(i).GUID case HCEVENT_SEIZURE obj.eventMarkers(i).IDStr = 'Seizure'; %disp(' Seizure event'); case HCEVENT_ANNOTATION obj.eventMarkers(i).IDStr = 'Annotation'; fseek(h,32,'cof'); % Skip Reserved5 array evtAnnotation = fread(h,evtTextLen,'uint16'); obj.eventMarkers(i).annotation = deblank(char(evtAnnotation).'); %disp(sprintf(' Annotation:%s',evtAnnotation)); case HCEVENT_FORMATCHANGE obj.eventMarkers(i).IDStr = 'Format change'; case HCEVENT_PHOTIC obj.eventMarkers(i).IDStr = 'Photic'; case HCEVENT_POSTHYPERVENT obj.eventMarkers(i).IDStr = 'Posthyperventilation'; case HCEVENT_REVIEWPROGRESS obj.eventMarkers(i).IDStr = 'Review progress'; case HCEVENT_EXAMSTART obj.eventMarkers(i).IDStr = 'Exam start'; case HCEVENT_HYPERVENTILATION obj.eventMarkers(i).IDStr = 'Hyperventilation'; case HCEVENT_IMPEDANCE obj.eventMarkers(i).IDStr = 'Impedance'; case HCEVENT_AMPLIFIERDISCONNECT obj.eventMarkers(i).IDStr = 'Amplifier Disconnect'; case HCEVENT_AMPLIFIERRECONNECT obj.eventMarkers(i).IDStr = 'Amplifier Reconnect'; case HCEVENT_PAUSED obj.eventMarkers(i).IDStr = 'Recording Paused'; otherwise obj.eventMarkers(i).IDStr = 'UNKNOWN'; end % Next packet offset = offset + pktLen; fseek(h,offset,'bof'); pktGUID = fread(h,16,'uint8'); pktLen = fread(h,1,'uint64'); end %% Get montage - Andrei Barborica, Dec 2015 % Derivation (montage) mtgIdx = Tags(find(strcmp({Tags.IDStr},'DERIVATIONGUID'),1)).index; indexIdx = find([obj.index.sectionIdx]==mtgIdx,1); fseek(h,obj.index(indexIdx(1)).offset + 40,'bof'); % Beginning of current montage name mtgName = deblank(char(fread(h,32,'uint16')).'); fseek(h,640,'cof'); % Number of traces in the montage numDerivations = fread(h,1,'uint32'); numDerivations2 = fread(h,1,'uint32'); obj.montage = struct(); for i = 1:numDerivations obj.montage(i).derivationName = deblank(char(fread(h,64,'uint16')).'); obj.montage(i).signalName1 = deblank(char(fread(h,32,'uint16')).'); obj.montage(i).signalName2 = deblank(char(fread(h,32,'uint16')).'); fseek(h,264,'cof'); % Skip additional info end % Display properties dispIdx = Tags(find(strcmp({Tags.IDStr},'DISPLAYGUID'),1)).index; indexIdx = find([obj.index.sectionIdx]==dispIdx,1); fseek(h,obj.index(indexIdx(1)).offset + 40,'bof'); % Beginning of current montage name displayName = deblank(char(fread(h,32,'uint16')).'); fseek(h,640,'cof'); % Number of traces in the montage numTraces = fread(h,1,'uint32'); numTraces2 = fread(h,1,'uint32'); if (numTraces == numDerivations) for i = 1:numTraces fseek(h,32,'cof'); obj.montage(i).color = fread(h,1,'uint32'); % Use typecast(uint32(montage(i).color),'uint8') to convert to RGB array fseek(h,136-4,'cof'); end else disp('Could not match montage derivations with display color table'); end % Close File fclose(h); end function out = getNrSamples(obj, segment) % GETNRSAMPLES Returns the number of samples per channel in segment. % % OUT = GETNRSAMPLES(OBJ, SEGMENT) returns a 1xn array of values % indicating the number of samples for each of the channels in the % associated SEGMENT, where SEGMENT is the index of the % OBJ.segments array. assert(length(obj.segments)>= segment, ... 'Incorrect SEGMENT argument; must be integer representing segment index.'); out = obj.segments(segment).samplingRate .* obj.segments(segment).duration; end function cSumSegs = getCSumSegs(obj, chI) % GETCSUMSEGS Returns the cumulative sum of a channels segments % % CSUMSEGS = GETCSUMSEGS(OBJ, CHI) returns a 1xn array of the % cumulative sum of the given channel, where chI is the channel % index. Different TsInfos can have different numbers of channels % and different sampling rates. samplingRates = zeros(1, length(obj.tsInfos)); for ts_i = 1:length(obj.tsInfos) if length(obj.tsInfos{ts_i}) < chI continue else samplingRates(ts_i) = obj.tsInfos{ts_i}(chI).dSamplingRate; end end cSumSegs = [0 cumsum(samplingRates.*[obj.segments.duration])]; end function out = getdata(obj, segment, range, chIdx) % GETDATA Returns data from Nicolet file. % % OUT = GETDATA(OBJ, SEGMENT, RANGE, CHIDX) returns data in an nxm array of % doubles where n is the number of datapoints and m is the number % of channels. RANGE is a 1x2 array with the [StartIndex EndIndex] % and CHIDX is a vector of channel indeces. % Assert range is 1x2 vector if isempty([obj.tsInfos{segment}.label]) warning('Segment %d has an empty tsInfo. Skipping', segment) out = 0; return end assert(length(range) == 2, 'Range is [firstIndex lastIndex]'); assert(length(segment) == 1, 'Segment must be single value.'); % Reopen .e file. h = fopen(obj.fileName,'r','ieee-le'); % Find sectionID for channels lChIdx = length(chIdx); sectionIdx = zeros(lChIdx,1); for i = 1:lChIdx tmp = find(strcmp(num2str(chIdx(i)-1),{obj.sections.tag}),1); sectionIdx(i) = obj.sections(tmp).index; end % Iterate over all requested channels and populate array. out = zeros(range(2) - range(1) + 1, lChIdx); for i = 1 : lChIdx % Get cumulative sum segments. cSumSegments = obj.getCSumSegs(chIdx(i)); % Get sampling rate for current channel curSF = obj.segments(segment).samplingRate(chIdx(i)); mult = obj.segments(segment).scale(chIdx(i)); % Find all sections allSectionIdx = obj.allIndexIDs == sectionIdx(i); allSections = find(allSectionIdx); % Find relevant sections sectionLengths = [obj.index(allSections).sectionL]./2; cSectionLengths = [0 cumsum(sectionLengths)]; skipValues = cSumSegments(segment); firstSectionForSegment = find(cSectionLengths > skipValues, 1) - 1 ; lastSectionForSegment = firstSectionForSegment + ... find(cSectionLengths > curSF*obj.segments(segment).duration,1) - 2 ; if isempty(lastSectionForSegment) lastSectionForSegment = length(cSectionLengths); end offsetSectionLengths = cSectionLengths - cSectionLengths(firstSectionForSegment); firstSection = find(offsetSectionLengths < range(1) ,1,'last'); lastSection = find(offsetSectionLengths >= range(2),1)-1; if isempty(lastSection) lastSection = length(offsetSectionLengths); end if lastSection > lastSectionForSegment error('Index out of range for current section: %i > %i, on channel: %i', ... range(2), cSectionLengths(lastSectionForSegment+1), chIdx(i)); end useSections = allSections(firstSection: lastSection) ; useSectionL = sectionLengths(firstSection: lastSection) ; % First Partial Segment curIdx = 1; curSec = obj.index(useSections(1)); fseek(h, curSec.offset,'bof'); firstOffset = range(1) - offsetSectionLengths(firstSection); lastOffset = min([range(2) useSectionL(1)]); lsec = lastOffset-firstOffset + 1; fseek(h, (firstOffset-1) * 2,'cof'); out(1 : lsec,i) = fread(h, lsec, 'int16') * mult; curIdx = curIdx + lsec; if length(useSections) > 1 % Full Segments for j = 2: (length(useSections)-1) curSec = obj.index(useSections(j)); fseek(h, curSec.offset,'bof'); out(curIdx : (curIdx + useSectionL(j) - 1),i) = ... fread(h, useSectionL(j), 'int16') * mult; curIdx = curIdx + useSectionL(j); end % Final Partial Segment curSec = obj.index(useSections(end)); fseek(h, curSec.offset,'bof'); out(curIdx : end,i) = fread(h, length(out)-curIdx + 1, 'int16') * mult; end end % Close the .e file. fclose(h); end function out = getdataQ(obj, segment, range, chIdx) % GETDATAQ Returns data from Nicolet file. This is a "QUICK" version of getdata, % that uses more memory but operates faster on large datasets by reading % a single block of data from disk that contains all data of interest. % % OUT = GETDATAQ(OBJ, SEGMENT, RANGE, CHIDX) returns data in an nxm array of % doubles where n is the number of datapoints and m is the number % of channels. RANGE is a 1x2 array with the [StartIndex EndIndex] % and CHIDX is a vector of channel indeces. % % Andrei Barborica, Dec 2015 % if isempty([obj.tsInfos{segment}.label]) warning('Segment %d has an empty tsInfo. Skipping', segment) out = 0; return end % Assert range is 1x2 vector assert(length(range) == 2, 'Range is [firstIndex lastIndex]'); assert(length(segment) == 1, 'Segment must be single value.'); % Reopen .e file. h = fopen(obj.fileName,'r','ieee-le'); % Find sectionID for channels lChIdx = length(chIdx); sectionIdx = zeros(lChIdx,1); for i = 1:lChIdx tmp = find(strcmp(num2str(chIdx(i)-1),{obj.sections.tag}),1); sectionIdx(i) = obj.sections(tmp).index; end usedIndexEntries = zeros(size([obj.index.offset])); % Iterate over all requested channels and populate array. out = zeros(range(2) - range(1) + 1, lChIdx); for i = 1 : lChIdx % Get cumulative sum segments. cSumSegments = obj.getCSumSegs(chIdx(i)); % Get sampling rate for current channel curSF = obj.segments(segment).samplingRate(chIdx(i)); mult = obj.segments(segment).scale(chIdx(i)); % Find all sections allSectionIdx = obj.allIndexIDs == sectionIdx(i); allSections = find(allSectionIdx); % Find relevant sections sectionLengths = [obj.index(allSections).sectionL]./2; cSectionLengths = [0 cumsum(sectionLengths)]; skipValues = cSumSegments(segment); firstSectionForSegment = find(cSectionLengths > skipValues, 1) - 1 ; lastSectionForSegment = firstSectionForSegment + ... find(cSectionLengths > curSF*obj.segments(segment).duration,1) - 2 ; if isempty(lastSectionForSegment) lastSectionForSegment = length(cSectionLengths); end offsetSectionLengths = cSectionLengths - cSectionLengths(firstSectionForSegment); firstSection = find(offsetSectionLengths < range(1) ,1,'last'); lastSection = find(offsetSectionLengths >= range(2),1)-1; if isempty(lastSection) lastSection = length(offsetSectionLengths); end if lastSection > lastSectionForSegment error('Index out of range for current section: %i > %i, on channel: %i', ... range(2), cSectionLengths(lastSectionForSegment+1), chIdx(i)); end useSections = allSections(firstSection: lastSection) ; useSectionL = sectionLengths(firstSection: lastSection) ; % First Partial Segment usedIndexEntries(useSections(1)) = 1; if length(useSections) > 1 % Full Segments for j = 2: (length(useSections)-1) usedIndexEntries(useSections(j)) = 1; end % Final Partial Segment usedIndexEntries(useSections(end)) = 1; end end % Read a big chunk of the file, containing data of interest. ix = find(usedIndexEntries); fseek(h, obj.index(ix(1)).offset,'bof'); dsize = obj.index(ix(end)).offset - obj.index(ix(1)).offset + obj.index(ix(end)).sectionL; tmp = fread(h,dsize/2,'int16').'; % Close the .e file. fclose(h); baseOffset = obj.index(ix(1)).offset; % Extract specified channels for i = 1 : lChIdx % Get cumulative sum segments. cSumSegments = obj.getCSumSegs(chIdx(i)); % Get sampling rate for current channel curSF = obj.segments(segment).samplingRate(chIdx(i)); mult = obj.segments(segment).scale(chIdx(i)); % Find all sections allSectionIdx = obj.allIndexIDs == sectionIdx(i); allSections = find(allSectionIdx); % Find relevant sections sectionLengths = [obj.index(allSections).sectionL]./2; cSectionLengths = [0 cumsum(sectionLengths)]; skipValues = cSumSegments(segment); firstSectionForSegment = find(cSectionLengths > skipValues, 1) - 1 ; lastSectionForSegment = firstSectionForSegment + ... find(cSectionLengths > curSF*obj.segments(segment).duration,1) - 2 ; if isempty(lastSectionForSegment) lastSectionForSegment = length(cSectionLengths); end offsetSectionLengths = cSectionLengths - cSectionLengths(firstSectionForSegment); firstSection = find(offsetSectionLengths < range(1) ,1,'last'); lastSection = find(offsetSectionLengths >= range(2),1)-1; if isempty(lastSection) lastSection = length(offsetSectionLengths); end if lastSection > lastSectionForSegment error('Index out of range for current section: %i > %i, on channel: %i', ... range(2), cSectionLengths(lastSectionForSegment+1), chIdx(i)); end useSections = allSections(firstSection: lastSection) ; useSectionL = sectionLengths(firstSection: lastSection) ; % First Partial Segment curIdx = 1; curSec = obj.index(useSections(1)); %fseek(h, curSec.offset,'bof'); firstOffset = range(1) - offsetSectionLengths(firstSection); lastOffset = min([range(2) useSectionL(1)]); lsec = lastOffset-firstOffset + 1; out(1 : lsec,i) = tmp( (curSec.offset - baseOffset)/2 + (firstOffset-1) + (1:lsec) ) * mult; curIdx = curIdx + lsec; if length(useSections) > 1 % Full Segments for j = 2: (length(useSections)-1) curSec = obj.index(useSections(j)); out(curIdx : (curIdx + useSectionL(j) - 1),i) = ... tmp( (curSec.offset - baseOffset)/2 + (1:useSectionL(j)) ) * mult; curIdx = curIdx + useSectionL(j); end % Final Partial Segment curSec = obj.index(useSections(end)); out(curIdx : end,i) = tmp( (curSec.offset - baseOffset)/2 + (1:(length(out)-curIdx + 1)) ) * mult; % length(out) ?????? end end end function labels = getlabels(obj,str) % Returns annotations containing specified string % % Cristian Donos, Dec 2015 % labels=[]; counter = 1; for i = 1:length(obj.eventMarkers) if strfind(lower(obj.eventMarkers(i).annotation),lower(str)) labels{counter,1} = obj.eventMarkers(i).annotation; % annotation string labels{counter,2} = i; % annotation index in obj.eventMarkers % identify segment time_vector = []; for j = 1:length(obj.segments) time_vector = [time_vector etime(datevec(obj.eventMarkers(i).dateStr),datevec(obj.segments(j).dateStr))]; end labels{counter,3}= find(time_vector==min(time_vector(time_vector>0))); % annotation part of this segment labels{counter,4}= min(time_vector(time_vector>0)); % annotation offset in seconds, relative to its segment start counter = counter+1; end end end end end obj = NicoletFile('janbrogger.e'); % save data out = getdata(obj, 1, [1 20000], 1:16); plot(out) |
【本文地址】
今日新闻 |
推荐新闻 |