首页
学习
活动
专区
圈层
工具
发布
首页
学习
活动
专区
圈层
工具
MCP广场
社区首页 >问答首页 >使用NetTopologySuite计算点之间的地理距离

使用NetTopologySuite计算点之间的地理距离
EN

Stack Overflow用户
提问于 2021-04-13 01:03:06
回答 2查看 248关注 0票数 2

我正在尝试使用NetTopologySuite计算两点之间的距离。因为我参考了微软的文档,所以我想出了以下GeometryExtension和GeometryHelper类:

代码语言:javascript
复制
        public static class GeometryExtensions
    {
        private static readonly CoordinateSystemServices _coordinateSystemServices = new CoordinateSystemServices(new Dictionary<int, string>
            {
                // Coordinate systems:

                [4326] = GeographicCoordinateSystem.WGS84.WKT,

                // CRS for Denmark ESPG 25832. Source: https://epsg.io/25832 and https://sdfe.dk/
                [25832] = @"PROJCS[""ETRS89 / UTM zone 32N"",
                                GEOGCS[""ETRS89"",
                                    DATUM[""European_Terrestrial_Reference_System_1989"",
                                        SPHEROID[""GRS 1980"",6378137,298.257222101,
                                            AUTHORITY[""EPSG"",""7019""]],
                                        TOWGS84[0,0,0,0,0,0,0],
                                        AUTHORITY[""EPSG"",""6258""]],
                                    PRIMEM[""Greenwich"",0,
                                        AUTHORITY[""EPSG"",""8901""]],
                                    UNIT[""degree"",0.0174532925199433,
                                        AUTHORITY[""EPSG"",""9122""]],
                                    AUTHORITY[""EPSG"",""4258""]],
                                PROJECTION[""Transverse_Mercator""],
                                PARAMETER[""latitude_of_origin"",0],
                                PARAMETER[""central_meridian"",9],
                                PARAMETER[""scale_factor"",0.9996],
                                PARAMETER[""false_easting"",500000],
                                PARAMETER[""false_northing"",0],
                                UNIT[""metre"",1,
                                    AUTHORITY[""EPSG"",""9001""]],
                                AXIS[""Easting"",EAST],
                                AXIS[""Northing"",NORTH],
                                AUTHORITY[""EPSG"",""25832""]]"
            }
        );
        


        /// <summary>
        /// Projects a geometry to another SRID
        /// </summary>
        /// <param name="geometry"></param>
        /// <param name="targetSrid"></param>
        /// <param name="defaultSourceSrid">If the geometry SRID has not been specified (i.e. equals 0) defaultSourceSrid is used</param>
        /// <returns></returns>
        public static Geometry ProjectTo(this Geometry geometry, int targetSrid, int? defaultSourceSrid = null)
        {
            if (geometry == null)
                throw new Exception("Geometry is null, cannot project");

            var sourceSrid = geometry.SRID == 0 && defaultSourceSrid.HasValue ? defaultSourceSrid.Value : geometry.SRID;
            var transformation = _coordinateSystemServices.CreateTransformation(sourceSrid, targetSrid);

            var result = geometry.Copy();
            result.Apply(new MathTransformFilter(transformation.MathTransform));

            return result;
        }
}




     public  class GeometryHelper
    {
        private static readonly int _longLatSRID = 4326;
        private static readonly int _targetSRID = 25832;

        /// <summary>
        /// In order to get the distance in meters, we need to project to an appropriate
        /// coordinate system. If no SRID is provided 25832, which covers Denmark is used.
        /// If the provided Points have no SRID, 4326 (longitude/latitude) is assumed.
        /// </summary>
        /// <param name="a"></param>
        /// <param name="b"></param>
        /// <param name="targetSrid"></param>
        /// <returns></returns>
        public static double DistanceInMeters(Point a, Point b, int? targetSrid = null)
        {
            targetSrid ??= _targetSRID;

            try
            {
                //If SRID is not set, assume long/lat, ie. 4326
                return a.ProjectTo(targetSrid.Value, _longLatSRID).Distance(b.ProjectTo(targetSrid.Value, _longLatSRID));
            }
            catch (Exception e)
            {
                throw new Exception("Failed to calculate distance", e);
            }
        }
}

为了测试我的距离计算是否正确,我打开Google Maps并选择了计算3对不同点之间的距离。然后,我复制粘贴了测试中的值,以查看它们是否匹配。我从我的代码中获得的值比我从Google Maps中获得的值大2倍。我做错了什么?以下是我的测试:

代码语言:javascript
复制
public class GeometryHelperTests
{
    [Theory]
    [InlineData(55.676518126293466, 12.567203858066554, 55.67645068023813, 12.56698376863065, 15.68)]//I am getting 36.1
    [InlineData(55.659368700924475, 12.546625254609248, 55.65940085355421, 12.546601114728679, 3.77)]//I am getting 6.2
    [InlineData(55.65896978705746, 12.546674114514795, 55.6596855501795, 12.547258269821455, 87.09)]//I am getting 173.7
    public void DistanceInMeters(double x1, double y1, double x2, double y2, double expected)
    {
        // setup
        Point point1 = new Point(x1, y1) { SRID = 4326 };
        Point point2 = new Point(x2, y2) { SRID = 4326 };

        //CoordinateSystemFactory csFact = new CoordinateSystemFactory();
        

        // act
        var result = GeometryHelper.DistanceInMeters(point1, point2);

        

        // assert
        result.Should().Be(expected);
    }

}
EN

回答 2

Stack Overflow用户

回答已采纳

发布于 2021-04-13 18:47:04

您需要计算great circle distance。NetTopologySuite Point.Distance方法返回cartesian distance

尝试以下操作:

代码语言:javascript
复制
public static double Radians(double x)
{
    return x * Math.PI / 180;
}

public static double GreatCircleDistance(double lon1, double lat1, double lon2, double lat2)
{
    double R = 6371e3; // m

    double sLat1 = Math.Sin(Radians(lat1));
    double sLat2 = Math.Sin(Radians(lat2));
    double cLat1 = Math.Cos(Radians(lat1));
    double cLat2 = Math.Cos(Radians(lat2));
    double cLon = Math.Cos(Radians(lon1) - Radians(lon2));

    double cosD = sLat1*sLat2 + cLat1*cLat2*cLon;

    double d = Math.Acos(cosD);

    double dist = R * d;

    return dist;
}

您还可以使用内置的GeoCoordinate.GetDistanceTo(),它实现了小距离更精确的Haversine公式。

票数 2
EN

Stack Overflow用户

发布于 2021-04-13 20:55:56

尽管saeedkazemi的答案是正确的,给出的结果更接近真实值,但我的代码给出的结果也相当接近,但我发现我必须反转谷歌地图给出的坐标。所以如果给我55.6765,12.5672,我需要给这个公式12.5672,55.6765。话虽如此,saeedkazemi的答案提供了更接近真实价值的结果。

票数 0
EN
页面原文内容由Stack Overflow提供。腾讯云小微IT领域专用引擎提供翻译支持
原文链接:

https://stackoverflow.com/questions/67062622

复制
相关文章

相似问题

领券
问题归档专栏文章快讯文章归档关键词归档开发者手册归档开发者手册 Section 归档