[GRASS-dev] latlon support

[Glynn Clements]

Markus Metz wrote:

> > IMHO, the optimal solution for lat/lon is to remove the region
> > constraints, and require the data to be "Euclidified" either on input
> > or by the low-level retrieval functions (e.g. G_get_raster_row etc).
> 
> As long as the input is ok, I don't see a problem with GRASS modules, in 
> theory the output should conform as well to Euclidean geometry as the 
> input. But here the problem starts. I get all sorts of dirty data, in 
> particular shapefiles, that I am supposed to work with, and cleaning 
> them up is sometimes really a pain... After cleaning I have 
> topologically correct vectors, but that doesn't guarantee that they are 
> "Euclidified". If I understand you right I would have to check for every 
> single line and every single boundary if every single vertex is correct. 
> In practice I can only do an educated guess if, to stick to my previous 
> example, it is correct that one point has lon 179 and the next point has 
> lon 181, i.e. crossing the datum border is right, or if it is correct to 
> go instead from lon 179 to lon -179, i.e. go over Greenwich meridian and 
> not over the datum border.

By "Euclidified", I mean that the data is in a form such that code
which process it can treat -179 and +181 as being distinct points 360
units apart.

E.g. for a global map 360 degrees wide and centred on Greenwich,
Euclidification would involve either copying or moving the Chukchi
Peninsula over to the left-hand side of the map, so that it doesn't
get chopped off by the right-hand edge.

Similarly, to use r.resamp.interp on a global DEM without getting a
gap at the 180th meridian, the longitude bounds would need to be set
to 360 degrees plus 4 cells, so that the cells around the 180th
meridan appear on both sides of the map, as the code will assume that
the columns either side of column n are simply n+1 and n-1.

[r.resamp.interp can perform the enlargement itself, but libgis has to
permit it, and not insist that -182..+182 means +178..+182.]

> > E.g. you could freely set the region to west=-360, east=+360, and get
> > two copies of the globe. The points 0,-360 and 0,+360 would be treated
> > as being separated by twice the equatorial circumference. And so on.
> 
> That can be intended but I am afraid not all people generating spatial 
> data think about this effect. It may be safer to convert coordinates to 
> the 180 deg lon and 90 deg lat limits and take the shortest way between 
> two points.

My point is that this needs to be done in the lowest-level libraries,
not in higher-level libraries or the modules.

High-level code must be free to assume that eastings increase
monotonically from left to right, that if a point has an easting of x,
then the point which is d units to its east has an easting of x+d, and
so on.

Requiring code to explicitly use specialised (spherical-aware)
functions in place of +,-,<,>,== etc is not a realistic strategy.

-- 
Glynn Clements <glynn at gclements.plus.com>