[GRASSLIST:5855] Re: Question regarding r.proj

[Glynn Clements]

Victor Wren wrote:

> I need a suggestion how to proceed.  I've got a large number of raster maps 
> that I need to project from NAD27 datum to NAD83 (both in UTM, grs80 
> ellipsoid, though I think NAD27 is supposed to be clark66, but corpscon 
> doesn't support that ellipsoid, so I couldn't check the coordinate shift 
> there).  I wrote a shell script to do this, but it failed dismally, because I 
> didn't have the correct region selected for each raster at the time of 
> import.
> 
> Given that I won't know the extents of a raster in the current map until I 
> import it, it seems a chicken-and-egg problem to decide what the current 
> region should be set to before importing.  If the region is too large, the 
> raster will be expanded on import to its limits (requiring a huge amount of 
> memory and conversion time).  If it's too small, the maps will end up being 
> trimmed (or will be skipped entirely, if they're outside the region.)  I 
> could set them huge, and r.patch them after, but that seems inelegant (and 
> I'm not sure how to automate that, either).

This shouldn't be the case; r.proj includes code to project the
borders of both the current region and the source raster to determine
the final borders.

What should happen is that r.proj projects the current region to the
source location's projection to determine the region of the source map
which needs to be read, then projects that region back to the
destination location to determine the extents of the resulting map.

> m.proj does not seem to support datum conversion (and m.datum.shift doesn't 
> support UTM coordinates), so I can't readily use that to find out where the 
> projected raster extents will end up.  Experiments with corpscon show that in 
> my area of interest, the raster will shift about 200m Northing and about 50m 
> Easting.

AFAICT, the CVS version of m.proj2 should handle datum conversion.

> A switch in [rsv].proj to override the current region and create a new raster 
> with the same geometric extents (projected) as the source raster would be a 
> real timesaver, but I'm not enough of a code hog to even know whether that 
> can be done at the level that r.proj works at.

Provided that the current region is initially larger than the
resulting map (and you don't use the -n flag), you should end up with
a map whose extents are "shrink wrapped" to the actual data. If that
doesn't happen, it's a bug.

> I'm using Grass 5.0.2.

5.0.2 hasn't been released yet; it's either 5.0.1, or one of the CVS
versions (which are currently named "5.0.2-cvs").

-- 
Glynn Clements <glynn.clements at virgin.net>