如何在PROJ库中调用proj_lp_dist

Question

我正在尝试计算libray 7.2.2项目中的测地距离。但是，无论我输入什么数据，proj_lp_dist函数都会返回无穷大。命令行工具可以正常工作，但它不会调用相同的API。我尝试反转long-lat值，但这不能解决问题。

c:\PROJ\build_vs2019\bin\Release>geod +ellps=WGS84 -I
28.472672560770082 -81.466635213256282 37.41439966801741 -122.07673673962702
-64d33'38.315"  92d36'51.975"   3892821.141

        TEST_METHOD(TestDistanceExperimentalReverse)
        {
            wchar_t message[2000];
            SetEnv();

            PJ_CONTEXT* C = proj_context_create();
            Assert::IsNotNull(C, L"Cannot create PROJ context.");

            PJ* P = proj_create(C,"EPSG:4326"/*WGS84*/);
            Assert::IsNotNull(P, L"Cannot create PROJ converter");

            PJ_COORD a, b;
            //Universal Studios, Orlando, Florida
            a = proj_coord(-81.466635213256282, 28.472672560770082, 0, 0);
            //Computer history Museum, Mountain View, California
            b = proj_coord(-122.07673673962702, 37.41439966801741, 0, 0);

            double dist = proj_lp_dist(P, a, b);

            swprintf_s(message, L"distance %g\n", dist);
            Logger::WriteMessage(message);

            {
                wchar_t msg[2000];

                int errnum = proj_context_errno(C);
                if (errnum != 0)
                {
                    const char* c = proj_errno_string(errnum);
                    size_t cSize = strlen(c) + 1;
                    mbstowcs_s(&cSize, msg, c, cSize);
                    swprintf_s(message, L"Context state %d %s", errnum, msg);
                    Logger::WriteMessage(message);
                }
            }
            {
                wchar_t msg[2000];
                int errnum = proj_errno(P);
                if (errnum != 0)
                {
                    const char* c = proj_errno_string(errnum);
                    size_t cSize = strlen(c) + 1;
                    mbstowcs_s(&cSize, msg, c, cSize);
                    swprintf_s(message, L"Projection state %d %s", errnum, msg);
                    Logger::WriteMessage(message);
                }
            }
            Assert::IsFalse(isnan(dist), L"Distance is NaN");
            Assert::IsFalse(isinf(dist), L"Distance is Infinity");
            Assert::AreEqual(3892821.141, dist);
        };

Answer 1

我也有同样的问题。我在项目邮件列表上发了帖子，发现它是个bug。对于某些创建方法，PJ对象未正确初始化：

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

完整讨论：http://osgeo-org.1560.x6.nabble.com/PROJ-Projection-type-required-for-proj-lp-dist-td5456628.html

在8.0.1中修复了这个问题：https://github.com/OSGeo/PROJ/releases/tag/8.0.1

使proj_lp_dist()和proj_geod()在PJ* CRS对象(2570) https://github.com/OSGeo/PROJ/pull/2570上工作

我无法验证它是否已经修复，因为我在Windows上使用OSGeo4W安装PROJ，它仍然使用Proj6。我还没有尝试自己构建PROJ。