# Reading periodic elements from a large file

I have a large binary data file (big endian) with 100+ million "rows" of 11 elements, combination of floats and integers.

This is the format:

{"Real32", "Real32", "Real32", "Real32", "Real32", "Real32", "Real32", "Real32", "Real32", "Integer32", "Integer32"}


This question: How to read data file quickly?, is related but not exactly the same.

I've been reading in the whole file like this:

str = OpenRead[filename, BinaryFormat -> True];
data = BinaryReadList[str, {"Real32", "Real32", "Real32", "Real32", "Real32", "Real32", "Real32", "Real32", "Real32", "Integer32", "Integer32"}, ByteOrdering -> +1];


This requires lots and lots of memory and in the end I throw away most of the data most of the time. Usually I am just interested in the 4th real32 and the 2ndint32, or each "row". I would like to read the only the 4th real32 and the 2nd int32 of each "row" if possible and skip over the rest.

I've tried to use Skip but the documentation isn't clear if it works with BinaryReadList. I get the error Skip::readf: Real32 is not a valid format specification. ".

The documentation doesn't describe that you can skip byte by byte, but you can...

str = OpenRead[name, BinaryFormat -> True];
count = FileByteCount[name]/(11*4);
Skip[str, Byte, 24];
{count}]; // AbsoluteTiming


edit: This code works now, but it is very slow, about a minute to load a file that takes only 15 seconds with BinaryReadList, however, the memory overhead is orders of magnitude lower.

edit2: Skip appears to be very slow, much slower than SetStreamPosition for some reason. So I wrote some new code that uses SetStreamPosition with a precomputed list of StreamPositions in bytes. It is about twice as fast as the Skip version, which is okay, but its still about 3x slower than BinaryReadList

pos = Range[12, FileByteCount[name], 11*4];
data = {SetStreamPosition[str, #];
SetStreamPosition[str, # + 28];
BinaryRead[str, "Integer32", ByteOrdering -> +1]} & /@ pos; // AbsoluteTiming


Hopefully, someone will have an idea how this can be improved. Memory usage is still low, as expected.

I'm willing to tolerate a slight slow down (maybe 2x but not 5-10x) if there is a considerable memory savings to be gained but it would be great if the process could be sped up as well.

I can't really easily provide a copy of my data file as they are 100s of megabytes. I tried to write some code that generates some random data and writes it to a file, however, BinaryWrite appears to be extremely slow... I'm on a fast machine with a solid state drive and its going only a few 100 kilobytes per second... Here is the code, regardless, maybe someone knows a faster way to make a random binary data file. This will make an ~40 MB file.

outputstr = OpenWrite["randomdata", BinaryFormat -> True]
reals = RandomReal[100, {10^6, 9}];
ints = RandomInteger[100, {10^6, 2}];
both = Flatten@Transpose@Join[Transpose@reals, Transpose@ints];
BinaryWrite[outputstr, both, {"Real32", "Real32", "Real32", "Real32",
"Real32", "Real32", "Real32", "Real32", "Real32", "Integer32",
"Integer32"}, ByteOrdering -> +1]
Close[outputstr]

-
if anyone wants to try the various approaches or find a better one I'd like to mention that every large file on your computer will serve just as well as the generated file with random content (as long as the reading routines don't fail when the number of entries doesn't exactly match...). Make sure to not overwrite existing files, though... –  Albert Retey Nov 26 '12 at 11:22
@AlbertRetey I'm afraid my code would crash on generic files since I assume that the file size is an integer multiple of the single row byte count, but I do not perform the checks...Not that they are hard to add, but I was not focusing on bullet-proofing the code. –  Leonid Shifrin Nov 26 '12 at 15:55

I was able to get 50x speedup w.r.t. your fastest code by using highly optimized Java buffered read functionality.

### The idea

The idea is quite simple: use buffered read to reduce the IO overhead, and use Java to reduce the symbolic Mathematica overhead.

### Implementation

You will have to run the Java reloader. Then, you call

JCompileLoad@"
import java.io.*;
import java.nio.ByteBuffer;
import java.nio.MappedByteBuffer;
import java.nio.channels.FileChannel;
import java.nio.channels.FileChannel.MapMode;
import java.util.Arrays;

public static int byteArrayToInt(byte[] b){
return   b[3] & 0xFF |
(b[2] & 0xFF) << 8 |
(b[1] & 0xFF) << 16 |
(b[0] & 0xFF) << 24;
}

public  static int[] getIntegerColumn(String filename, int rowByteCount,
int skipBefore, int skipAfter, int rowChunkSize)
throws FileNotFoundException,IOException{
File fl = new File(filename);
FileInputStream str = new FileInputStream(fl);
FileChannel ch = str.getChannel( );
MappedByteBuffer mb = ch.map( FileChannel.MapMode.READ_ONLY, 0L, ch.size( ) );
final int buffrows = rowChunkSize;
final int buffSize = buffrows * rowByteCount;
byte[] buffer = new byte[buffSize];
int rows  = (int)(fl.length()/rowByteCount);
int[] result = new int[rows];
int cycles  = (int)(rows/buffrows);
int remaining = rows % buffrows;
byte[] remBuffer = new byte[remaining * rowByteCount];
int ctr=0;
try{
for(int j=0;j<cycles;j++){
int bctr = 0;
mb.get(buffer);
for(int i=0;i < buffrows;i++){
bctr+=skipBefore;
result[ctr++] = byteArrayToInt(Arrays.copyOfRange(buffer,bctr,bctr+4));
bctr+=4+skipAfter;
}
}
int bctr = 0;
mb.get(remBuffer);
for(int i=0; i < remaining;i++){
bctr+=skipBefore;
result[ctr++] = byteArrayToInt(Arrays.copyOfRange(remBuffer,bctr,bctr+4));
bctr+=4+skipAfter;
}
} finally{
str.close();
}
return result;
}

public  static float[] getFloatColumn(String filename, int rowByteCount, int skipBefore, int skipAfter, int rowChunkSize)
throws FileNotFoundException, IOException{
File fl = new File(filename);
FileInputStream str = new FileInputStream(fl);
FileChannel ch = str.getChannel();
MappedByteBuffer mb = ch.map( FileChannel.MapMode.READ_ONLY, 0L, ch.size( ));
final int buffrows = rowChunkSize;
final int buffSize = buffrows * rowByteCount;
byte[] buffer = new byte[buffSize];
int rows  = (int)(fl.length()/rowByteCount);
float[] result = new float[rows];
byte[] intermediate = new byte[4*rows];
int cycles  = (int)(rows/buffrows);
int remaining = rows % buffrows;
byte[] remBuffer = new byte[remaining * rowByteCount];
int ctr=0;
try{
for(int j=0;j<cycles;j++){
int bctr = 0;
mb.get(buffer);
for(int i=0;i < buffrows;i++){
bctr+=skipBefore;
System.arraycopy(buffer, bctr,intermediate,4*ctr++,4);
bctr+=4+skipAfter;
}
}
int bctr = 0;
mb.get(remBuffer);
for(int i=0; i < remaining;i++){
bctr+=skipBefore;
System.arraycopy(remBuffer, bctr,intermediate,4*ctr++,4);
bctr+=4+skipAfter;
}
ByteBuffer buf2 = ByteBuffer.wrap(intermediate);
for(int i=0;i<rows;i++){
result[i]=buf2.getFloat();
}
} finally{
str.close();
}
return result;
}
}"


### Usage

There are 2 static methods, to extract a single column, of integer or floating point numbers. Both take the same set of 5 parameters: file name, total bytes in one row, bytes to skip before reading the element in one row, bytes to skip after, and the number of rows in a buffer for buffered read.

### Benchmarks

Using your code to produce the 40Mb file, I get then:

(jdataInt = TableReadergetIntegerColumn[name,11*4,10*4,0,100])
//Length//AbsoluteTiming

(*  {0.0898438,1000000}  *)

//Length//AbsoluteTiming

(* {0.0839844,1000000} *)


while your code on my machine gives

str = OpenRead[name, BinaryFormat -> True];
pos = Range[12, FileByteCount[name], 11*4];
data = {
SetStreamPosition[str, #];
SetStreamPosition[str, # + 28];
} & /@  pos; // AbsoluteTiming
Close[str];

(* {9.1044922,Null} *)


And we can verify:

Flatten[data[[All,1]]] == jdataFl

(* True *)

Flatten[data[[All,2]]]==jdataInt

(* True *)


### Conclusions

By using buffered reads implemented in a compiled language with a highly optimized native IO routines (the latter implemented by much brighter people than me :)), I was able to gain two orders of magnitude speedup. I suspect that the running time for the Java code is dominated by data transfer, so the read itself is quite a bit faster still. I also think that going to C one can gain further significant speedups, although probably not as dramatic as here.

Note that I was quite sloppy with the Java-side error-handling, partly intentionally to optimize for speed, partly because at this point I did not really care. If someone decides to use this in production code, more care must be taken however.

-
Sorry, if this is silly, but where is the endianness specified? –  s0rce Nov 26 '12 at 23:56
@s0rce It is big endian, which is implicit in both the code for byteArrayToInt method and the getFloat() method of the ByteBuffer class that I use for the floats. Thanks for the accept,b.t.w. I hope this will help you. For really large files, some changes may be needed - one may have to launch JRE (Java) with more memory. –  Leonid Shifrin Nov 27 '12 at 9:44
@s0rce But why do you need to modify it at all? It has already a method getFloatColumn, and above there is an example of how to use it. All you have to do is to compute the offsets and total byte count of a row - and then just use it. –  Leonid Shifrin Mar 6 '13 at 22:23
Thanks for your quick response. I just realized I'm dumb. I wanted to read the whole file but I thought about it wrong and I put in a value for the rowSize but really that should just be 4 with 0 offset before and after, obviously. However, I'm still getting memory errors on OSX, might have to do with how much memory java is allocated? –  s0rce Mar 6 '13 at 22:31

Leonid has shown how to get a big improvement in speed with Java and has mentioned that by using a library function written e.g. in C you could probably get further improvement as you can then avoid additional copies/data transfer.

I just wanted to add an all Mathematica solution as it probably indicates how similar problems can be tackled where for whatever reason something like the Java approach isn't possible or necessary. I've used a similar approach to e.g. read only small parts of huge XML-files. This isn't faster as the memory intensive read-all-at-once but at least isn't slower und is much more memory efficient.

It basically reads the file in chunks of adjustable size, extracts the part of interest of each chunk and collects everything. Here is the code:

AbsoluteTiming[
str = OpenRead[filename, BinaryFormat -> True];
chunksize = 500;
data = {};
res = {1};
While[Or[data == {},
And[res[[-1, -1]] =!= EndOfFile, Length[res] > 0]],
data = {data,
res =
str, {"Real32", "Real32", "Real32", "Real32", "Real32",
"Real32", "Real32", "Real32", "Real32", "Integer32",
"Integer32"}, chunksize, ByteOrdering -> +1][[All, {4, 11}]]
};
];
Close[str];
data = Partition[Flatten[data], 2];
]


A slightly more elegant way using Sow and Reap with localized variables and packed as a function is this:

readFile[filename_, chunksize_Integer: 500] := Module[{
str = OpenRead[filename, BinaryFormat -> True],
res = {{1}},
data
},
data = Reap[
While[And[res[[-1, -1]] =!= EndOfFile, Length[res] > 0],
Sow[res =
str, {"Real32", "Real32", "Real32", "Real32", "Real32",
"Real32", "Real32", "Real32", "Real32", "Integer32",
"Integer32"}, chunksize, ByteOrdering -> +1][[
All, {4, 11}]]];
]
][[2]];
Close[str];
Flatten[data[[1]], 1]
]


I did expect this to use about twice the memory of the final result, as we need one copy and depending on the chunksize the additional memory for the intermediate data (that is all columns except for 4 and 11) which is "thrown away" should become neglectable for small chunk sizes.

I have now done some measuring and one interesting thing is that memory usage doesn't necessarily become much larger than the final total datasize, which indicates that Mathematica doesn't make a full copy for the rearranged final result. Another interesting feature is that the memory usage of this approach for small chunksizes seems to decrease with increasing chunksizes. Only for very large chunk sizes (larger than about 100 000 for this example with a 600MB file on my computer) the behavior is "as expected": the memory consumption will then increase as the chunk size increases. My interpretation of this is that Mathematica is doing some internal optimization/buffering for BinaryReadList which mainly influences the situation for small chunk sizes. The result is that very small chunksizes are slow and use more memory than larger values, presumably a side effect of the large number of iterations and thus evaluations of the loop-body. Then for junksizes larger than about 50 neither runtime nor memory usage seem to vary very much until chunk sizes get relatively large, presumably comparable to the internal buffer sizes. The precise explanation for all details of this behaviour is beyond my current knowledge, all in all I see a great decrease in memory usage for this approach, but max. memory usage for large files is always in the order of sizeOfFinalResult + O[100MB] for medium sized chunks and can be kept at roughly that value without much cost in run time even for very large files. With this approach it is possible to read files too large to be read all-at-once. Your mileage may vary depending on data, chunksize and hardware, so some testing-measuring when using this approach is certainly a good idea...

-
+1. I was thinking in this direction as well, actually my Java code implements the same idea. But since I started with Java, I wanted to get as fast as I could :) –  Leonid Shifrin Nov 26 '12 at 11:06
@LeonidShifrin: your answer of course makes mine look somewhat lame, but I was almost having it ready when I saw yours so I thought I'd post it anyway :-). I think it's at least an alternative which might be easier to understand and adopt for some of us not-so-java-fluents –  Albert Retey Nov 26 '12 at 11:12
Actually having an all-Mathematica solution is very good. Besides, I am not sure how memory-efficient my solution is. It should be ok, but I did not really measure that. It may very well be that yours is better in this respect. I personally was more interested in pushing the speed envelope. Once I tried RandomAccessFile and FileInputStream vs BinaryReadList and found that they are on par. Now I know that the buffered read via Java can still offer a speed advantage. And, another use case for my Java reloader :) –  Leonid Shifrin Nov 26 '12 at 11:16

You write

This requires lots and lots of memory and in the end I throw away most of the data most of the time.

If your only concern is the memory which is used by Mathematica after the values are extracted from the file, then there is a simple solution which I describe in the following. If your concern is the memory usage during reading the file, which means, you have not enough RAM to read the whole file even temporarily in memory, than you should look at the solution of Leonid.

The only thing which is necessary is, that you don't store the whole read data, but you access the wanted columns instantly. Then you can use BinaryReadList as expected and you don't waste memory:

name = "tmp/randomdata.dat";
MemoryInUse[]/2.^20


@s0rce Hmm, that wasn't clear from the question. I suspected, that (since you stored data completely) you were worried about the memory used by the kernel after you have read all data. –  halirutan Nov 26 '12 at 15:36