點到折線最短距離所在點距離折線起點的累積距離
阿新 • • 發佈:2018-12-24
點到折線最短距離所在點 距離 折線起點 的累積距離
using System; using System.Collections.Generic; using System.Linq; using System.Text; using ESRI.ArcGIS.Geometry; using RGeos.Geometry namespace RGeos.Geometry { public class CulmulateDistance { /// <summary> /// 點到折線最短距離處距離折線起點的累積距離 /// </summary> /// <param name="P">任意一點</param> /// <param name="polyline">折線</param> /// <returns></returns> public static double CulmulateDist_Point_to_Polyline(IPoint P, IPolyline polyline) { ISegmentCollection segs = polyline as ISegmentCollection; double min = double.MaxValue; int segIndex = -1;//最短距離所在檔的索引 for (int i = 0; i < segs.SegmentCount; i++) { //點到每條線段的最短距離 double dis = Dist_Point_to_Segment(P, segs.get_Segment(i)); if (dis < min) { min = dis; segIndex = i;//取出最小的一個 } } double culmulateDis = 0; for (int i = 0; i < segs.SegmentCount; i++) { if (segIndex != i) { culmulateDis += segs.get_Segment(i).Length; } else { ISegment current = segs.get_Segment(i); Vector3d v = new Vector3d(); v.X = current.ToPoint.X - current.FromPoint.X; v.Y = current.ToPoint.Y - current.FromPoint.Y; Vector3d w = new Vector3d(); w.X = P.X - current.FromPoint.X; w.Y = P.Y - current.FromPoint.Y; double c1 = dot(w, v);//投影長度 if (c1 <= 0)//這種情況最短距離在該線段的起點處 { break; } double c2 = dot(v, v); if (c2 <= c1)//這種情況最短距離在該線段的終點處 { culmulateDis += Math.Sqrt(c2); break; } double b = c1 / c2; culmulateDis += Math.Sqrt(c1); IPoint Pb = new PointClass(); Pb.X = current.FromPoint.X + b * v.X; Pb.Y = current.FromPoint.Y + b * v.Y; break; } } return culmulateDis; } public static IPoint GetCentrePoint(IPolygon polygon) { IArea pArea = polygon as IArea; IPoint pt = new PointClass(); pt.X = pArea.Centroid.X; pt.Y = pArea.Centroid.Y; return pt; } public static double Dist_Point_to_Segment(IPoint P, ISegment S) { Vector3d v = new Vector3d(); v.X = S.ToPoint.X - S.FromPoint.X; v.Y = S.ToPoint.Y - S.FromPoint.Y; Vector3d w = new Vector3d(); w.X = P.X - S.FromPoint.X; w.Y = P.Y - S.FromPoint.Y; double c1 = dot(w, v); if (c1 <= 0) return d(P, S.FromPoint); double c2 = dot(v, v); if (c2 <= c1) return d(P, S.ToPoint); double b = c1 / c2; IPoint Pb = new PointClass(); Pb.X = S.FromPoint.X + b * v.X; Pb.Y = S.FromPoint.Y + b * v.Y; return d(P, Pb); } /// <summary> /// 向量的模 /// </summary> /// <param name="v"></param> /// <returns></returns> public static double norm(Vector3d v) { return Math.Sqrt(dot2(v, v)); // norm = length of vector } /// <summary> /// 2D數量積,點乘 /// </summary> /// <param name="u"></param> /// <param name="v"></param> /// <returns></returns> public static double dot2(Vector3d u, Vector3d v) { return ((u).X * (v).X + (u).Y * (v).Y); } public static double dot(Vector3d u, Vector3d v) { return ((u).X * (v).X + (u).Y * (v).Y + (u).Z * (v).Z); } public static double d(IPoint P, IPoint P1) { return Math.Sqrt((P1.X - P.X) * (P1.X - P.X) + (P1.Y - P.Y) * (P1.Y - P.Y)); } /// <param name="x">增量X</param> /// <param name="y">增量Y</param> /// <returns>象限角</returns> public static double GetQuadrantAngle(double x, double y) { double theta = Math.Atan(y / x); if (x > 0 && y == 0) return 0; if (x == 0 && y > 0) return Math.PI / 2; if (x < 0 && y == 0) return Math.PI; if (x == 0 && y < 0) return 3 * Math.PI / 2; if (x > 0 && y > 0) return theta; if (x > 0 && y < 0) return Math.PI * 2 + theta; if (x < 0 && y > 0) return theta + Math.PI; if (x < 0 && y < 0) return theta + Math.PI; return theta; } } }