[Geotiff] Get one part out of the GeoTIFF Image

[Frank Warmerdam]

StefanReisser at t-online.de wrote:
> Hello!
> 
> At first I want to thank everbody for their quick response...
> 
> However, I have another question. The following is the output of the
> gdalinfo tool of both the GeoTIFF and the DEM.
> How can I use the corner coordinates of the DEM to compute the
> corresponding raster coordinates of the GeoTIFF?
> And if the projection parameters are not the same, how can I overcome
> this problem?

Stefan,

First, both images are UTM Zone 10N / NAD27.  That is good.  Otherwise we
would have a much harder problem to resolve.  But now we just need to
concentrate on the following information:

> [root at localhost include]# gdalinfo /mnt/cdrom/L7_P46R26S00.TIF
> Size is 19242, 17442
> Coordinate System is:
...
> Origin = (509375.000000,5522290.000000)
> Pixel Size = (12.50000000,-12.50000000)


> [root at localhost include]# gdalinfo /root/q1826.dem
> Size is 974, 1406
...
> Origin = (585735.000000,5191825.000000)
> Pixel Size = (10.00000000,-10.00000000)

So the top left corner of the DEM is (585735.000000,5191825.000000).
We can translate this to a location on the image file by subtracting the
offset to it's top left corner and then dividing by the pixel size.

So:

pixel_img = (geo_dem - origin_x_img) / x_pixel_size_img
           = (585735 - 509375) / 12.5
           = 6108.8

line_img = (geo_dem - origin_y_img) / y_pixel_size_img
          = (5191825 - 5522290) / -12.5
          = 26437.2

Hmm, this is bad.  It means top left corner of the DEM isn't even on
the image file. In fact, a quick inspection of the corner coordinates of
the two files shows they don't even overlap:

Corner Coordinates:
Upper Left  (  509375.000, 5522290.000)
Lower Left  (  509375.000, 5304265.000)
Upper Right (  749900.000, 5522290.000)
Lower Right (  749900.000, 5304265.000)
Center      (  629637.500, 5413277.500)

Corner Coordinates:
Upper Left  (  585735.000, 5191825.000)
Lower Left  (  585735.000, 5177765.000)
Upper Right (  595475.000, 5191825.000)
Lower Right (  595475.000, 5177765.000)
Center      (  590605.000, 5184795.000)

However, in general, you would use the top left corner offset and
pixel size on each file to convert locations back and forth between
pixel/line and georeferenced.

The offset and pixel size information is accessable via the GDALGetGeoTransform()
function in C and corresponding calls in Python and C++.

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