postgis - Understand gdalwarp and gdalinfo on raster files

I am using:

• PostgreSQL 9.5


• PostGIS 2.2.1


• GDAL 1.11.3

I have two raster files:

-rw-rw-r-- 1 jlandercy jlandercy 1,3G Aug 25 21:02 contour_raster.dem
-rw-rw-r-- 1 jlandercy jlandercy 2,0G Sep  7 17:42 contours.dem

The second is the result of:

gdalwarp contour_raster.dem contours.dem

I would like to understand what it has done, files are different:

• in size;



• in structure;


• and behave differently when I import them into PostgreSQL.

If I compare their gdalinfo using diff:

2,3c2,3
< Files: contour_raster.dem
< Size is 22807, 22715
---
> Files: contours.dem
> Size is 22810, 22718

24,26c24,27
< GeoTransform =
<   138103.531187929, 1.000005446, -0.000117871259
<   181238.971718448, -0.000116195632, -1.000007394
---
> Origin = (138100.853742280829465,181238.971718448010506)
> Pixel Size = (1.000000000000000,-1.000000000000000)
> Metadata:
>   AREA_OR_POINT=Area
28d28

<   COMPRESSION=LZW
31,36c31,36
< Upper Left  (  138103.531,  181238.972) (  4d11'53.62"E, 50d56'30.71"N)
< Lower Left  (  138100.854,  158523.804) (  4d11'56.13"E, 50d44'15.61"N)
< Upper Right (  160910.655,  181236.322) (  4d31'21.79"E, 50d56'30.69"N)
< Lower Right (  160907.978,  158521.154) (  4d31'19.23"E, 50d44'15.60"N)
< Center      (  149505.755,  169880.063) (  4d21'37.69"E, 50d50'23.56"N)
< Band 1 Block=128x128 Type=Float32, ColorInterp=Gray
---
> Upper Left  (  138100.854,  181238.972) (  4d11'53.49"E, 50d56'30.71"N)

> Lower Left  (  138100.854,  158520.972) (  4d11'56.13"E, 50d44'15.52"N)
> Upper Right (  160910.854,  181238.972) (  4d31'21.80"E, 50d56'30.78"N)
> Lower Right (  160910.854,  158520.972) (  4d31'19.38"E, 50d44'15.59"N)
> Center      (  149505.854,  169879.972) (  4d21'37.70"E, 50d50'23.56"N)
> Band 1 Block=22810x1 Type=Float32, ColorInterp=Gray

The original raster has GeoTransform the other is in pixel (I don't know exactly what it means). Bounding Boxes are different, and Band is not encoded in the same way.

Further more if I use raster2pgpsql -C -t 1000 1000 to square rasters and insert tiles, the first file lacks constraints that cannot be set:

• enforce_same_alignment_rast;



• enforce_scalex_rast;


• enforce_scalex_rast.

My questions are:

• Why files do not have the same size (is it because of COMPRESSION=LZW)?


• What is the difference between the two formats (alignment problem due to coordinate system)?


• Why bounding boxes are different?


• Why constraints fail in the first case?

Answer

The chapter about Affine GeoTransform in http://www.gdal.org/gdal_datamodel.html gives the first hint

GDAL datasets have two ways of describing the relationship between raster positions (in pixel/line coordinates) and georeferenced coordinates. The first, and most commonly used is the affine transform (the other is GCPs). The affine transform consists of six coefficients returned by GDALDataset::GetGeoTransform() which map pixel/line coordinates into georeferenced space using the following relationship:

Xgeo = GT(0) + Xpixel*GT(1) + Yline*GT(2)
Ygeo = GT(3) + Xpixel*GT(4) + Yline*GT(5)

In case of north up images, the GT(2) and GT(4) coefficients are zero, and the GT(1) is pixel width, and GT(5) is pixel height. The (GT(0),GT(3)) position is the top left corner of the top left pixel of the raster.

Note that the pixel/line coordinates in the above are from (0.0,0.0) at the top left corner of the top left pixel to (width_in_pixels,height_in_pixels) at the bottom right corner of the bottom right pixel. The pixel/line location of the center of the top left pixel would therefore be (0.5,0.5).

From http://www.gdal.org/gdal_tutorial.html you can read also

In the particular, but common, case of a "north up" image without any rotation or shearing, the georeferencing transform takes the following form :

adfGeoTransform[0] /* top left x */

adfGeoTransform[1] /* w-e pixel resolution */

adfGeoTransform[2] /* 0 */

adfGeoTransform[3] /* top left y */

adfGeoTransform[4] /* 0 */

adfGeoTransform[5] /* n-s pixel resolution (negative value) *

Your first image is not of this particular, but commom, sort. The image is slightly rotated and therefore it needs all the six parameters.

What gdalwarp does is that it rotates the image so that it becomes a north-up image. After this process four parameters are enough for defining the geotransformation. Knowing pixel width, height, and pixel count from left to right and top to bottom are enough for georeferencing the image.

Rotation that gdalwarp does makes the image wider and higher. Your image is almost north-up even before warping and therefore the difference is just three pixels in both directions.

Difference in file size is mostly due to LZW compression. As uncompressed the warped image is only slightly bigger because of the three additional rows and columns.

Raster2pgsql seems to require that images which have the simple geotransformation and also that pixels are square with same width and height.