Non-looped fread of binary file with multiple value types

2012 9 17
4 个回答
6 次查看（30 天）

0 个投票

Hi all,
This is my first post here and I am having quite the trouble with this little program. I am trying to read data out of a proprietary binary generated by a photon counting card. I've got it reading the data correctly but it is super slow. I'm guessing this is because of my implementation of fread, however I can't figure out how to do this without a loop. Any suggestions from you matlab gurus out there? Currently it takes about nine seconds on my computer to read in a 452kB file. I've attached my program with a set of example data. For those of you who don't care to download anything, here is the code excerpt in suspect:
% This reads the TTTR event records
Ofltime = 0;
Counts = 0;
outfile = ['t3r.out'];
fpout = fopen(outfile,'W');
fprintf(1,'Writing data to %s\n', outfile);
fprintf(1,'This may take a while...');
% fprintf(fpout,' # timetag time/s channel route\n\n');
for i=1:NumberOfRecords
TTTRRecord = fread(fid, 1, 'uint32');
TimeTag = bitand(TTTRRecord,65535); %the lowest 16 bits
Channel = bitand(bitshift(TTTRRecord,-16),4095); %the next 12 bits
Route = bitand(bitshift(TTTRRecord,-28),3); %the next 2 bits
Valid = bitand(bitshift(TTTRRecord,-30),1); %the next bit
Truetime = (Ofltime + TimeTag) * TTTRGlobclock * 1e-9;
if Valid
fprintf(fpout,'%7u %7u %11.7f %5u %2u', Counts, TimeTag, Truetime, Channel, Route);
Counts = Counts+1;
else %this means we have a special record
if bitand(Channel,2048)
Ofltime = Ofltime+65536;
end;
if bitand (Channel,7)
fprintf(fpout,'Marker=%u\n',bitand(Channel,7));
end;
end;
if Valid
fprintf(fpout,'\n');
end;
end;
Thanks,
Josh
Program zip
https://docs.google.com/open?id=0B84phdetR0DkYktXT2pxWGt2YUE

请先登录，再回答此问题。

回答（4 个）

Josh Parks
Josh Parks 2012-9-27

1 个投票

Hi all,
in case anyone was wondering, this is what I came up with:
%% Data Processing
% This reads the TTTR event records
outfile = [filename(1:length(filename)-4),'.txt'];
fprintf(1,'Analyzing data from %s\n', outfile);
fprintf(1,'This may take a few seconds...');
Ofltime = 0;
Counts = 0;
TTTRData = fread(fid, [1 NumberOfRecords], 'uint32');
TTTRProcessed = zeros(NumberOfRecords,5);
%read values out of bitrecord for each photon and assimilate
%accordinlgy
position = ftell(fid)
ValidA = bitand(bitshift(TTTRData(:),-30),ones(length(TTTRData),1)) == 1;
CountsA = cumsum(ValidA);
ChannelA = bitand(bitshift(TTTRData(:),-16),ones(length(TTTRData),1).*4095);
%special record
OfltimeA = cumsum(bitand(ChannelA(:),ones(length(TTTRData),1).*2048)./2048.*65536.*(~ValidA));
TimeTagA = bitand(TTTRData(:),ones(length(TTTRData),1).*65535);
TrueTimeA = (OfltimeA + TimeTagA).*TTTRGlobclock.*1e-9;
RouteA = bitand(bitshift(TTTRData(:),-28),ones(length(TTTRData),1).*3);
%zero unimportant data tags reduction
ChannelA = ChannelA.*ValidA;
TimeTagA = TimeTagA.*ValidA;
RouteA = RouteA.*ValidA;
TrueTimeA = TrueTimeA.*ValidA;
CountsA = CountsA.*ValidA;
%concatenate values for single matrix printing/saving
TTTRProcessed = horzcat(CountsA,TimeTagA,TrueTimeA,ChannelA,RouteA);
%removes non-data rows
TTTRProcessed = TTTRProcessed(any(TTTRProcessed,2),:);
if printTXT == 1
fprintf(1,'Writing data to %s\n', outfile);
fpout = fopen(outfile,'W');
fprintf(fpout,'%7u\t%7u\t%11.7f\t%5u\t%2u\n',TTTRProcessed.');
fclose(fid);
fclose(fpout);
fprintf(1,'\n%u photon records written to %s\n\n',max(CountsA),outfile);
end
Walter Roberson
Walter Roberson 2012-9-18

0 个投票

  • you can read the entire file into memory and process it once it is in memory
  • if Valid is not true, you can skip calculation of Route and so on.
  • You can extract Valid with a single bitand() without doing a shift, as you only care if it is zero or non-zero
  • assign the result of bitand(Channel,7) to a variable so you do not need to recalculate it for printing purposes
Josh Parks
Josh Parks 2012-9-20

0 个投票

Thanks walter, but is there a way to completely vectorize the data read in? This is what I've come up with so far, but it still takes 2 minutes to run on large files because of the loop
% This reads the TTTR event records
Ofltime = 0;
Counts = 0;
outfile = ['t3r.out'];
fpout = fopen(outfile,'W');
fprintf(1,'Writing data to %s\n', outfile);
fprintf(1,'This may take a while...');
TTTRData = fread(fid, [1 NumberOfRecords], 'uint32');
TTTRProcessed = zeros(NumberOfRecords,5);
% fprintf(fpout,' # timetag time/s channel route\n\n');
for i=1:NumberOfRecords
TTTRRecord = TTTRData(i);
Valid = bitand(bitshift(TTTRRecord,-30),1); %the next bit
Channel = bitand(bitshift(TTTRRecord,-16),4095); %the next 12 bits
if Valid
TimeTag = bitand(TTTRRecord,65535); %the lowest 16 bits
Route = bitand(bitshift(TTTRRecord,-28),3); %the next 2 bits
Truetime = (Ofltime + TimeTag) * TTTRGlobclock * 1e-9;
TTTRProcessed(i,1) = Counts;
TTTRProcessed(i,2) = TimeTag;
TTTRProcessed(i,3) = Truetime;
TTTRProcessed(i,4) = Channel;
TTTRProcessed(i,5) = Route;
% fprintf(fpout,'%7u %7u %11.7f %5u %2u', Counts, TimeTag, Truetime, Channel, Route);
Counts = Counts+1;
else %this means we have a special record
if bitand(Channel,2048)
Ofltime = Ofltime+65536;
end;
if bitand(Channel,7)
TTTRProcessed(i,1) = -1;
TTTRProcessed(i,2) = -1;
TTTRProcessed(i,3) = -1;
TTTRProcessed(i,4) = -1;
TTTRProcessed(i,5) = -1;
end;
end;
end;
%removes non-data rows
TTTRProcessed = TTTRProcessed(any(TTTRProcessed,2),:);
fprintf(fpout,'%7u\t%7u\t%11.7f\t%5u\t%2u\n',TTTRProcessed.');
fclose(fid);
fclose(fpout);
fprintf(1,'\n%u photon records written to %s\n',Counts,outfile);

2 个评论
显示 无 隐藏 无

Walter Roberson
Walter Roberson 2012-9-20
bitand() and bitshift() both work on arrays, so you can do all those extractions at once (possibly including TimeTag and Route). After that you could use logical indexing to select the Valid entries for printing.
Jan
Jan 2012-9-20
Try if a multiplication is much faster than bitshift(). This has been the case at least in older Matlab versions, as far as I remember in R2009a.

请先登录，再进行评论。

Jan
Jan 2012-9-20
编辑：Jan 2012-9-20

0 个投票

Small but easy improvement:
TTTRRecord = TTTRData(i);
Valid = bitand(bitshift(TTTRRecord,-30),1); %the next bit
Channel = bitand(bitshift(TTTRRecord,-16),4095); %the next 12 bits
if Valid
...
To:
TTTRRecord = TTTRData(i);
Valid = bitand(TTTRRecord), 1073741824); % 2^30
if Valid
Channel = bitand(bitshift(TTTRRecord,-16),4095); %the next 12 bits
...
Perhaps:
Channel = bitand(TTTRRecord/65536, 4095); %the next 12 bits
Another standard procedure ist to move all repeated calculations out of the loop, e.g.:
c1 = TTTRGlobclock * 1e-9;
...
Truetime = (Ofltime + TimeTag) * c1;
Another improvement: Instead of 5 single calls use:
TTTRProcessed(i,1:5) = -1;
c1 = TTTRGlobclock * 1e-9;
TTTRProcessed = zeros(NumberOfRecords,5);
for i=1:NumberOfRecords
TTTRRecord = TTTRData(i);
Valid = bitand(TTTRRecord), 1073741824); % 2^30
if Valid
Channel = bitand(TTTRRecord / 65536, 4095); % next 12 bits
TimeTag = bitand(TTTRRecord, 65535); % lowest 16 bits
Route = bitand(TTTRRecord / 268435456), 3); % next 2 bits
Truetime = (Ofltime + TimeTag) * c1;
TTTRProcessed(i,1) = Counts;
TTTRProcessed(i,2) = TimeTag;
TTTRProcessed(i,3) = Truetime;
TTTRProcessed(i,4) = Channel;
TTTRProcessed(i,5) = Route;
Counts = Counts + 1;
else %this means we have a special record
if bitand(Channel,2048)
Ofltime = Ofltime + 65536;
end
if bitand(Channel, 7)
TTTRProcessed(i, 1:5) = -1;
end
end
end

类别

在 帮助中心 和 File Exchange 中查找有关 HDF5 Files 的更多信息

产品

标签

提问:

Josh Parks
Josh Parks
2012-9-17

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by