Speeding up fread for random position in big file

dear community,

at the moment i need to read one or multiple image(s) from a big data file. since I build a live data viewer and the file is too big (20-80 GB) too read all at once i read it one by one. at the moment I use the following command

% this only once
fid=fopen('myfile.mraw','r');

% and this each time I need a specific or multiple image(s)
first_frame=251674 % exemplary frame number
numOfFrames=1 % can be 1 or 3, user dependent (3 only if i apply a temporal median filter of length 3)
bitOrder = 'b';
color = 1; % just gray in my case and not color
N=[320 384];
pixels=320*384;
ColorBit=12;

I=zeros(Pixels*numOfFrames,1,'uint16');

start = (first_frame-1)*pixels*ColorBit/8;
fseek(fid,start,'bof');
I = fread(fid,Pixels*numOfFrames,'ubit12=>uint16',bitOrder);

% convert to image matrix
A(:,:,:)=permute(reshape(I,[N numOfFrames]),[2 1 3 4]); 

is there an easy way to speed up the fread call? compared to the the further displaying etc. this line takes about 97% of the time.

when reading one frame per call I dont get more than 18--20 frames out per second. at least 50 fps would be nice

unfortunately i cannot supply the mraw file due to its size.

reading more frames at once is not possible here, since first_frame can change e.g. by 1000 between calls.

The file contains high speed recordings at 20k fps, together with other data.

best regards

Jonas

[dpb]

I = fread(fid,Pixels+numOfFrames,'ubit12=>uint16',bitOrder);

I'd guess the 'ubit12=>uint16' translation is the bottleneck here, not the reading itself.

The '+' sign there looks suspicious to me, should that not be '*' instead?

While by your statement it isn't much of a slowdown, the initializtion of I to the zeros array is superluous here since are doing a complete array assignment; I'd retime with it commented out just to see...

[Jonas]

the plus was just a error from me here, in the actual code i wrote * ;-) -> corrected!

i doubt that commenting out the will solve that, the profiler did not show much time consumption there, but I will try

edit: commenting out did not have any effect.

[Bjorn Gustavsson]

An uggly workaround might be to make a 16-bit version of your original data-file (yeah, it will be a bit bigger). If the conversion-on-read is the bottle-neck during multiple viewings/processings it might be worthwhile to do the conversion once, and pay the conversion-time once, and have a smoother experience ever-after.

[Jonas]

@Bjorn Gustavsson the problem is, that there are not only images in this file, also up to 10 other channels with non image data like sound or vibration. Beside that, I have many files of that kind, together more multiple terabyte actually. And editing and replaceing the files may cause errors to the native software which uses this file.

can one of the following be adapted to my case? i was not able to do it

https://www.mathworks.com/matlabcentral/answers/28178-how-to-read-12bit-raw-data-faster

https://www.mathworks.com/matlabcentral/answers/276069-how-do-i-read-in-a-binary-file-that-has-a-very-unique-data-structure

[dpb]

We knew it wouldn't make any big difference in overall time to not execute the initialization -- but there's no point in doing something that's just going to be overwritten later, anyways.

I would guess @Bjorn Gustavsson's suggestion is the only way you'll gain any throughput in the end application; the 12-->16bit conversion is going to be what it is; you don't have any control over the internals other than try to roll your own in mex.

Native machine-width storage might be somewhat faster than the 16-bit type that is not word-aligned in storage, but with large images you may not be able to afford the extra storage/memory required -- Catch-22. Efficiency requires larger memory but memory is still limited despite modern GB RAM.

[Bjorn Gustavsson]

Maybe you can read the data with ubit12 and then convert the frame to uint16. That might be faster. I certainly did not mean that you would replace the original files with the preprocessed ones, just that it might be a better user-interface to first separately process/conver the file-of-interest into 16-bit frames. How to handle the other data in the original file I don't know, and that might depend on what your ui-viewing-function is expected to handle. If reading that format is swifter it might still be worth weasling around the 12-bit-reading-and-converting during the viewing-phase.

[dpb]

"can one of the following be adapted to my case?"

That depends on how the data are actually packed in the file -- in the example @Jan gave of maksing off the upper nyble, the data were actually written as 16-bit integers, just other stuff in the upper 4 bits. Sounds like maybe they actually did pack your data, but we don't know; you've not shared that crucial piece of information.

..." editing and replaceing the files may cause errors to the native software..."

If there are applications available to read these files, it might be better to wrap those in mex functions to use them -- or, if the vendor supplies the C or other source code, turn it directly into mex.

Again, you've left us in the dark without all the necessary information to have any real insight other than direct 12-->16 bit translation takes time (which is obvious it must).

[dpb]

"I certainly did not mean that you would replace the original files with the preprocessed ones,..."

One wouldn't replace the original files, of course not. It might, however, well be worth preprocessing them and build sets of files that contain the data in more i/o efficient rather than storage-efficient format. That could be done offline/overnight in batch processes so the end user never sees it. Just keep a database of the auxiliary files associated with the vendor file.

Can the vendor-supplied codes process images much more rapidly? If so, then they're playing some sort of tricks; at the root level MATLAB fread devolves to the system library i/o routines and they're all the same or essentially the same for everybody's compiler; the MATLAB version is somewhat different in that TMW has vectorized it for MATLAB, but under the hood the disk transfer is still what it is; and almost all applications are built with the same vendor toolsets which use the same standard libraries.

[dpb]

I = fread(fid,Pixels*numOfFrames,'ubit12=>uint16',bitOrder);

has a syntax error, too. The bitOrder parameter is in the position of the skip parameter; this probably is another typo as it likely wouldn't run at all as is, but if fread is actually interpreting the 'b' as a number, that's a problem.

I = fread(fid,Pixels*numOfFrames,'ubit12=>uint16',0,bitOrder);

would appear to be what is wanted/needed.

[Jonas]

the syntax is alright, this way it also used in the doc of fread. the resulting image(s) are also right, the command works as intended and the images are correct

Answer 1

[Walter Roberson]

In my experience, you can be more efficient on conversion of 12 bit data than what is done by fread()

If memory serves

raw = uint16(fread(fid, [3 size], '*uint8'));
firsts = First12(raw(1,:),raw(2,:));
seconds = Second12(raw(2,:),raw(3,:));

where First12 and Second12 are lookup tables such that First12(A,B) = A*256 + bitand(B, 240)/16) and Second12(A,B) = bitand(A,15)*256 + B

One might think that it is more efficient to just do those numeric computations for the whole array of raw data, without using lookup tables -- and if you were to write a bit of C code to do the conversion, it might well be more efficient to do that. But my memory is that in MATLAB, it turns out to be more efficient to pre-compute the values and use array lookup -- since the pre-computation is over a relatively small array compared to the input values, and table lookup is comparatively fast since it only involves address arithmetic and pulling values out of memory. Applying those bitand() and multiplications over the whole large array is, if my memory is correct, slower at the MATLAB level (but again, some C code might well rip through it.)

[Jonas]

so if i understand correctly, the code reads 3 bytes as uint8, which contain two 12 bit values. after that byte 1 and 2 are used to generate the first 12 bit value and then byte 2 and 3 are used for the second 12 bit value.

following that, 'size' in your code is ceil(numOfPixels*numOfFrames/2)

which output class will then be 'firsts' and 'seconds', with the given formula? uint16?

[Walter Roberson]

Yes, they would be uint16. But if you needed them as double that would be an easy change.

[Jonas]

first, the brace in the first formular is missing, it is (A*256 + bitand(B, 240))/16 if I am not mistaken

alright, update on this: actually this IS much faster than the previous syntax. I ran into a problem because

1) in the presented form I need a 2d lookup table

2) the input to the lookup table is a vector

applying arrayfun to solve 2) leads to much poorer performance and 2) interfers with 1)

i solved this by decoupling the lookup tables:

LUT1Part1=zeros(2^16,1,'uint16');
    for A=0:2^8-1
        LUT1Part1(A+1)= A*256/16;
    end

LUT1Part2=zeros(2^16,1,'uint16');
    for B=0:2^8-1
        LUT1Part2(B+1)= bitand(B, 240)/16;
    end
    
LUT2=zeros(2^16,1,'uint16');
    for A=0:2^8-1
        LUT2(A+1)= bitand(A,15)*256;
    end
    

later, i read and convert the tables as follows:

I use intlut from Matlab to lookup the values, this is also the reason why the lookup tables need 2^16 entries although 2^8 would be enough

    raw = uint16(fread(fid, [3 ceil(obj.Pixels/2)], 'uint8',bitOrder));
            
    firsts = intlut(raw(1,:),LUT1Part1) + intlut(raw(2,:),LUT1Part2);
    seconds = intlut(raw(2,:),LUT2) + raw(3,:); % here we can find the second part of Walter's conversion formular for second parts

    I=[firsts;seconds];
    
    % reshape also takes the elements columns wise, so values of firsts
    % and seconds do alternate
    A(:,:)=permute(reshape(I,N),[2 1 3]);  

the elapsed time to load single 1000 random images is about

29.660486 seconds for the old read and processing method

to 7.204475 seconds for this method

[Bjorn Gustavsson]

You shouldn't need to have 2-D lookup tables. Matlab allows you to do thing like this:

first = zeros(3,4);
second = zeros(3,4);
first(:) = randn(12,1);

So if you know the image dimension you should be able to utilize this pattern.

[Walter Roberson]

When I was designing the code with 2D lookup tables, I could have used linear indexing by multiplying one input by a constant and adding the other -- but that would have required an explicit multiply and add. I figured that it would probably be faster to use 2D indices, with the implicit lookup, as the Execution Engine would be able to convert that into low-level multiply-and-add instead of explicit. It might even possibly be able to convert to base register + offset machine code.

Answer 2

[Steven Lord]

Since your data is Big Data (too large to fit in memory all at once) you may want to investigate creating a datastore for your files and using that datastore to create a tall array on which to operate. See this section of the documentation for more information about the tools and techniques you can use for working with tall arrays and datastore objects.

[dpb]

My playing with the datastore and tall arrays is while they are convenient from coding standpoint; that convenience comes at a pretty noticeable performance hit. Not sure it'll help here where throughput is the key....but I'll admit I've not tried again very recently so maybe performance had gotten better.

[Jonas]

is this computationally competitive compared to fread which also does not load the data into memory?

[dpb]

"...compared to fread which also does not load the data into memory?"

Say What!!!??? Loading data into memory is precisely what fread does (with translation in your case as well, of course, which is probably a major fraction of the actual time being spent -- in the one case previous, it was almost a 25X speedup to read 16-bits and mask off upper four as compared to the implicit conversion of 12 bits to 16 inside fread.

Again, we need to see what the actual format of the data file is -- is it really packed in 12-bit chunks, etc. ...

The (quite minimal, only playing with it in conjunction with other Q? here in forum) experience I have with the datastore and tall arrays was that they did have a fair amount of performance overhead when it came to actually addressing the portions not presently in memory -- a straight fread of sections if could process in independent chunks not requiring all data to be processed to get result was noticeably faster albeit more coding.

[Bjorn Gustavsson]

Here's 100 SEK saying that the suggestion I've made above is the one you will have to resort to eventually. It is a wellknown programming pattern of "divide and conquer" (not in the q-sort way but the label still fits.)

[Jonas]

@dpb of course the specified portion is loaded into the memory, but not he whole file. that is what i wanted to say

how could i show you how the actual data looks like whithout the need to upload the whole file?

[dpb]

"how could i show you how the actual data looks like without the need to upload the whole file?"

You stated above this is from a vendor product and there are other coding tools for the files -- is there not documentation that gives the actual file structure detail? --would be the most convenient way. Tell us what it is you're working with; give us links to documentation -- or post the pertinent details from the doc here.

[Jonas]

i am working with high speed data which is managed by the Photron Fastcam viewer. Version is 3 if I am not mistaken. the best i could find about file format is here

for reading outside of this viewer i found two functions on fileexchange and a python bariant on github

mraw description from the documentation

An uncompressed binary data file with no header information. The color image data saved with PFV is interleaved (RGBRGBRGB...) and output is supported for each color plane, RGB, at 8-bits, 16-bits, or the sensor's bit depth (color output capable).

[dpb]

OK, that pretty-much confirms that, presuming your camera has 12-bit depth the output is, indeed, packed 12-bit. Looking at one of the CIH files would be interesting exercise.

It looks like it also will save 16-bit files that, while bigger, may well load much faster than translating 12 bits into 16 -- why not start with 16 to begin with????

[Jonas]

actually that was not my decision. I just got those 12bit packed files to work with and need to make the best out of it.

Answer 3

[Image Analyst]

There is also memmapfile - a function older than the datastore family of functions.

Maybe @Steven Lord can explain the differences and when to use each. 🙂🤞