On Fri, Oct 10, 2014 at 1:44 PM, Glynn Clements <gl...@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 <gl...@gclements.plus.com> >
_______________________________________________ grass-dev mailing list grass-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/grass-dev