How do I convert coordinates to raster (bsq, bil) lookup?

I have a .bil file, which has interleave = bsq in the hdr file.

I'm having trouble finding any insight into how I would essentially do a lookup on the file.

I have no problem reading a specific chunk of bytes randomly from the file, but am unsure of how I would go about converting a lat/lon pair to a specific offset within the file. I've only seen solutions in R, but ideally I'd like to be able to do this in python or javascript, using proj4 or something to that effect.

Any advice?

Answer

BSQ is band sequential. The GDAL doc is here

Extract from ESRI whitepaper: Band sequential format stores information for the image one band at a time. In other words, data for all pixels for band one is stored first, then data for all pixels for band two, and so on. in http://downloads.esri.com/support/whitepapers/other_/eximgav.pdf on page 8.

If there's only one band there is no difference between BIL and BSQ.

If the NBITS isn't specified then it can be calculated by the size of file (in bytes) / (NROWS * NCOLS). This should give you either 8, 16 or 32. If it's 32 it may be float or 32 bit integer (long), you will need this in bytes 1, 2 or 4 respectively.

The hook point is usually upper left, so you find the difference from your coordinate to the reference point (in cells):

Row = int( (ULYMAP - Your.Y) / YDIM )
Col = int( (Your.X - ULXMAP) / XDIM )

That is assuming your point is in the same spatial reference as the image, if not project the point - or if all of your requests are going to be in the same spatial reference that's different to the raster consider projecting the raster.

With the row and column known calculating the location in file is (Col * NCOLS * bytesPerCell) + Row * bytesPerCell