MATLAB Answers

Non-looped fread of binary file with multiple value types

3 views (last 30 days)
Josh Parks
Josh Parks on 17 Sep 2012
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

Answers (4)

Josh Parks
Josh Parks on 27 Sep 2012
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 on 18 Sep 2012
  • 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 on 20 Sep 2012
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 Comments
Jan
Jan on 20 Sep 2012
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.

Sign in to comment.


Jan
Jan on 20 Sep 2012
Edited: Jan on 20 Sep 2012
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

Community Treasure Hunt

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

Start Hunting!