[Proj] vertical projection

Thomas Knudsen knudsen.thomas at gmail.com
Wed Oct 12 14:12:24 EST 2016


Well, the pipeline code is under reconstruction: Having had to build a
large amount of infrastructure while implementing the pipeline code, the
pull request ended up somewhat messy, and I closed it,

I am now in the process of refactoring in order to first build the
infrastructure, then the solution on top of it. The first part appeared as
pull request #431 yesterday (cf. https://github.com/OSGeo/proj.4/pull/431
 -  I would be thankful for any comments you may have to offer).

Next steps will probably be to reintroduce the +proj=pipeline "pseudo
projection", followed by the elementary operations (Helmert, Horner,
Cartesian/Ellipsoidal)., and finally the 3D-transformation program "tran",
and some documentation.

The intention is to make pipelines work transparently also for the original
proj.4 2D API calls, and the more recently added 3D calls
(pj_fwd3d/pj_inv3d): All the ugly stuff is taken care of inside the
pipeline-pseudo-projection driver, so given the transformation intended
makes sense in 2D it should work once I get around to polish it up -
although the main intention of the work is to facilitate full 3D geodetic
transformations in all their glorious convolution :-)

So Vincent, as Even points out: The pipeline stuff may be what you need -
but it is not ready for prime time yet.

I do, however, believe that the refactoring of the code will make it much
more comprehensible, and probably also more stable. But stability is
correlated with usage as well, so I will highly appreciate any test efforts
and/or guidance you may be able to provide as the work proceeds.

Best regards
Thomas



2016-10-05 12:18 GMT+02:00 Even Rouault <even.rouault at spatialys.com>:

> Le mercredi 05 octobre 2016 10:21:31, Vincent Mora a écrit :
> > Hi,
> >
> > I need to draw terrain sections for several projects (hydrology,
> geology).
> >
> > The vertical section is defined by a linestring in the (x,y) plane (e.g.
> > +proj=lcc). Data in the vicinity of the section are ortho-projected on
> > the folded plane with coordinates (s,z), where s is the point location
> > on the line (meters or feet).
> >
> > For the moment the CRS part of the QGis plugin we are developping is a
> > bit of a hack, and I'd like to introduce that kind of projection in
> > proj4 since I believe it belongs there.
> >
> > A projection on a folded plane is rather common, but in this case:
> > - the folds are sharp,
> > - the projection is not well defined for points that lie on the bisector
> > of two line-segments
> > - the projection is meaninigless for points that are not in the plane
> > vicinity.
> >
> > The proj4 string should contain:
> > - the line definition (e.g. +wkt=LINESTRING(...))
> > - the classical definition of the CRS of the line,
> > - and a thickness (vicinity of the plane)
> >
> > Questions:
> > - has such a CRS it's place in proj4 ?
> > - do you see obvious errors/misconceptions ?
>
> Vincent,
>
> The fact that the line definition can be defined in any CRS handled by
> proj.4
> would make it quite a odd object for proj.4 if we used current mechanisms.
> You'd need to have a definition like :
> +proj=vertplane +wkt=LINESTRING(...) +wkt_proj=lcc +wkt_proj_lat_1= ...
>
> Probably a better fit would be to rely on Thomas Knudsen's ongoing work on
> transformation pipelines (
> http://osgeo-org.1560.x6.nabble.com/Transformation-pipelines-your-opinion-
> td5269960.html ,
> http://thomasknudsen.net/tran_users_guide.html ).
> I'm not sure this has yet landed into master. Perhaps Thomas can confirm ?
>
> With that framework, I think you could define a 2 stage pipeline :
> - first converting from geographical coordinates to LCC projected ones
> - then doing the projection onto the vertical plane
>
> +proj=pipeline +step=lcc +lat_1=...  +step +proj=vertplane
> +wkt=LINESTRING(...) +thickness=
>
> Regarding QGIS integration, this would require some changes, due to
> pipeline support
> requiring a new dedicated API in proj.4 "TRIPLET pj_apply_projection
> (TRIPLET point, int direction, PJ *P)" since you cannot use the
> existing int pj_transform( projPJ srcdefn, projPJ dstdefn, long
> point_count, int point_offset, double *x, double *y, double *z )
>
> By the way, Thomas, wouldn't that make sense to have a pipeline
> transformation
> API that would accept several points at once ? The default implementation
> could
> just loop over them, but that would be a provision for potential later
> optimizations
> where you could potentially vectorize computations.
>
> Even
>
> --
> Spatialys - Geospatial professional services
> http://www.spatialys.com
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://lists.maptools.org/pipermail/proj/attachments/20161012/65f65324/attachment.htm 


More information about the Proj mailing list