geotiff tiff - Read Elevation from 1/3 Arc-Second NED .tif or .img

I'm trying to read multiple elevations from a USGS 1/3 arc-second IMG / GeoTIFF file.

A sample of the files is available here (5 MB): http://tdds3.cr.usgs.gov/Ortho9/ned/ned_13/img/n43w108.zip

I can query the .img file and get the elevation using GDAL: gdallocationinfo -wgs84 imgn60w145_13.img -144.267361111152354 59.9981944444444082

This returns: Report: Location: (7918P,25L) Band 1: Value: 1.2593731880188

Because a GeoTIFF can be read quite a bit faster than an IMG file (in my benchmarks) I convert the IMG to a WGS84 GeoTIFF: gdal_translate -a_srs WGS84 -of GTiff imgn60w145_13.img imgn60w145_13.tif

gdallocationinfo only supports querying one point at a time. I need to get the elevation every 10 meters between 2 points (along a straight path), so I'm trying to read the file directly, with the following PHP code (but the approach is the same in C or Python):

$STRIPOFFSETS = 230; // 0xe6 $LEN_OFFSET = 4; $DATA_VOID = 0x8000; // data void ( = signed int -32768) $LEN_DATA = 4; // the number of bytes containing each item of elevation data // ( = BitsPerSample tag value / 8) $fp = fopen("imgn60w145_13.tif", 'rb'); if ($fp === false) { echo "Could not open the file\n"; } // first data offset fseek($fp, $STRIPOFFSETS); // find the location of the required data row in the StripOffsets data $dataOffset = $STRIPOFFSETS + ($row * $LEN_OFFSET); fseek($fp, $dataOffset); $dataBytes = fread($fp, $LEN_OFFSET); $data = unpack('VdataOffset', $dataBytes); echo print_r($data, true); // this is the offset of the 1st column in the required data row $firstColOffset = $data['dataOffset']; // now work out the required column offset relative to the 1st column $requiredColOffset = $col * $LEN_DATA; // combine the two and read the elevation data at that address fseek($fp, $firstColOffset + $requiredColOffset); $dataBytes = fread($fp, $LEN_DATA); echo "1: " . $firstColOffset . " + " . $requiredColOffset . "\n"; echo "2: " . $LEN_DATA . "\n"; echo "3: " . $dataBytes . "\n"; $data = unpack('velevation', $dataBytes); $elevation = $data['elevation']; if ($elevation == $DATA_VOID) { $elevation = 0; } echo $elevation . "\n";

I don't know that I have the right values for $STRIPOFFSETS, $LEN_OFFSET, or $LEN_DATA. The code above is a limited snippet from SRTMGeoTIFFReader (http://www.osola.org.uk/elevations/index.htm), which I'm not able to get to work either.

I'd really prefer to be able to do this using a direct approach similar to the above (no external libraries, etc), but I can live with a Python / GDAL approach as well.

Thanks!

Answer

The solution for me was to convert the files to BIL files using GDAL:

gdal_translate -a_srs NAD83 -of EHdr file.img file.bil

This created .bil and .hdr files. The .hdr file contains the details you need to determine where to read the file. I was then able to get the elevation using the PHP unpack function. Performance is lightning fast.