public class Logic
{
private static readonly double eps = 1e-10;
public static Point FindFermat(List<Point> poly, double step = 10)
{
if (poly == null || poly.Count == 0) return null;
if (poly.Count == 1) return poly[0];
if (poly.Count == 2) return new Point() { x = (poly[0].x + poly[1].x) / 2, y = (poly[0].y + poly[1].y) / 2 };
Point basePoint = new Point() { x = poly[0].x, y = poly[0].y };
double ans = AllDistance(poly[0].x, poly[0].y, poly);
double a = poly[0].x, b = poly[0].y;
double ta = a, tb = b, temp = 0;
while (step > eps)
{
bool flag = true;
while (flag)
{
flag = false;
temp = AllDistance(a + step, b, poly);
if (temp < ans)
{
flag = true;
ans = temp;
ta = a + step;
tb = b;
}
temp = AllDistance(a - step, b, poly);
if (temp < ans)
{
flag = true;
ans = temp;
ta = a - step;
tb = b;
}
temp = AllDistance(a, b + step, poly);
if (temp < ans)
{
flag = true;
ans = temp;
ta = a;
tb = b + step;
}
temp = AllDistance(a, b - step, poly);
if (temp < ans)
{
flag = true;
ans = temp;
ta = a;
tb = b - step;
}
a = ta;
b = tb;
}
step *= 0.5;
}
return new Point() { x = a, y = b };
}
private Point Move(Point p, Point direct)
{
return new Point()
{
x = p.x + direct.x,
y = p.y + direct.y
};
}
private static double AllDistance(double tx, double ty, List<Point> pList)
{
double ans = 0;
foreach (var point in pList)
{
ans += Distance(new Point() { x = tx, y = ty }, point);
}
return ans;
}
private static double Distance(Point a, Point b)
{
return Math.Sqrt((a.x - b.x) * (a.x - b.x) + (a.y - b.y) * (a.y - b.y));
}
private const double EARTH_RADIUS = 6378.137;
private static double rad(double d)
{
return d * Math.PI / 180.0;
}
/// <summary>
/// 点到点间的距离
/// </summary>
/// <param name="lat1"></param>
/// <param name="lng1"></param>
/// <param name="lat2"></param>
/// <param name="lng2"></param>
/// <returns></returns>
public static double GetDistance(double lat1, double lng1, double lat2, double lng2)
{
double radLat1 = rad(lat1);
double radLat2 = rad(lat2);
double a = radLat1 - radLat2;
double b = rad(lng1) - rad(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;
}
}