Найти координаты пересечения двух кругов на земле? - PullRequest
0 голосов
/ 15 ноября 2018

Я пытаюсь найти вторую точку пересечения двух кругов. Одна из точек, которые я уже знаю, использовалась для расчета расстояния, а затем использовалась в качестве радиуса круга ( например ). Проблема в том, что я не знаю точку, я получаю две новые координаты, даже если они похожи. Возможно, проблема связана с искривлением Земли, но я искал какое-то решение и ничего не нашел.

Радиус окружностей рассчитывается с учетом кривизны земли. И вот код, который у меня есть:

function GET_coordinates_of_circles(position1,r1, position2,r2) {
  var deg2rad = function (deg) { return deg * (Math.PI / 180); };
  x1=position1.lng;
  y1=position1.lat;
  x2=position2.lng;
  y2=position2.lat;
  var centerdx = deg2rad(x1 - x2); 
  var centerdy = deg2rad(y1 - y2); 
  var R = Math.sqrt(centerdx * centerdx + centerdy * centerdy);

  if (!(Math.abs(r1 - r2) <= R && R <= r1 + r2)) { // no intersection
    console.log("nope");
    return []; // empty list of results
  }

  // intersection(s) should exist
  var R2 = R*R;
  var R4 = R2*R2;
  var a = (r1*r1 - r2*r2) / (2 * R2);
  var r2r2 = (r1*r1 - r2*r2);
  var c = Math.sqrt(2 * (r1*r1 + r2*r2) / R2 - (r2r2 * r2r2) / R4 - 1);

  var fx = (x1+x2) / 2 + a * (x2 - x1);
  var gx = c * (y2 - y1) / 2;
  var ix1 = fx + gx;
  var ix2 = fx - gx;

  var fy = (y1+y2) / 2 + a * (y2 - y1);
  var gy = c * (x1 - x2) / 2;
  var iy1 = fy + gy;
  var iy2 = fy - gy;
  // note if gy == 0 and gx == 0 then the circles are tangent and there is only one solution
  // but that one solution will just be duplicated as the code is currently written
  return [[iy1, ix1], [iy2, ix2]];
}

Переменная deg2rad, которая должна корректировать другие расчеты с кривизной земли.

Спасибо за любую помощь.

1 Ответ

0 голосов
/ 16 ноября 2018

Ваши вычисления для R и т. Д. Неверны, потому что плоская формула Пифагора не работает для сферической тригонометрии (например - у нас может быть треугольник со всеми тремя прямыми углами на сфере!).Вместо этого мы должны использовать специальные формулы.Некоторые из них взяты из этой страницы .

Сначала найдите большие круговые дуги в радианах для обоих радиусов, используя R = Earth radius = 6,371km

a1 = r1 / R
a2 = r2 / R

И расстояние (снова дугав радианах) между центром круга, используя формулу haversine

var R = 6371e3; // metres
var φ1 = lat1.toRadians();
var φ2 = lat2.toRadians();
var Δφ = (lat2-lat1).toRadians();
var Δλ = (lon2-lon1).toRadians();

var a = Math.sin(Δφ/2) * Math.sin(Δφ/2) +
        Math.cos(φ1) * Math.cos(φ2) *
        Math.sin(Δλ/2) * Math.sin(Δλ/2);
var ad = 2 * Math.atan2(Math.sqrt(a), Math.sqrt(1-a));

и ориентиром из положения 1 в положение 2:

 //where    φ1,λ1 is the start point, φ2,λ2 the end point 
 //(Δλ is the difference in longitude)
var y = Math.sin(λ2-λ1) * Math.cos(φ2);
var x = Math.cos(φ1)*Math.sin(φ2) -
        Math.sin(φ1)*Math.cos(φ2)*Math.cos(λ2-λ1);
var brng = Math.atan2(y, x);

Теперь посмотрите на картинку из моего ответа с учетом случая равных радиусов.(Здесь радиусы окружности могут быть разными, и мы должны использовать другой подход, чтобы найти необходимые дуги)

enter image description here

У нас есть сферические прямоугольные треугольники ACB и FCB (аналогично плоскому случаю BD перпендикулярен AF в точке C, а угол BCA является прямым).
Сферическая теорема Пифагора (из книги по sph. trig) говорит, что

 cos(AB) = cos(BC) * cos(AC)
 cos(FB) = cos(BC) * cos(FC)

или (используя x дляAC, y для BC и (ad-x) для FC)

 cos(a1) = cos(y) * cos(x)
 cos(a2) = cos(y) * cos(ad-x)

делим уравнения для исключения cos (y)

 cos(a1)*cos(ad-x) = cos(a2) * cos(x)
 cos(a1)*(cos(ad)*cos(x) + sin(ad)*sin(x)) = cos(a2) * cos(x)
 cos(ad)*cos(x) + sin(ad)*sin(x) = cos(a2) * cos(x) / cos(a1)
 sin(ad)*sin(x) = cos(a2) * cos(x) / cos(a1) - cos(ad)*cos(x)
 sin(ad)*sin(x) = cos(x) * (cos(a2) / cos(a1) - cos(ad))
 TAC = tg(x) = (cos(a2) / cos(a1) - cos(ad)) / sin(ad)

Имея гипотенузу и катет треугольника ACB, мы можем найтиугол между направлениями AC и AB (правила Нейпира для правильных сферических треугольников) - обратите внимание, что мы уже знаем TAC = tg(AC) и a1 = AB

cos(CAB)= tg(AC) * ctg(AB)
CAB = Math.acos(TAC * ctg(a1))

Теперь мы можем вычислить точки пересечения - они лежат на расстоянии дуги a1 от позиции 1вдоль подшипников brng-CAB и brng+CAB

B_bearing = brng - CAB
D_bearing = brng + CAB

Координаты точек пересечения:

var latB = Math.asin( Math.sin(lat1)*Math.cos(a1) + 
              Math.cos(lat1)*Math.sin(a1)*Math.cos(B_bearing) );
var lonB = lon1.toRad() + Math.atan2(Math.sin(B_bearing)*Math.sin(a1)*Math.cos(lat1), 
                     Math.cos(a1)-Math.sin(lat1)*Math.sin(lat2));

и то же для D_bearing

широта, lB - в радианах

...