[GRASS-user] bug in r.proj or just artifact of weird (unrealistic) reprojection ?

Glynn Clements glynn at gclements.plus.com
Tue Sep 25 03:16:20 PDT 2012


Moritz Lennert wrote:

> - A location in EPSG 31370 (Belgian Lambert 1972)
> - Objectif: project the elevation raster layer from the nc_spm_08 
> dataset into that Belgian projection dataset
> - A first run of r.proj to get the coordinates of the layer in the new 
> projection system. Here's the result:
> 
> n=2370907.92051779 s=2378797.24912944 w=-6087464.12326059 
> e=-6068211.64408602 rows=1350 cols=1500
> 
> As you can see s > n. g.region bails out on that.
> 
> But, when you just switch n and s, then everything works fine and the 
> elevation dataset gets projected to the correct spot.
> 
> Now, I know that this is not the typical reprojection scenario, but I'm 
> wondering whether this indicates an underlying issue in r.proj's 
> boundary calculation algorithm.

This has nothing to do with r.proj's boundary calculation algorithm,
which isn't used for -p/-g. Instead, it simply projects the south-east
and north-west corners of the input map, and assumes that the
projected points are also the south-east and north-west corners.

Clearly, this will be inaccurate for anything other than projecting
between cylindrical projections and/or lat-lon.

The region-cropping optimisation which is used for projecting (when
the -n isn't given) projects the entire boundary, sampled at the
source resolution, then uses the bounding box of the projected
boundary. AFAICT, this works for all but pathological cases (e.g. 
a projection which is tangential to the region being projected).

Can you try the attached patch?

-- 
Glynn Clements <glynn at gclements.plus.com>

-------------- next part --------------
An embedded and charset-unspecified text was scrubbed...
Name: r.proj.patch
URL: <http://lists.osgeo.org/pipermail/grass-user/attachments/20120925/1f7f26b4/attachment.ksh>


More information about the grass-user mailing list