[PROJ] Projection type required for proj_lp_dist

Even Rouault even.rouault at spatialys.com
Wed Mar 10 10:07:10 PST 2021


Fix pending in https://github.com/OSGeo/PROJ/pull/2570

Documentation of those functions would deserve enhancements though. See 
my comment in the PR

Le 10/03/2021 à 16:39, Kristian Evers a écrit :
> Sorry, I thought your problem was getting set up before the call to proj_lp_dist.
> I’ve tried to reproduce your problem and I believe you’ve come across a bug.
> Below is a bit of test code I wrote to demonstrate the behaviour. As far as I can tell
> you get similar results.
>
>
> #include <stdio.h>
> #include "proj.h"
>
> int main() {
>
>      PJ_CONTEXT *C = proj_context_create();
>      PJ *P;
>      //P = proj_create(C, "EPSG:4258");                                  // returns inf
>      //P = proj_get_ellipsoid(C, proj_create(C, "EPSG:4258"));           // returns inf
>      //P = proj_create(C, "+proj=latlong +ellps=GRS80");                 // returns 12836.740931
>      //P = proj_create(C, "+proj=latlong +ellps=GRS80 +type=crs");       // returns inf
>      //P = proj_create(C, "+proj=utm +zone=32 +ellps=GRS80");            // returns 12836.740931
>      //P = proj_create_crs_to_crs (C,"EPSG:4326", "EPSG:4258", NULL);    // returns 12836.740931
>      P = proj_create_crs_to_crs (C,"EPSG:4326", "EPSG:32632", NULL);     // returns 12836.740931
>
>      PJ_COORD c1 = proj_coord(proj_torad(12.0), proj_torad(55.0), 0, 0);
>      PJ_COORD c2 = proj_coord(proj_torad(12.1), proj_torad(55.1), 0, 0);
>
>      double d = proj_lp_dist(P, c1, c2);
>
>      printf("Distance: %f", d);
>
>      proj_destroy(P);
>
>      return 0;
>
> }
>
>
> I would expect all of the above (perhaps with the exception of the proj_get_ellipsoid
> line) to return the same distance.
>
> This is most likely because internally P->geod is not initialized when creating a CRS object.
> I haven’t been able to pinpoint exactly where but someone more familiar with that part of the
> code might be able to do so quickly.
>
>
>
>> On 10 Mar 2021, at 15:27, Jason Dick <JasonD at TheTranstecGroup.com> wrote:
>>
>> On 3/9/2021 5:06 PM, Kristian Evers wrote:
>>> Try taking a look at this function from the testsuite:
>>> https://github.com/OSGeo/PROJ/blob/f278b5bee641dd94245ed54f4eb75e29c3d4d993/test/unit/test_c_api.cpp#L166-L203
>> Thanks for the pointer, but it's not getting through to me. From what I gather, Infinity is the correct answer when using EPSG:4326, but it doesn't explain why. Why does distance require an operation PJ and not a CRS? I've looked at the source code for the distance function before, but it wasn't helpful either. Keep in mind, I have no geo background. I just use libraries for transforming coordinates. I've learned a lot, but PROJ is much more "technical" than what I've used in the past.
>>
> _______________________________________________
> PROJ mailing list
> PROJ at lists.osgeo.org
> https://lists.osgeo.org/mailman/listinfo/proj

-- 
http://www.spatialys.com
My software is free, but my time generally not.



More information about the PROJ mailing list