[Gdal-dev] Find elevation data at specific point

[Frank Warmerdam]

Gerald Brandt wrote:
> Thanks!  What would the pixelx, liney formula be if the image is not
> north up?  
> 
> I'm working with DEM and DTED and have not yet seen anything that isn't
> north up, but that doesn't mean there aren't any!

Gerald,

Ummm, well it gets more complicated.

I normally invert the geotransform matrix, and then use it in format fashion.

The following code will invert the geotransform[]:

/************************************************************************/
/*                          InvGeoTransform()                           */
/*                                                                      */
/*      Invert a standard 3x2 "GeoTransform" style matrix with an       */
/*      implicit [1 0 0] final row.                                     */
/************************************************************************/

int InvGeoTransform( double *gt_in, double *gt_out )

{
     double	det, inv_det;

     /* we assume a 3rd row that is [1 0 0] */

     /* Compute determinate */

     det = gt_in[1] * gt_in[5] - gt_in[2] * gt_in[4];

     if( fabs(det) < 0.000000000000001 )
         return 0;

     inv_det = 1.0 / det;

     /* compute adjoint, and devide by determinate */

     gt_out[1] =  gt_in[5] * inv_det;
     gt_out[4] = -gt_in[4] * inv_det;

     gt_out[2] = -gt_in[2] * inv_det;
     gt_out[5] =  gt_in[1] * inv_det;

     gt_out[0] = ( gt_in[2] * gt_in[3] - gt_in[0] * gt_in[5]) * inv_det;
     gt_out[3] = (-gt_in[1] * gt_in[3] + gt_in[0] * gt_in[4]) * inv_det;

     return 1;
}

Then you could compute the pixel/line from geo_x, geo_y like this:

pixel = gt[0] + geo_x * gt[1] + geo_y * gt[2]
line = gt[3] + geo_x * gt[4] + geo_y * gt[5]

Best regards,
-- 
---------------------------------------+--------------------------------------
I set the clouds in motion - turn up   | Frank Warmerdam, warmerdam at pobox.com
light and sound - activate the windows | http://pobox.com/~warmerdam
and watch the world go round - Rush    | Geospatial Programmer for Rent