[GRASS-user] r.proj - Input raster map is outside current region

[Hamish]

Kirk wrote:
> I am trying to re-project a raster
> into a location, and I am getting an error about the region
> of the raster being outside the region of the new location.
> 
> 
> r.proj input=foresttype location=mnwimi_usfs
> mapset=PERMANENT output=foresttype
> Input Projection Parameters:  +proj=aea +lat_1=29.5
> +lat_2=45.5 +lat_0=23 +lon_0=-96 +x_0=0 +y_0=0 +no_defs
> +a=6378137 +rf=298.257222101 +towgs84=0,0,0,0,0,0,0
> Input Unit Factor: 1
> Output Projection Parameters:  +proj=aea
> +lat_1=42.122774 +lat_2=49.01518 +lat_0=45.568977
> +lon_0=-84.455955 +x_0=1000000 +y_0=1000000 +no_defs
> +a=6378137 +rf=298.257222101 +towgs84=0.000,0.000,0.000
> Output Unit Factor: 1
> ERROR: Input raster map is outside current region
> 
> However, running g.region on both the raster to be
> reprojected and the new locations shows the regions to be
> the same:
> 
> GRASS 6.4.1 (greatlakesalbers):~ > g.region -p
> projection: 99 (Albers Equal Area)
> zone:       0
> datum:      nad83
> ellipsoid:  grs80
> north:      2960848
> south:      2087348
> west:       -122110
> east:       1141890
> nsres:      250
> ewres:      250
> rows:       3494
> cols:       5056
> cells:      17665664
> 
> GRASS 6.4.1 (mnwimi_usfs):~ > g.region -p
> rast=foresttype
> projection: 99 (Albers Equal Area)
> zone:       0
> datum:      nad83
> ellipsoid:  grs80
> north:      2960848
> south:      2087348
> west:       -122110
> east:       1141890
> nsres:      250
> ewres:      250
> rows:       3494
> cols:       5056
> cells:      17665664
> 
> I am most perplexed by the differences in the lat_1 lat_2,
> etc parameters reported by r.proj, and the region values
> reported by g.region. Any suggestions would be much
> appreciated.

if lat_1 and lat_2 are different, then they won't have the same bounds,
just that the same family of math formulas will be used with different
coefficients. (so you get different answers)  I notice the false eastings
and northings are also different in the above.

With recent GRASS versions the method to use is to first run 'r.proj -p'
from the target location to show the bounds of source map reprojected into
the current (target) CRS. If that looks ok, re-run with 'r.proj -g' to get
the same formatted in a string you can feed directly to g.region in the
target location.

once that is done run 'g.region -p' and chose an appropriate resolution
(often a little rounding is needed). then align the region to that
chosen resolution like:
  'g.region -a res=250 -p'

only once you are happy with that, then run 'r.proj' without those flags
to do the actual map reprojection. The resolution must be set manually,
and the number of rows and columns should be reviewed by hand too (for maps
which will be rotated up to 45 deg you'll want to increase that from the
original) so it is not appropriate to have that all happen automatically
for you as it leads to enhanced data/generational loss. There's a bug
report for i.rectify in the old RT bug tracker which explains the issue
of automatic bounds setting in a bit more detail if you're interested.


Hamish