[GRASS-dev] shift raster when projecting

[Paulo van Breugel]

Dear all,

I am trying to project a DSM raster layer from Belgium (
https://overheid.vlaanderen.be/en/producten-diensten/digital-elevation-model;
EPSG 31370 to EPSG 28992 (Amersfoort / RD New).

The result is a DEM raster layer with the features shifted roughly 100 m as
compared to the same features in the AHN (Dutch DEM,
https://www.ahn.nl/ahn-4). See screenshot and code below.

Reprojecting the same in QGIS using the warp function (gdalwarp) results in
the Belgium DSM layer to be fairly closely aligned to the Dutch one.

I tried to set the bounds before, align to the target raster, and other
suggestions in the help file of r.proj, to no avail. Any idea what I might
be doing wrong?

Cheers,

Paulo

[image: image795.png]

AHN (Netherlands)

EPSG 28992
projection: 99 (Amersfoort / RD New)
zone:       0
datum:      rd18
ellipsoid:  bessel
north:      397562.15979713
south:      373497.15979713
west:       107708.49748491
east:       133913.49748491
nsres:      5
ewres:      5
rows:       4813
cols:       5241
cells:      25224933

Belgium DSM

EPSG 31370
projection: 99 (BD72 / Belgian Lambert 72)
zone:       0
datum:      bel72
ellipsoid:  international
north:      250000
south:      238000
west:       162000
east:       194000
nsres:      5
ewres:      5
rows:       2400
cols:       6400
cells:      15360000


Reproject

(Thu Jun 15 17:11:08 2023)

r.proj location=Merkske_Belgium mapset=PERMANENT input=DHMVIIDSMRAS5m_k02
output=test memory=3000 resolution=5
Selected PROJ pipeline:
+proj=pipeline +step +inv +proj=sterea +lat_0=52.1561605555556
+lon_0=5.38763888888889 +k=0.9999079 +x_0=155000 +y_0=463000 +ellps=bessel
+step +proj=lcc +lat_0=90 +lon_0=4.36748666666667 +lat_1=51.1666672333333
+lat_2=49.8333339 +x_0=150000.013 +y_0=5400088.438 +ellps=intl
************************
Selected PROJ pipeline:
+proj=pipeline +step +inv +proj=lcc +lat_0=90 +lon_0=4.36748666666667
+lat_1=51.1666672333333 +lat_2=49.8333339 +x_0=150000.013 +y_0=5400088.438
+ellps=intl +step +proj=sterea +lat_0=52.1561605555556
+lon_0=5.38763888888889 +k=0.9999079 +x_0=155000 +y_0=463000 +ellps=bessel
************************

Input:
Cols: 4114 (original: 6400)
Rows: 2400 (original: 2400)
North: 250000.000000 (original: 250000.000000)
South: 238000.000000 (original: 238000.000000)
West: 173430.000000 (original: 162000.000000)
East: 194000.000000 (original: 194000.000000)
EW-res: 5.000000
NS-res: 5.000000

Output:
Cols: 4108 (original: 5241)
Rows: 2456 (original: 4813)
North: 396852.159797 (original: 397562.159797)
South: 384572.159797 (original: 373497.159797)
West: 107708.497485 (original: 107708.497485)
East: 128248.497485 (original: 133913.497485)
EW-res: 5.000000
NS-res: 5.000000

Allocating memory and reading input raster map...
Selected PROJ pipeline:
+proj=pipeline +step +inv +proj=sterea +lat_0=52.1561605555556
+lon_0=5.38763888888889 +k=0.9999079 +x_0=155000 +y_0=463000 +ellps=bessel
+step +proj=lcc +lat_0=90 +lon_0=4.36748666666667 +lat_1=51.1666672333333
+lat_2=49.8333339 +x_0=150000.013 +y_0=5400088.438 +ellps=intl
************************
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/grass-dev/attachments/20230616/1dff220e/attachment-0001.htm>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: image795.png
Type: image/png
Size: 351793 bytes
Desc: not available
URL: <http://lists.osgeo.org/pipermail/grass-dev/attachments/20230616/1dff220e/attachment-0001.png>