GMap.net 地图标绘实现(一)-程序员宅基地

  • GeometryUtil帮助类:主要用于空间参考(坐标)转换,距离计算等等。
public class GeometryUtil
    {
        private const double EARTH_RADIUS = 6378.137;//地球半径

        /// <summary>
        /// 获取多边形几何中心点
        /// </summary>
        /// <param name="mapPolygon"></param>
        /// <returns></returns>
        public static PointLatLng GetCenterPoint(GMapPolygon mapPolygon)
        {
            if (mapPolygon == null)
            {
                return new PointLatLng();
            }
            List<GeoAPI.Geometries.Coordinate> coordinates = new List<GeoAPI.Geometries.Coordinate>();
            foreach (var value in mapPolygon.Points)
            {
                GeoAPI.Geometries.Coordinate coordinate = new GeoAPI.Geometries.Coordinate(value.Lng, value.Lat);
                coordinates.Add(coordinate);
            }
            if (mapPolygon.Points.Count == 2)
            {
                NetTopologySuite.Geometries.LineString lineString = new NetTopologySuite.Geometries.LineString(coordinates.ToArray());
                IPoint point = lineString.Centroid;
                return new PointLatLng(point.Y, point.X);
            }
            else
            {
                if (mapPolygon.Points[0].Lng != mapPolygon.Points[mapPolygon.Points.Count - 1].Lng ||
                mapPolygon.Points[0].Lat != mapPolygon.Points[mapPolygon.Points.Count - 1].Lat)
                {
                    coordinates.Add(new GeoAPI.Geometries.Coordinate(mapPolygon.Points[0].Lng, mapPolygon.Points[0].Lat));
                }

                NetTopologySuite.Geometries.LinearRing linearRing = new NetTopologySuite.Geometries.LinearRing(coordinates.ToArray());
                NetTopologySuite.Geometries.Polygon polygon = new NetTopologySuite.Geometries.Polygon(linearRing);
                IPoint point = polygon.Centroid;
                return new PointLatLng(point.Y, point.X);
            }
        }

        public static List<Coordinate> WGS84ToWebMercator(List<PointLatLng> pointLatLngs)
        {
            if (pointLatLngs == null)
                return null;

            IList<Coordinate> coors = PointsToCoords(pointLatLngs);

            NetTopologySuite.Geometries.PrecisionModel precisionModel = new NetTopologySuite.Geometries.PrecisionModel(GeoAPI.Geometries.PrecisionModels.Floating);
            GeoAPI.CoordinateSystems.ICoordinateSystem wgs84 = ProjNet.CoordinateSystems.GeographicCoordinateSystem.WGS84;
            GeoAPI.CoordinateSystems.ICoordinateSystem webMercator = ProjNet.CoordinateSystems.ProjectedCoordinateSystem.WebMercator;


            int SRID_wgs84 = Convert.ToInt32(wgs84.AuthorityCode);    //WGS84 SRID
            int SRID_webMercator = Convert.ToInt32(webMercator.AuthorityCode);    //WebMercator SRID

            ProjNet.CoordinateSystems.Transformations.CoordinateTransformationFactory ctFact = new ProjNet.CoordinateSystems.Transformations.CoordinateTransformationFactory();
            GeoAPI.CoordinateSystems.Transformations.ICoordinateTransformation transformation = ctFact.CreateFromCoordinateSystems(wgs84, webMercator);

            NetTopologySuite.Geometries.GeometryFactory factory_wgs84 = new NetTopologySuite.Geometries.GeometryFactory(precisionModel, SRID_wgs84);

            IList<Coordinate> coords = transformation.MathTransform.TransformList(coors);

            //var gcgs2000 = CreateCoordinateSystemFromWkt("GEOGCS[\"China Geodetic Coordinate System 2000\",DATUM[\"China_2000\",SPHEROID[\"CGCS2000\",6378137,298.257222101,AUTHORITY[\"EPSG\",\"1024\"]],AUTHORITY[\"EPSG\",\"1043\"]],PRIMEM[\"Greenwich\",0,AUTHORITY[\"EPSG\",\"8901\"]],UNIT[\"degree\",0.0174532925199433,AUTHORITY[\"EPSG\",\"9122\"]],AUTHORITY[\"EPSG\",\"4490\"]]");
            //var gcgs4508 = CreateCoordinateSystemFromWkt("PROJCS[\"CGCS2000 / Gauss-Kruger CM 111E\",GEOGCS[\"China Geodetic Coordinate System 2000\",DATUM[\"China_2000\",SPHEROID[\"CGCS2000\",6378137,298.257222101,AUTHORITY[\"EPSG\",\"1024\"]],AUTHORITY[\"EPSG\",\"1043\"]],PRIMEM[\"Greenwich\",0,AUTHORITY[\"EPSG\",\"8901\"]],UNIT[\"degree\",0.0174532925199433,AUTHORITY[\"EPSG\",\"9122\"]],AUTHORITY[\"EPSG\",\"4490\"]],PROJECTION[\"Transverse_Mercator\"],PARAMETER[\"latitude_of_origin\",0],PARAMETER[\"central_meridian\",111],PARAMETER[\"scale_factor\",1],PARAMETER[\"false_easting\",500000],PARAMETER[\"false_northing\",0],UNIT[\"metre\",1,AUTHORITY[\"EPSG\",\"9001\"]],AUTHORITY[\"EPSG\",\"4508\"]]");
            //IList<Coordinate> temps = ConvertCoordinates(coors.ToList(), gcgs2000,gcgs4508);
            //IList<Coordinate> temps1= ConvertCoordinates(temps.ToList(), gcgs4508,gcgs2000);


            return coords.ToList();
        }

        public static List<PointLatLng> WGS84ToWebMercator2(List<PointLatLng> pointLatLngs)
        {
            List<Coordinate> temps = WGS84ToWebMercator(pointLatLngs);
            List<PointLatLng> points = CoordsToPoints(temps);
            return points;
        }

        public static List<PointLatLng> WebMercatorToWGS84(List<PointLatLng> pointLatLngs)
        {
            return ConvertCoordinates(pointLatLngs, ProjNet.CoordinateSystems.ProjectedCoordinateSystem.WebMercator, ProjNet.CoordinateSystems.GeographicCoordinateSystem.WGS84);
        }

        public static List<PointLatLng> CoordsToPoints(List<Coordinate> coords)
        {
            List<PointLatLng> points = new List<PointLatLng>();
            foreach (var entity in coords)
            {
                points.Add(new PointLatLng(entity.Y, entity.X));
            }
            return points;
        }

        public static List<Coordinate> PointsToCoords(List<PointLatLng> points)
        {
            List<Coordinate> coords = new List<Coordinate>();
            foreach(var entity in points)
            {
                coords.Add(new Coordinate(entity.Lng, entity.Lat));
            }
            return coords;
        }

        /// <summary>
        /// 默认创建EPSG:4490 国家2000大地坐标系
        /// </summary>
        /// <param name="wkt"></param>
        /// <returns></returns>
        public static GeoAPI.CoordinateSystems.ICoordinateSystem CreateCoordinateSystemFromWkt(string wkt)
        {
            if (string.IsNullOrEmpty(wkt))
                wkt = "GEOGCS[\"China Geodetic Coordinate System 2000\",DATUM[\"China_2000\",SPHEROID[\"CGCS2000\",6378137,298.257222101,AUTHORITY[\"EPSG\",\"1024\"]],AUTHORITY[\"EPSG\",\"1043\"]],PRIMEM[\"Greenwich\",0,AUTHORITY[\"EPSG\",\"8901\"]],UNIT[\"degree\",0.0174532925199433,AUTHORITY[\"EPSG\",\"9122\"]],AUTHORITY[\"EPSG\",\"4490\"]]";
            try
            {
                ICoordinateSystemFactory coordinateSystemFactory = new ProjNet.CoordinateSystems.CoordinateSystemFactory();
                GeoAPI.CoordinateSystems.ICoordinateSystem coordinateSystem = coordinateSystemFactory.CreateFromWkt(wkt);
                return coordinateSystem;
            }
            catch
            {
                return null;
            }
        }

        public static List<Coordinate> ConvertCoordinates(List<Coordinate> coordinates, ICoordinateSystem sourceCoordinateSystem, ICoordinateSystem targetCoordinateSystem)
        {
            if (coordinates == null || sourceCoordinateSystem == null || targetCoordinateSystem==null)
                throw new Exception("请输入正确的参数");


            NetTopologySuite.Geometries.PrecisionModel precisionModel = new NetTopologySuite.Geometries.PrecisionModel(GeoAPI.Geometries.PrecisionModels.Floating);


            int SRID_source = Convert.ToInt32(sourceCoordinateSystem.AuthorityCode);    //source SRID
            int SRID_target = Convert.ToInt32(targetCoordinateSystem.AuthorityCode);    //target SRID

            ProjNet.CoordinateSystems.Transformations.CoordinateTransformationFactory transformationFactory = new ProjNet.CoordinateSystems.Transformations.CoordinateTransformationFactory();
            GeoAPI.CoordinateSystems.Transformations.ICoordinateTransformation transformation = transformationFactory.CreateFromCoordinateSystems(sourceCoordinateSystem, targetCoordinateSystem);
            
            //NetTopologySuite.Geometries.GeometryFactory factory_source = new NetTopologySuite.Geometries.GeometryFactory(precisionModel, SRID_source);

            IList<Coordinate> coords = transformation.MathTransform.TransformList(coordinates);
            return coords.ToList();

        }

        public static List<PointLatLng> ConvertCoordinates(List<PointLatLng> coordinates, ICoordinateSystem sourceCoordinateSystem, ICoordinateSystem targetCoordinateSystem)
        {
            List<Coordinate> tempCoordinates = new List<Coordinate>();
            foreach(var entity in coordinates)
            {
                tempCoordinates.Add(new Coordinate(x: entity.Lng, y: entity.Lat));
            }


            if (coordinates == null || sourceCoordinateSystem == null || targetCoordinateSystem == null)
                throw new Exception("请输入正确的参数");


            NetTopologySuite.Geometries.PrecisionModel precisionModel = new NetTopologySuite.Geometries.PrecisionModel(GeoAPI.Geometries.PrecisionModels.Floating);


            int SRID_source = Convert.ToInt32(sourceCoordinateSystem.AuthorityCode);    //source SRID
            int SRID_target = Convert.ToInt32(targetCoordinateSystem.AuthorityCode);    //target SRID

            ProjNet.CoordinateSystems.Transformations.CoordinateTransformationFactory transformationFactory = new ProjNet.CoordinateSystems.Transformations.CoordinateTransformationFactory();
            GeoAPI.CoordinateSystems.Transformations.ICoordinateTransformation transformation = transformationFactory.CreateFromCoordinateSystems(sourceCoordinateSystem, targetCoordinateSystem);

            //NetTopologySuite.Geometries.GeometryFactory factory_source = new NetTopologySuite.Geometries.GeometryFactory(precisionModel, SRID_source);

            IList<Coordinate> coords = transformation.MathTransform.TransformList(tempCoordinates);

            List<PointLatLng> points = new List<PointLatLng>();
            foreach(var entity in coords)
            {
                if (double.IsNaN(entity.X) || double.IsNaN(entity.Y))
                    continue;
                points.Add(new PointLatLng(entity.Y, entity.X));
            }
            return points;

        }

        public static double GetArea(List<PointLatLng> pointLatLngs)
        {
            if (pointLatLngs == null)
            {
                return 0;
            }
            if (pointLatLngs.Count <= 2)
            {
                return 0;
            }
            
            List<GeoAPI.Geometries.Coordinate> coordinates = WGS84ToWebMercator(pointLatLngs);

            if (coordinates[0].X != coordinates[coordinates.Count - 1].X ||
                coordinates[0].Y != coordinates[coordinates.Count - 1].Y)
            {
                coordinates.Add(new GeoAPI.Geometries.Coordinate(coordinates[0].X, coordinates[0].Y));
            }

            NetTopologySuite.Geometries.LinearRing linearRing = new NetTopologySuite.Geometries.LinearRing(coordinates.ToArray());
            NetTopologySuite.Geometries.Polygon polygon = new NetTopologySuite.Geometries.Polygon(linearRing);
            IPoint point = polygon.Centroid;
            return polygon.Area;
        }

        public static double GetArea(GMapPolygon mapPolygon)
        {
            if (mapPolygon == null)
            {
                return 0;
            }
            return GetArea(mapPolygon.Points);
        }

        public static RectLatLng GetRegionMaxRect(GMapPolygon polygon)
        {
            double latMin = 90;
            double latMax = -90;
            double lngMin = 180;
            double lngMax = -180;
            foreach (var point in polygon.Points)
            {
                if (point.Lat < latMin)
                {
                    latMin = point.Lat;
                }
                if (point.Lat > latMax)
                {
                    latMax = point.Lat;
                }
                if (point.Lng < lngMin)
                {
                    lngMin = point.Lng;
                }
                if (point.Lng > lngMax)
                {
                    lngMax = point.Lng;
                }
            }

            return new RectLatLng(latMax, lngMin, lngMax - lngMin, latMax - latMin);
        }

        public static double GetDistance(double lat1, double lng1, double lat2, double lng2)
        {
            PointLatLng first = new PointLatLng(lat2,lng1);
            PointLatLng second = new PointLatLng(lat2,lng2);
            List<PointLatLng> pointLatLngs = new List<PointLatLng>();
            pointLatLngs.Add(first);
            pointLatLngs.Add(second);

            var coordinates= WGS84ToWebMercator(pointLatLngs);
            NetTopologySuite.Geometries.LineString lineString = new NetTopologySuite.Geometries.LineString(coordinates.ToArray());
            return lineString.Length;
        }

        public static double GetDistance2(double lat1, double lng1, double lat2, double lng2)
        {
            double radLat1 = DegreesToRadians(lat1);
            double radLat2 = DegreesToRadians(lat2);
            double a = radLat1 - radLat2;
            double b = DegreesToRadians(lng1) - DegreesToRadians(lng2);

            double s = 2 * Math.Asin(Math.Sqrt(Math.Pow(Math.Sin(a / 2), 2) +
             Math.Cos(radLat1) * Math.Cos(radLat2) * Math.Pow(Math.Sin(b / 2), 2)));
            s = s * EARTH_RADIUS;
            s = Math.Round(s * 10000) / 10000;
            return s;
        }

        /// <summary>
        /// 求B点经纬度
        /// </summary>
        /// <param name="A">已知点的经纬度</param>
        /// <param name="distance">AB两地的距离  单位km</param>
        /// <param name="angle">AB连线与正北方向的夹角(0~360)</param>
        /// <returns>B点的经纬度</returns>
        public static CustomLatLng GetCustomLatLng(CustomLatLng A, double distance, double angle)
        {
            double dx = distance * 1000 * Math.Sin(DegreesToRadians(angle));
            double dy = distance * 1000 * Math.Cos(DegreesToRadians(angle));

            double bjd = (dx / A.Ed + A.m_RadLo) * 180 / Math.PI;
            double bwd = (dy / A.Ec + A.m_RadLa) * 180 / Math.PI;
            return new CustomLatLng(bjd, bwd);
        }

        //度 转换成 弧度
        public static double DegreesToRadians(double degrees)
        {
            const double degToRadFactor = Math.PI / 180;
            return degrees * degToRadFactor;
        }

        /// <summary>
        /// 弧度 转换成 度
        /// </summary>
        /// <param name="radians"></param>
        /// <returns></returns>
        public static double RadiansToDegrees(double radians)
        {
            const double radToDegFactor = 180 / Math.PI;
            return radians * radToDegFactor;
        }

        /// <summary>
        /// 描述:以centerP为圆心,绘制半径为radius的圆
        /// </summary>
        /// <param name="gMapControl">overlay:图层</param>
        /// <param name="overlay">图层</param>
        /// <param name="centerP">圆心点</param>
        /// <param name="radius">圆半径(单位: km)</param>
        /// <param name="name">多边形id</param>
        public static GMapPolygon DrawEllipse(GMapControl gMapControl, GMapOverlay overlay, PointLatLng centerP, int radius, string name)
        {
            try
            {
                if (radius <= 0)
                    return null;

                List<PointLatLng> latLngs = new List<PointLatLng>();
                CustomLatLng centerLatLng = new CustomLatLng(centerP.Lng, centerP.Lat);

                // 0 - 360度 寻找半径为radius,圆心为centerP的圆上点的经纬度
                for (int i = 0; i < 360; i++)
                {
                    //获取目标经纬度
                    CustomLatLng tempLatLng = GetCustomLatLng(centerLatLng, radius, i);
                    //将自定义的经纬度类 转换成 标准经纬度类
                    PointLatLng p = new PointLatLng(tempLatLng.m_Latitude, tempLatLng.m_Longitude);
                    latLngs.Add(p);
                }

                //安全性检查
                if (latLngs.Count < 20)
                {
                    return null;
                }

                //通过绘制多边形的方式绘制圆
                GMapPolygon gpol = new GMapPolygon(latLngs, name);
                gpol.Stroke = new Pen(Color.Red, 1.0f);
                gpol.Fill = new SolidBrush(Color.FromArgb(20, Color.Red));
                gpol.IsHitTestVisible = true;
                overlay.Polygons.Add(gpol);
                return gpol;
            }
            catch (Exception ex)
            {
                return null;
            }

        }
    }

    //自定义经纬度类
    public class CustomLatLng
    {
        public double Rc = 6378137;     //赤道半径
        public double Rj = 6356725;     //极半径
        public double m_LoDeg, m_LoMin, m_LoSec;
        public double m_LaDeg, m_LaMin, m_LaSec;
        public double m_Longitude, m_Latitude;
        public double m_RadLo, m_RadLa;
        public double Ec;
        public double Ed;
        public CustomLatLng(double longitude, double latitude)
        {
            m_LoDeg = (int)longitude;
            m_LoMin = (int)((longitude - m_LoDeg) * 60);
            m_LoSec = (longitude - m_LoDeg - m_LoMin / 60) * 3600;

            m_LaDeg = (int)latitude;
            m_LaMin = (int)((latitude - m_LaDeg) * 60);
            m_LaSec = (latitude - m_LaDeg - m_LaMin / 60) * 3600;

            m_Longitude = longitude;
            m_Latitude = latitude;
            m_RadLo = longitude * Math.PI / 180;
            m_RadLa = latitude * Math.PI / 180;
            Ec = Rj + (Rc - Rj) * (90 - m_Latitude) / 90;
            Ed = Ec * Math.Cos(m_RadLa);
        }
    }
  • 用于绘制的接口定义IPlot
public interface IPlot
    {
        /// <summary>
        /// 是否为图元
        /// </summary>
        /// <returns></returns>
       bool IsPlot();

        /// <summary>
        /// 获取当前图元的点集数量
        /// </summary>
        /// <param name="value"></param>
        void SetPoints(PointLatLng[] value);

        /// <summary>
        /// 获取当前图元的点集
        /// </summary>
        /// <returns>图元的点集</returns>
        PointLatLng[] GetPoints();

        /// <summary>
        /// 获取当前图元的点集数量
        /// </summary>
        /// <returns>图元的点集的数量</returns>
        int GetPointCount();

        /// <summary>
        /// 更新某个索引的点
        /// </summary>
        /// <param name="point">点</param>
        /// <param name="index">位置</param>
        void UpdatePoint(PointLatLng point, int index);

        /// <summary>
        /// 更新最后一个点
        /// </summary>
        /// <param name="point">{ol.Coordinate}</param>
        void UpdateLastPoint(PointLatLng point);

        /// <summary>
        /// 图元绘制逻辑.各个图元用来覆盖
        /// </summary>
        void Generate();

        /// <summary>
        /// 图元结束绘制回调
        /// </summary>
        void FinishDrawing();
    }
  • Constants:常量定义
public class Constants
    {
        public static readonly double TWO_PI = Math.PI*2;
        public static readonly double HALF_PI = Math.PI / 2;
        public static readonly int FITTING_COUNT = 100;
        public static readonly double ZERO_TOLERANCE = 0.0001;
    }
  • PlotTypes:标绘图元枚举类型
public enum PlotTypes
    {
        [Description("点")]
        MARKER =1,
        [Description("折线")]
        POLYLINE =2,
        [Description("多边形")]
        POLYGON =3,
        [Description("圆")]
        CIRCLE =4,
        [Description("椭圆")]
        ELLIPSE=5,
        [Description("矩形")]
        RECTANGLE =6,
        [Description("弧线")]
        ARC =7,
        [Description("进攻箭头")]
        ATTACK_ARROW =8,
        [Description("曲线面")]
        CLOSED_CURVE =9,
        [Description("曲线")]
        CURVE =10,
        [Description("双箭头")]
        DOUBLE_ARROW =11,
        [Description("细直箭头")]
        FINE_ARROW =12,
        [Description("突击方向")]
        ASSAULT_DIRECTION =13,
        [Description("自由线")]
        FREEHAND_LINE =14,
        [Description("自由面")]
        FREEHAND_POLYGON =15,
        [Description("聚集地")]
        GATHERING_PLACE =16,
        [Description("弓形")]
        LUNE =17,
        [Description("扇形")]
        SECTOR =18,
        [Description("分队战斗行动")]
        SQUAD_COMBAT=19,
        [Description("直箭头")]
        STRAIGHT_ARROW =20,
        [Description("进攻箭头(尾)")]
        TAILED_ATTACK_ARROW =21,
        [Description("分队战斗行动(尾)")]
        TAILED_SQUAD_COMBAT =22,

    }
  • PlotUtil:这个是核心类库,包括箭头算法、贝塞尔曲线和B样条曲线等。这儿的箭头计算坐标是EPSG:4326下的坐标
public class PlotUtil
    {
        public static double ConvertDegreesToRadians(double degrees)
        {
            double radians = (Math.PI / 180) * degrees;
            return (radians);
        }

        public static double Distance(PointLatLng pnt1,PointLatLng pnt2)
        {
            return Math.Sqrt(Math.Pow((pnt1.Lng - pnt2.Lng), 2) + Math.Pow((pnt1.Lat - pnt2.Lat), 2));
        }

        public static double WholeDistance(List<PointLatLng> points)
        {
            double _distance = 0.0;
            for(int i = 0; i < points.Count-1; i++)
            {
                _distance += Distance(points[i], points[i + 1]);
            }
            return _distance;
        }

        public static double GetBaseLength(List<PointLatLng> points)
        {
            return Math.Pow(WholeDistance(points), 0.99);
        }

        public static PointLatLng Mid(PointLatLng pnt1,PointLatLng pnt2)
        {
            return new PointLatLng((pnt1.Lat + pnt2.Lat) / 2, (pnt1.Lng + pnt2.Lng) / 2);
        }

        public static PointLatLng GetCircleCenterOfThreePoints(PointLatLng pnt1, PointLatLng pnt2, PointLatLng pnt3)
        {
            var pntA = new PointLatLng((pnt1.Lat + pnt2.Lat) / 2, (pnt1.Lng + pnt2.Lng) / 2);
            var pntB = new PointLatLng(pntA.Lat + pnt1.Lng - pnt2.Lng, pntA.Lng - pnt1.Lat + pnt2.Lat);
            var pntC = new PointLatLng((pnt1.Lat + pnt3.Lat) / 2, (pnt1.Lng + pnt3.Lng) / 2);
            var pntD = new PointLatLng(pntC.Lat + pnt1.Lng - pnt3.Lng, pntC.Lng - pnt1.Lat + pnt3.Lat);
            return GetIntersectPoint(pntA, pntB, pntC, pntD);
        }

        public static PointLatLng GetIntersectPoint(PointLatLng pntA, PointLatLng pntB, PointLatLng pntC, PointLatLng pntD)
        {
            if (pntA.Lat == pntB.Lat)
            {
                var f = (pntD.Lng - pntC.Lng) / (pntD.Lat - pntC.Lat);
                var x = f * (pntA.Lat - pntC.Lat) + pntC.Lng;
                var y = pntA.Lat;
                return new PointLatLng(y, x);
            }
            if (pntC.Lat == pntD.Lat)
            {
                var e = (pntB.Lng - pntA.Lng) / (pntB.Lat - pntA.Lat);
                var x = e * (pntC.Lat - pntA.Lat) + pntA.Lng;
                var y = pntC.Lat;
                return new PointLatLng(y, x);
            }
            else
            {
                var e = (pntB.Lng - pntA.Lng) / (pntB.Lat - pntA.Lat);
                var f = (pntD.Lng - pntC.Lng) / (pntD.Lat - pntC.Lat);
                var y = (e * pntA.Lat - pntA.Lng - f * pntC.Lat + pntC.Lng) / (e - f);
                var x = e * y - e * pntA.Lat + pntA.Lng;
                return new PointLatLng(y, x);
            }
        }

        /// <summary>
        /// 获取方位角
        /// </summary>
        /// <param name="startPnt">起点</param>
        /// <param name="endPnt">终点</param>
        public static double GetAzimuth(PointLatLng startPnt,PointLatLng endPnt)
        {
            double azimuth = 0.0;
            var angle = Math.Asin(Math.Abs(endPnt.Lat - startPnt.Lat) / Distance(startPnt, endPnt));
            if (endPnt.Lat >= startPnt.Lat && endPnt.Lng >= startPnt.Lng)
                azimuth = angle + Math.PI;
            else if (endPnt.Lat >= startPnt.Lat && endPnt.Lng < startPnt.Lng)
                azimuth = Constants.TWO_PI - angle;
            else if (endPnt.Lat < startPnt.Lat && endPnt.Lng < startPnt.Lng)
                azimuth = angle;
            else if (endPnt.Lat < startPnt.Lat && endPnt.Lng >= startPnt.Lng)
                azimuth = Math.PI - angle;
            return azimuth;
        }

        public static double GetAngleOfThreePoints(PointLatLng pntA, PointLatLng pntB, PointLatLng pntC)
        {
            var angle = GetAzimuth(pntB, pntA) - GetAzimuth(pntB, pntC);
            return (angle < 0 ? angle + Constants.TWO_PI : angle);
        }

        public static bool IsClockWise(PointLatLng pnt1, PointLatLng pnt2, PointLatLng pnt3)
        {
            return ((pnt3.Lat - pnt1.Lat) * (pnt2.Lng - pnt1.Lng) > (pnt2.Lat - pnt1.Lat) * (pnt3.Lng - pnt1.Lng));
        }

        public static PointLatLng GetPointOnLine(double t,PointLatLng startPnt,PointLatLng endPnt)
        {
            var x = startPnt.Lng + (t * (endPnt.Lng - startPnt.Lng));
            var y = startPnt.Lat + (t * (endPnt.Lat - startPnt.Lat));
            return new PointLatLng(y,x);
        }

        public static PointLatLng GetCubicValue(double t, PointLatLng startPnt, PointLatLng cPnt1, PointLatLng  cPnt2, PointLatLng endPnt)
        {
            t = Math.Max(Math.Min(t, 1), 0);
            var tp = 1 - t;
            var t2 = t * t;
            var t3 = t2 * t;
            var tp2 = tp * tp;
            var tp3 = tp2 * tp;
            var x = (tp3 * startPnt.Lng) + (3 * tp2 * t * cPnt1.Lng) + (3 * tp * t2 * cPnt2.Lng) + (t3 * endPnt.Lng);
            var y = (tp3 * startPnt.Lat) + (3 * tp2 * t * cPnt1.Lat) + (3 * tp * t2 * cPnt2.Lat) + (t3 * endPnt.Lat);
            return new PointLatLng(y,x);
        }

        public static PointLatLng GetThirdPoint(PointLatLng startPnt, PointLatLng endPnt, double angle,double distance,bool clockWise)
        {
            var azimuth = GetAzimuth(startPnt, endPnt);
            var alpha = clockWise ? azimuth + angle : azimuth - angle;
            var dx = distance * Math.Cos(alpha);
            var dy = distance * Math.Sin(alpha);
            return new PointLatLng(endPnt.Lat + dy,endPnt.Lng + dx);
        }

        public static List<PointLatLng> GetArcPoints(PointLatLng center, double radius, double startAngle, double endAngle)
        {
            List<PointLatLng> pnts = new List<PointLatLng>();
            var angleDiff = endAngle - startAngle;
            angleDiff = angleDiff < 0 ? angleDiff + Constants.TWO_PI : angleDiff;
            for (var i = 0; i <= Constants.FITTING_COUNT; i++)
            {
                var angle = startAngle + angleDiff * i / Constants.FITTING_COUNT;
                double x = center.Lng + radius * Math.Cos(angle);
                double y = center.Lat + radius * Math.Sin(angle);
                pnts.Add(new PointLatLng(y,x));
            }
            return pnts;
        }

        public static PointLatLng[] GetBisectorNormals(double t, PointLatLng pnt1, PointLatLng pnt2, PointLatLng pnt3)
        {
            var normal = GetNormal(pnt1, pnt2, pnt3);
            var dist = Math.Sqrt(normal.Lng * normal.Lng + normal.Lat * normal.Lat);
            var uX = normal.Lng / dist;
            var uY = normal.Lat / dist;
            var d1 = Distance(pnt1, pnt2);
            var d2 = Distance(pnt2, pnt3);
            PointLatLng bisectorNormalRight;
            PointLatLng bisectorNormalLeft;

            if (dist > Constants.ZERO_TOLERANCE)
            {
                if (IsClockWise(pnt1, pnt2, pnt3))
                {
                    var dt = t * d1;
                    var x = pnt2.Lng - dt * uY;
                    var y = pnt2.Lat + dt * uX;
                    bisectorNormalRight = new PointLatLng(y, x);
                    dt = t * d2;
                    x = pnt2.Lng + dt * uY;
                    y = pnt2.Lat - dt * uX;
                    bisectorNormalLeft = new PointLatLng(y, x);
                }
                else
                {
                    var dt = t * d1;
                    var x = pnt2.Lng + dt * uY;
                    var y = pnt2.Lat - dt * uX;
                    bisectorNormalRight = new PointLatLng(y, x);
                    dt = t * d2;
                    x = pnt2.Lng - dt * uY;
                    y = pnt2.Lat + dt * uX;
                    bisectorNormalLeft = new PointLatLng(y, x);
                }
            }
            else
            {
                var x = pnt2.Lng + t * (pnt1.Lng - pnt2.Lng);
                var y = pnt2.Lat + t * (pnt1.Lat - pnt2.Lat);
                bisectorNormalRight = new PointLatLng(y, x);
                x = pnt2.Lng + t * (pnt3.Lng - pnt2.Lng);
                y = pnt2.Lat + t * (pnt3.Lat - pnt2.Lat);
                bisectorNormalLeft = new PointLatLng(y, x);
            }
            return new PointLatLng[] { bisectorNormalRight, bisectorNormalLeft };
        }

        public static PointLatLng GetNormal(PointLatLng pnt1, PointLatLng pnt2, PointLatLng pnt3)
    {
        var dX1 = pnt1.Lng - pnt2.Lng;
        var dY1 = pnt1.Lat - pnt2.Lat;
        var d1 = Math.Sqrt(dX1 * dX1 + dY1 * dY1);
        dX1 /= d1;
        dY1 /= d1;

        var dX2 = pnt3.Lng - pnt2.Lng;
        var dY2 = pnt3.Lat- pnt2.Lat;
        var d2 = Math.Sqrt(dX2 * dX2 + dY2 * dY2);
        dX2 /= d2;
        dY2 /= d2;

        var uX = dX1 + dX2;
        var uY = dY1 + dY2;
        return new PointLatLng(uY,uX);
    }

        public static List<PointLatLng> GetCurvePoints(double t, PointLatLng[] controlPoints)
        {
            var leftControl = GetLeftMostControlPoint(t, controlPoints);
            List<PointLatLng> normals = new List<PointLatLng> { leftControl };
            for (int i = 0; i < controlPoints.Length - 2; i++)
            {
                var pnt1 = controlPoints[i];
                var pnt2 = controlPoints[i + 1];
                var pnt3 = controlPoints[i + 2];
                var normalPoints = GetBisectorNormals(t, pnt1, pnt2, pnt3);
                normals = normals.Concat(normalPoints).ToList();
            }
            var rightControl = GetRightMostControlPoint(t,controlPoints);
            normals.Add(rightControl);
            List<PointLatLng> points = new List<PointLatLng>();
            for (int i = 0; i < controlPoints.Length - 1; i++)
            {
                var pnt1 = controlPoints[i];
                var pnt2 = controlPoints[i + 1];
                points.Add(pnt1);
                for(var k = 0.0; k < Constants.FITTING_COUNT; k++)
                {
                    PointLatLng pnt = GetCubicValue(k / Constants.FITTING_COUNT, pnt1, normals[i * 2], normals[i * 2 + 1], pnt2);
                    points.Add(pnt);
                }
                points.Add(pnt2);
            }
            return points;
        }

        public static PointLatLng GetLeftMostControlPoint(double t, PointLatLng[] controlPoints)
        {
            var pnt1 = controlPoints[0];
            var pnt2 = controlPoints[1];
            var pnt3 = controlPoints[2];
            var pnts = GetBisectorNormals(0, pnt1, pnt2, pnt3);
            var normalRight = pnts[0];
            var normal = GetNormal(pnt1, pnt2, pnt3);
            var dist = Math.Sqrt(normal.Lng * normal.Lng + normal.Lat * normal.Lat);
            double controlX = 0;
            double controlY = 0;
            if (dist > Constants.ZERO_TOLERANCE)
            {
                var arr_mid = Mid(pnt1, pnt2);
                var pX = pnt1.Lng - arr_mid.Lng;
                var pY = pnt1.Lat - arr_mid.Lat;

                var d1 = Distance(pnt1, pnt2);
                // normal at midpoint
                var n = 2.0 / d1;
                var nX = -n * pY;
                var nY = n * pX;

                // upper triangle of symmetric transform matrix
                var a11 = nX * nX - nY * nY;
                var a12 = 2 * nX * nY;
                var a22 = nY * nY - nX * nX;

                var dX = normalRight.Lng - arr_mid.Lng;
                var dY = normalRight.Lat - arr_mid.Lat;

                // coordinates of reflected vector
                controlX = arr_mid.Lng + a11 * dX + a12 * dY;
                controlY = arr_mid.Lat + a12 * dX + a22 * dY;
            }
            else
            {
                controlX = pnt1.Lng + t * (pnt2.Lng - pnt1.Lng);
                controlY = pnt1.Lat + t * (pnt2.Lat - pnt1.Lat);
            }
            return new PointLatLng(controlY, controlX);
        }

        public static PointLatLng GetRightMostControlPoint(double t, PointLatLng[] controlPoints)
        {
            var count = controlPoints.Length;
            var pnt1 = controlPoints[count - 3];
            var pnt2 = controlPoints[count - 2];
            var pnt3 = controlPoints[count - 1];
            var pnts = GetBisectorNormals(0, pnt1, pnt2, pnt3);
            var normalLeft = pnts[1];
            var normal = GetNormal(pnt1, pnt2, pnt3);
            var dist = Math.Sqrt(normal.Lng * normal.Lng + normal.Lat * normal.Lat);
            double controlX = 0;
            double controlY = 0;
            if (dist > Constants.ZERO_TOLERANCE)
            {
                var arr_mid = Mid(pnt2, pnt3);
                var pX = pnt3.Lng - arr_mid.Lng;
                var pY = pnt3.Lat - arr_mid.Lat;

                var d1 = Distance(pnt2, pnt3);
                // normal at midpoint
                var n = 2.0 / d1;
                var nX = -n * pY;
                var nY = n * pX;

                // upper triangle of symmetric transform matrix
                var a11 = nX * nX - nY * nY;
                var a12 = 2 * nX * nY;
                var a22 = nY * nY - nX * nX;

                var dX = normalLeft.Lng - arr_mid.Lng;
                var dY = normalLeft.Lat - arr_mid.Lat;

                // coordinates of reflected vector
                controlX = arr_mid.Lng + a11 * dX + a12 * dY;
                controlY = arr_mid.Lat + a12 * dX + a22 * dY;
            }
            else
            {
                controlX = pnt3.Lng + t * (pnt2.Lng - pnt3.Lng);
                controlY = pnt3.Lat + t * (pnt2.Lat - pnt3.Lat);
            }
            return new PointLatLng(controlY, controlX);
        }

        public static List<PointLatLng> GetBezierPoints(PointLatLng[] points)
        {
            if (points.Length <= 2)
                return points.ToList();

            List<PointLatLng> bezierPoints = new List<PointLatLng>();
            var n = points.Length - 1;
            for (var t = 0.0; t <= 1; t += 0.01)
            {
                double x = 0.0;
                double y = 0.0;
                for (var index = 0; index <= n; index++)
                {
                    var factor = GetBinomialFactor(n, index);
                    var a = Math.Pow(t, index);
                    var b = Math.Pow((1 - t), (n - index));
                    x += factor * a * b * points[index].Lng;
                    y += factor * a * b * points[index].Lat;
                }
                bezierPoints.Add(new PointLatLng(y, x));
            }
            bezierPoints.Add(points[n]);
            return bezierPoints;
        }

        public static double GetBinomialFactor(double n, double index)
        {
            return GetFactorial(n) / (GetFactorial(index) * GetFactorial(n - index));
        }

        public static double GetFactorial(double n)
        {
            if (n <= 1)
                return 1;
            if (n == 2)
                return 2;
            if (n == 3)
                return 6;
            if (n == 4)
                return 24;
            if (n == 5)
                return 120;
            var result = 1;
            for (var i = 1; i <= n; i++)
                result *= i;
            return result;
        }

        public static List<PointLatLng> GetQBSplinePoints(PointLatLng[] points)
        {
            if (points.Length <= 2)
                return points.ToList();

            int n = 2;

            List<PointLatLng> bSplinePoints = new List<PointLatLng>();
            int m = points.Length - n - 1;
            bSplinePoints.Add(points[0]);
            for (int i = 0; i <= m; i++)
            {
                for (var t = 0.0; t <= 1; t += 0.05)
                {
                    var x = 0.0;
                    var y = 0.0;
                    for (var k = 0; k <= n; k++)
                    {
                        var factor = GetQuadricBSplineFactor(k, t);
                        x += factor * points[i + k].Lng;
                        y += factor * points[i + k].Lat;
                    }
                    bSplinePoints.Add(new PointLatLng(y, x));
                }
            }
            bSplinePoints.Add(points[points.Length - 1]);
            return bSplinePoints;
        }

        public static double GetQuadricBSplineFactor(int k, double t)
        {
            if (k == 0)
                return Math.Pow(t - 1, 2) / 2;
            if (k == 1)
                return (-2 * Math.Pow(t, 2) + 2 * t + 1) / 2;
            if (k == 2)
                return Math.Pow(t, 2) / 2;
            return 0;
        }
    }
  • 全部接口和帮助类全部从js到C#重写完毕,重写设计以后不合理的地方应该蛮多。有能力的可以自己设计重写。下一篇将开始写对应图元构建类。
  • 刚开始写博客,存在太多问题,希望多多包涵。

 

 

 

版权声明:本文为博主原创文章,遵循 CC 4.0 BY-SA 版权协议,转载请附上原文出处链接和本声明。
本文链接:https://blog.csdn.net/qq_17371831/article/details/107881654

智能推荐

Linux/Ubuntu(三)proxychains代理软件_linux下的代理软件-程序员宅基地

文章浏览阅读432次。配置代理安装proxychainssudo apt install proxychainsProxyChains 的配置文件位于/etc/proxychains.conf,打开后你在末尾添加代理,sock4改为sock5。例如:socks5 127.0.0.1 1080在需要代理的命令前加上proxychains4或启动bash$ proxychains ..._linux下的代理软件

HashMap的时间复杂度是O(1)_6.hashmap中的get方法执行情况,解释get是o1复杂度-程序员宅基地

文章浏览阅读8k次。今天在面试的时候说到HashMap,面试官问了这么一个问题:你说HashMap的get迭代了一个链表,那怎么保证HashMap的时间复杂度O(1)?链表的查找的时间复杂度又是多少? 在这之前我是阅读过HashMap的源码的:Java7源码浅析——对HashMap的理解 由上一个博客..._6.hashmap中的get方法执行情况,解释get是o1复杂度

分类算法————sklearn转换器和估计器_怎么用sklearn的ridge类构造估计器-程序员宅基地

文章浏览阅读831次。目录3.2.1 转换器-特征工程的父类3.2.2 估计器(sklearn机器学习算法的实现)1 转换器-特征工程的父类把特征工程的接口称之为转换器 fit_transform() 两个函数的封装 fit 做计算 transform 进行最终的转换2 估计器(sklearn机器学习算法的实现)1、用于分类的估计器: sklearn.neighbors k-近邻算法 sklearn.naive_bayes 贝叶斯 s..._怎么用sklearn的ridge类构造估计器

python机器学习--数据预处理API_机器学习api要做数据处理吗-程序员宅基地

文章浏览阅读333次。均值移除(标准差)由于一个样本的不同特征值差异较大,不利于使用现有机器学习算法进行样本处理。均值移除可以让样本矩阵中的每一列的平均值为0,标准差为1。#均值移除APIimport sklearn.preprocessing as spimport numpy as np# scale函数用于对函数进行预处理,实现均值移除。# array为原数组,返回A为均值移除后的结果。A = sp..._机器学习api要做数据处理吗

美国选举问题/完全背包/Knapsack-程序员宅基地

文章浏览阅读160次。using System;using System.Collections.Generic;using System.Linq;using System.Text;namespace Knapsack{ /// <summary> /// Project:Knapsack /// Author: Andy Zeng...

Cocos2d-X 移植 Android 平台详解_cocos2d-x-android-程序员宅基地

文章浏览阅读534次。在 Windows 下使用 VS 开发 Win32 平台的游戏程序后,需要把它们移植到不同的平台上。在这里首先介绍 Android 平台的移植工作。因为 Windows 和 手机等平台的 CPU 指令不同、架构不同,在 Windows 下编写的程序要想在手机上运行,首先要进行编译。在编译之前需要准备以下软件: Android NDK Android SDK Ap..._cocos2d-x-android

随便推点

anaconda安装所有库代码集总_anaconda代码下载安装所有库-程序员宅基地

文章浏览阅读1.3k次。1.安装jupyterconda install jupyter_anaconda代码下载安装所有库

hadoop2.2.0 mapred-site.xml--i/o properties_hadoop 2.10.0 mapred.site.xml配置-程序员宅基地

文章浏览阅读628次。property> name>mapreduce.task.io.sort.factorname> value>10value> description>The number of streams to merge at once while sorting files. This determines the number of open file handles.de_hadoop 2.10.0 mapred.site.xml配置

NodeJS如何免费发送邮件,使用Nodemailer发送电子邮件_nodemailer 怎收费-程序员宅基地

文章浏览阅读900次。案例:成功提交表单后,您需要发送确认电子邮件,或在每次购买时提供收据或订单详细信息。使用Nodemailer可以在 nodejs 本身中免费完成。_nodemailer 怎收费

漫谈程序员系列:一张图道尽程序员的出路_98tang.net-程序员宅基地

文章浏览阅读10w+次,点赞151次,收藏353次。今天的选择,决定了明天的路……_98tang.net

[VM] 如何给ESXi 上的VM缩小硬盘(VMDK)_esxi 减小vmdk-程序员宅基地

文章浏览阅读1.3w次。问题描述:最近遇到了一个需要给ESXi host上的VM缩小硬盘的问题,本以为和单机版vmware类似,把guest os里面的disk shrink出来,然后在vmware里面减小size就行了。但是:it’s not supported.给一个VMware VM缩小硬盘的官方方法是使用VMware vCenter Converter。ShrinkingVirtual disk shrinking is supported when using VMware Converter conve_esxi 减小vmdk

梯度反传-程序员宅基地

文章浏览阅读220次。pytorch 梯度反传详细说明:import torchw = torch.tensor([1.], requires_grad=True)x = torch.tensor([2.], requires_grad=True)a = torch.add(w, x)b = torch.add(w, 1)y = torch.mul(a, b)print(w)print(x)print(a)print(b)print(y) y.backward(retain_graph=True_梯度反传