[Proj] Data-driven projections Re: Understanding 2D Helmert

[Luke Bergmann]

Dear all,

As a fan of the proj[.4] project, as a lurker on these lists, and as hopefully a future contributor to proj, I thought I should jump in here with my first post. Please forgive me my mistakes in protocol and my ignorance! This thread's discussion has come close to functionality that a colleague of mine (David O’Sullivan) and I have long thought would really transform the capabilities of proj and projects downstream (e.g., QGIS), making them much more useful to geography, to ‘analytical cartography’, and to the digital humanities, among other pursuits. 

In particular, we would like to suggest (and we are willing to contribute work to see that) proj implement transformations that are data-driven, or as Waldo Tobler called them, ‘empirical projections’. The easiest example would be the transformation associated with a cartogram. Of course, the math/heuristic to calculate a cartogram is clearly beyond the scope of a library like proj. However, each existing cartogram embodies a transformation, and having the ability to represent that transformation in an approximate/empirical/data-driven form within a proj string would enable your average GIS to suddenly work with a vastly expanded world of coordinate systems. Similarly, there is an explosion of creative transformations out there these days (think about all the heuristics unfolding the globe these days) that may themselves be data-driven and/or may not have been implementable in proj for various reasons—this would enable everyday GIS users to use them. 

We recognize there is already a ‘grid’ functionality for datum transformations within proj, but having a proj projection that interprets and implements arbitrary projections encoded within proj strings is both closer to the spirit of understanding these transformations as projections (as Tobler and others have advocated) and would seem to be  important from a practical perspective, as there are already ways in GIS such as QGIS to input proj strings and then call proj. As academics, we want this functionality to be usable by students in the introductory GIS course. There are other benefits to knowledge as well. The humanities, which see space as more complex, often imagined, would benefit from these transformations. And they would be key to proj enabling the ‘analytical' use of cartography that Waldo Tobler advocated with his ‘analytical cartography’.

How might this work? We gave a paper at SIGSPATIAL in 2017 that started to sketch out one idea. I’ve linked the paper below—I can email you a copy but the proj list doesn’t accept big attachments and I haven’t been able to set up the open sharing function with ACM yet, sorry! Note that the paper goes beyond what I’m able to describe here. Just as this list has been discussing transformations with hundreds of triangles, so too were we. To start, we thought that a proj string could encode how a set of points in one space map into coordinates within another space. Here would be a generic notional proj.4 CRS string:

+proj=ptobler +interpolation=interptype +mapping=mapstring 

where:

- ptobler is a string literal that represents a possible shorthand name for the class of generic projections (‘Ptobler’ is a somewhat droll honorary offered to Waldo Tobler’s cartographic ideas which also references the cartographic achievements of Ptolemy. It first emerged out of conversations around the early Michigan Inter-University Community of Mathematical Geographers (MICMOG))

- interptype is a string representing the name of an interpolation method, e.g., polynomial, bicubic, thin plate spline, TIN. The interpolation method is used to estimate how arbitrary points in (φ,λ) coordinates map into (x,y) space.

- mapstring is a string serialization of the list of known, specified mappings between points in an absolute space such as latitude and longitude (φ, λ) and the projected space (x,y). A naive serialization replacing mapstring in an actual proj string might proceed as follows: ′(φ1,λ1),(x1,y1);(φ2,λ2),(x2,y2);...(φN,λN),(xN,yN);′ All φi, λi, xi, and yi are likewise to be replaced by the appropriate real numbers.

To determine how an arbitrary point would be transformed, as suggested above, different interpolation methods could be defined. One such interpolation would involve having first computed the Delaunay triangulation of the known source points. The arbitrary point would be located within one of those triangles and the associated known transformation of that triangle into the output coordinate system would be used to transform the arbitrary point.

We know that proj strings can’t be infinite in practice, for a variety of reasons, but those limitations still enable a tremendous amount of functionality. We were thinking that one could fit the needed code within a single classical .c file used to implement a projection in proj if needed, though if there were other places for useful code, such as for triangulation or spatial indexes or queries, that would be at least as helpful. 

All this being said, we are new to proj as contributors, and indeed, are academics not professional software developers, so we wanted to come to the community first before going further. Your thoughts are most appreciated. Further, we want to acknowledge and thank you for all your varied efforts to create, maintain, and extend this key node in the open source geo world!

Best wishes,
Luke


Happy to send you a copy of:

Luke R. Bergmann and David O’Sullivan. 2017. Computing with many spaces: Generalizing projections for the digital geohumanities and GIScience. In Proceedings of GeoHumanities’17: 1st ACM SIGSPATIAL Workshop on Geospatial Humanities, Los Angeles Area, CA, USA, November 7–10, 2017 (GeoHumanities’17), 8 pages.
https://doi.org/10.1145/3149858.3149866

https://dl.acm.org/citation.cfm?id=3149866




> On Oct 28, 2018, at 9:46 AM, Even Rouault <even.rouault at spatialys.com> wrote:
> 
> On dimanche 28 octobre 2018 15:02:05 CET Kristian Evers wrote:
>> Why is that? To me it seems somewhat similar to polynomial mappings
>> (Horner), in the way that many coefficients are likely to be stored for a
>> given transformation. Given a clever way of defining the triangles and
>> their corresponding Helmet parameters I think it should be possible to make
>> these transformations work. At least as PROJ strings.
> 
> You would need to store hundreds of triangles and their parameters, and have 
> PROJ figure out in which triangle a point it. Certainly doable, but not 
> directly in the scope of my current work.  And avoid the database to grow too 
> big with very particular transforms. That said, 600 triangles with 6 double 
> parameters each + the x,y coordinates is just 38.4 KB if you store that as a 
> compact binary blob. Some balance to know if it is something acceptable in-db 
> (EPSG style transformations have at most 18 parameters each) or out-db 
> (grids). The suggestion for the grid approach was to make it more immediatly 
> usable.
> 
> -- 
> Spatialys - Geospatial professional services
> http://www.spatialys.com
> _______________________________________________
> Proj mailing list
> Proj at lists.maptools.org
> http://lists.maptools.org/mailman/listinfo/proj