Opened 3 years ago

Closed 3 years ago

Last modified 3 years ago

#2911 closed defect (fixed)

[raster] ST_Rescale changes extent

Reported by: strk Owned by: strk
Priority: blocker Milestone: PostGIS 2.1.4
Component: raster Version: 2.1.x
Keywords: history Cc:

Description

According to documentation ( http://postgis.net/docs/RT_ST_Rescale.html ), ST_Rescale should only change the pixel scale, by retaining (or at most slightly enlarging) the extent of a raster.

While it mostly work that way, I've found this input to be handled differently, by retaining the same scale and halving the extent instead (!!)

0100000100000000000000F03F000000000000F0BF0000000000000000000000000000000000000000000000000000000000000000000000000A000A00640000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000

Its metadata:

upperleftx | 0
upperlefty | 0
width      | 10
height     | 10
scalex     | 1
scaley     | -1
skewx      | 0
skewy      | 0
srid       | 0
numbands   | 1

Can be created with:

SELECT ST_AddBand(
 ST_MakeEmptyRaster(10, 10, 0, 0, 1, -1, 0, 0, 0),
 1, '8BUI', 0, 0
);                                                

Attachments (1)

warp_with_default_geotransform.patch (1.8 KB) - added by strk 3 years ago.
I'm attaching the current patch I wrote to fix the reported case, but it needs more love

Download all attachments as: .zip

Change History (17)

comment:1 Changed 3 years ago by strk

In other words, ST_Rescale does not change the scale, but just the size, although the parameters are floating point numbers. Noticeable are the "default geotransform" notices:

pg21=# with data as ( 
   SELECT ST_AddBand(
          ST_MakeEmptyRaster(10, 10, 0, 0, 1, -1, 0, 0, 0),
          1, '8BUI', 0, 0
        ) r 
 ), 
 meta as ( 
   SELECT st_metadata(st_rescale(r,2.0,-2.0)) m from data 
 ) select (m).* from meta;
NOTICE:  Raster has default geotransform. Adjusting metadata for use of GDAL Warp API
CONTEXT:  SQL function "st_rescale" statement 1
-[ RECORD 1 ]--
upperleftx | 0
upperlefty | 0
width      | 5
height     | 5
scalex     | 1
scaley     | -1
skewx      | 0
skewy      | 0
srid       | 0
numbands   | 1

comment:2 Changed 3 years ago by strk

Another view on the matter, showing how origin 0,scaley*height is threated specially:

pg21=# select x, y, m.scalex, m.scaley from generate_series(-1, 1) x, generate_series(-1, 1) y, ST_AddBand(ST_MakeEmptyRaster(10, 10, x, y, 1, -1, 0, 0, 0), 1, '8BUI', 0, 0) r, st_metadata(st_rescale(r,2.0,-2.0)) m;
NOTICE:  Raster has default geotransform. Adjusting metadata for use of GDAL Warp API
CONTEXT:  SQL function "st_rescale" statement 1
 x  | y  | scalex | scaley 
----+----+--------+--------
 -1 | -1 |      2 |     -2
  0 | -1 |      2 |     -2
  1 | -1 |      2 |     -2
 -1 |  0 |      2 |     -2
  0 |  0 |      1 |     -1
  1 |  0 |      2 |     -2
 -1 |  1 |      2 |     -2
  0 |  1 |      2 |     -2
  1 |  1 |      2 |     -2
(9 rows)

comment:3 Changed 3 years ago by strk

r11215 may be the cause of this

comment:4 Changed 3 years ago by strk

Owner: changed from dustymugs to strk
Status: newassigned

I confirm the problem is with re-setting geotransform to the original one as the last step of warping when the default geotransform was changed to an arbitrary projection.

If a scale change was requested, the reset should take into account that change. This is likely to also affect rotation...

comment:5 Changed 3 years ago by strk

I've a fix for this, but it breaks an existing test which I belive is broken, the one introduced here: http://trac.osgeo.org/postgis/changeset/11215

That test pretends that the output scale of a ST_Rescale with parameter 2,2 is 1, -1 as the input. Shouldn't it be 2, -2 instead ?

Changed 3 years ago by strk

I'm attaching the current patch I wrote to fix the reported case, but it needs more love

comment:6 Changed 3 years ago by strk

r12928 fixes this in trunk. It does not make any attempt at handling translation or rotation of so called "default georeferenced rasters", is that supposed to be possible with _st_gdalwarp ? I'd guess so. It'd be better verified before backporting.

comment:7 Changed 3 years ago by strk

Priority: highblocker

comment:8 Changed 3 years ago by strk

Warping with a rotation still results in wrong output:

select x, y, m.width, m.height, m.skewx, m.skewy from generate_series(-1, 1) x, generate_series(-1, 1) y, ST_AddBand(ST_MakeEmptyRaster(10, 10, x, y, 1, -1, 0, 0, 0), 1, '8BUI', 0, 0) r, st_metadata(_st_gdalwarp(r,'NearestNeighbor', 0.125,NULL,NULL,NULL,NULL,NULL,1,1)) m;
 x  | y  | width | height | skewx | skewy 
----+----+-------+--------+-------+-------
 -1 | -1 |    10 |     15 |     1 |     1
  0 | -1 |    10 |     15 |     1 |     1
  1 | -1 |    10 |     15 |     1 |     1
 -1 |  0 |    10 |     15 |     1 |     1
  0 |  0 |     1 |     12 |     1 |     1
  1 |  0 |    10 |     15 |     1 |     1
 -1 |  1 |    10 |     15 |     1 |     1
  0 |  1 |    10 |     15 |     1 |     1
  1 |  1 |    10 |     15 |     1 |     1
(9 rows)

Trying to fix that I got the failures reported in #2912

comment:9 Changed 3 years ago by strk

Ok #2912 has nothing to do with the changes in output skew computation, so I'm continuing with tests of the fix for that.

comment:10 Changed 3 years ago by strk

For the record, this bug affects any raster having georeference of ul:0,0 scale:1,-1 skew:0,0 no matter having or not a SRID. I mention it because Bborie told on IRC that only SRID-less rasters were affected instead (that's not the case).

comment:11 Changed 3 years ago by strk

NOTE: a webmercator tile for a raster of 1 meter per pixel may very well be common with upper-left corner 0,0

comment:12 Changed 3 years ago by strk

I've noticed that ST_Width(rast) returns 10 and ST_Height(rast) returns 15, for all tiles (sounds good). So what is ST_Metadata.width and ST_Metadata.height ? Why are them different ?

comment:13 Changed 3 years ago by strk

Please forget teh ST_Width vs. ST_Metadata.Width comment, it was my mistake, they do both return the same (bogus) answer for the warped raster.

comment:14 Changed 3 years ago by strk

Summary: ST_Rescale changes extent[raster] ST_Rescale changes extent

comment:15 Changed 3 years ago by strk

Resolution: fixed
Status: assignedclosed

Well I backported r12928 to the 2.1 branch with r12951 and I'm closing this. Will file another one for the skew issue.

comment:16 Changed 3 years ago by robe

Keywords: history added
Note: See TracTickets for help on using tickets.