[GRASS5] raster inaccuracy

[Roger Miller]

On Wed, 10 Nov 2004 14:49:24 +0000, glynn wrote

> First, I suspect that any discrepancy between vectors and rasters is
> more likely to be in the vector code.

The output from d.where, d.what.vect and d.what.rast consistently show that
where there is a difference the vector map is in about the right place and the
raster map is in the wrong place.
 
> Also, bear in mind that GRASS' raster coordinate system is
> "grid-registered" rather than "node-registered". E.g. the north-east
> edge of the region is the location of the north-east corner of the
> north-east cell, not of the centre of the north-east cell. If you are
> given the region coordinates using a node-registered convention, you
> have to expand the region by half a cell in each direction for use
> with GRASS.

This difference doesn't effect my current problem, but I'm aware that it does
need to be kept in mind.

> For similar reasons, display modules which draw lines tend to be
> problematic, as the display architecture uses integer coordinates, so
> you can only place vertices on the corners of cells, not on their
> centres.

The cells in this case are thousands of units across, so the difference
between integer and floating point coordinates should be negligible.

> In any case, in order to look into this any further, I would need
> specific cases, i.e. a specific sequence of commands which can be
> performed on data which I either have (e.g. spearfish, global) or can
> generate.

I put together several examples that can be had from

http://members.spinn.net/~roger/example.tgz

The tar file produces a GRASS XY location named "example."  The examples are
in the PERMANENT map set.

Here is the step-by-step procedure:

copy the tar file to our Grass data file hierarchy and unarchive the file

start Grass, using location "example" and mapset "PERMANENT"

"g.list vect" should list 6 vector files.  Each file contains a simple grid
with square cells each 5280x5280.

Start with 2x2example.  The procedure for all of the examples is identical.

"v.info 2x2example" to get the limits of the grid.

Use "g.region" to set the region to match the grid.  In the case of 2x2example
the limits are n=10560 s=0 w=0 e=10560 nsres=5280 ewres=5280.  In all examples
the resolution is 5280.  The region so defined should have two rows and two
columns. 

Ideally it seems like you should be able to set the region with "g.region
vect=2x2example."  That may work, but in most examples it does not.  If the
region is not set exactly to the dimensions of the grid then the next step
will produce unwanted results.

"v.to.rast i=2x2example o=2x2example"
"d.mon start=x0"
"d.rast -o 2x2example"
"d.vect 2x2example"

This should produce an accurate map, with the lines of the vector grid exactly
outlining the colored cells in the raster map.

"d.pan"  and click on the image to recenter (e.g., click your pointer in the
center of the northwest cell).  Repeat this step several times with different
center locations.  With some centers the map is accurate.  With other centers
the map is drawn with the raster file offset by half a cell width in both the
x and y directions.  When there is an offset you can use d.where, d.what.rast
or d.what.vect to confirm that the vector map is correctly located (southwest
corner at x=0, y=0) and the raster map is offset.

As an aside, d.what.vect in my 5.3 installation does not work correctly.  It
identifies the two cells on the top row correctly and reports consistent
measurements.  It misidentifies the two cells on the lower row and reports -0
for all dimensions.

The additional examples illustrate different effects.

"3x3example" uses a 3x3 grid.  Grids with odd-numbered dimensions don't have
the same problem with inaccuracy.

"10x10example and 11x11example" are similar to 2x2example and 3x3example and
are just there to show that the same effect occurs on slightly large grids.

"14x11example" has 14 columns and 11 rows.  In this case the grid and the
raster map always agree in the y direction (odd number of rows) but often
disagree in the x direction (even number of columns)

"shiftexample" has 14 columns and 11 rows, but the origin is shifted by 880 in
both x and y.  The shift does not alter the results.

If the region is not set correctly when the raster file is created (for
instance, if "g.region vect=..." is used to set the region) then the raster
map will have additional problems and will not correspond accurately to the
vector grid.  That seems to be a different problem.


Roger Miller