[GRASS-dev] error v.what.rast

[Paulo van Breugel]

On Fri, Oct 10, 2014 at 1:44 PM, Glynn Clements <glynn at gclements.plus.com>
wrote:

>
> Paulo van Breugel wrote:
>
> > > Perhaps a vector point just hits the extents of the computational
> region
> > > or map?
> >
> > You are right that the extends do not match completely. But I was under
> the
> > impression that points outside the region simply got skipped?
>
> I believe that *ought* to happen. The question is whether the fix
> should go into v.what.rast or libraster.
>
> > For example, if I reduce the region to a small area, and run
> > v.what.rast, I get the warning message "WARNING: 2287133 points
> > outside current region were skipped", after which the it happily
> > continues to run updating the points that are within the region
> > bounds, as per below.
>
> The behaviour of Rast_get_row() etc is that requesting a row which is
> outside of the map's bounds returns a row of nulls, but requesting a
> row which is outside of the current region is a fatal error.
>
> The former is just how raster maps are conceived. Conceptually, a
> raster map has infinite extent, but all non-null values are contained
> within a finite region. A raster map's stored bounds provide a loose
> bound on the region containing the non-null values, as well as
> avoiding requiring infinite disk space for each map.
>

Yes, the way grass handles points outside the region makes sense to me; it
follows the general philosophy of grass to ignore what is outside the
region.

>
> The latter indicates a programming error. Any module which accesses a
> raster map first specifies (implicitly or explicitly) a finite grid of
> sample locations. Subsequent Rast_get_row() etc calls return a
> one-dimensional array of values obtained by sampling the map at the
> locations for one row of the grid. The error occurs when the program
> requests data for a row which doesn't exist in the sampling grid which
> it specified.
>

What would make sense to me is if for those points to update the table with
a NULL value, as that is what it really is (considering that all values
within the defined region are valid).


> It would be fairly straightfoward to work around the issue within
> libraster, by making compute_window_row() just return 0 instead of
> raising a fatal error. But this could mask more significant issues,
> resulting in modules silently returning wrong results.
>
> Personally, I feel that v.what.rast should be fixed. Specifically, the
> code
>
>         row = Rast_northing_to_row(Points->y[0], &window);
>         col = Rast_easting_to_col(Points->x[0], &window);
>
> should be followed by:
>
>         if (col < 0 || col >= window.cols || row < 0 || row >= window.rows)
>             continue;
>
> The Vect_point_in_box() call just above that code doesn't appear to be
> sufficient in the case where the point lies either on the boundary of
> the region or just inside.
>
> This is partly because Vect_point_in_box() uses >= and <= rather than
> > and < (i.e. a point on the boundary is considered inside), and
> partly because Rast_northing_to_row() etc are affected by rounding
> error (window.north - window.south is not necessarily equal to
> window.ns_res * window.rows, similarly for the horizontal direction).
>
> --
> Glynn Clements <glynn at gclements.plus.com>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/grass-dev/attachments/20141010/4c02854c/attachment.html>