Текущий архив: 2005.12.04;
Скачать: CL | DM;

Вниз

Получить координаты точки по удалению и азимуту   Найти похожие ветки 

← →
demon_god ©   (2005-11-18 01:10) [0]

Доброе время суток!

Задача такая: Есть некая точка с координатами Х с.ш. У в.д. Как узнать широту/долготу точки, расположенной на расстоянии R в направлении азимута A?

Подскажите алгоритм, пожалуйста или исходник.

---

← →
PKT   (2005-11-18 02:55) [1]

У меня была обратная задача:
есть точки с координатами Phi1, Lambda1
и Phi2, Lambda2 (широта-долгота).
По ним найти расстояние и азимут направления.
Может быть текст чем-то поможет.
На приблуды типа GetGMSStr(...) не обращай внимания -
это только перевод в текстовый вид.


procedure TForm1.Calculation;
const
 R = 6371.118;     // средний радиус Земли (по Красовскому)
type
 TCoordPoint = record
   X: Double;
   Y: Double;
   Z: Double
 end;

 function PlaneDegrees(A, B, C, D, E, F: TCoordPoint): Double;

   procedure Determinant(P1, P2, P3: TCoordPoint; var A, B, C: Double);
   begin
     A := (P2.Y-P1.Y) * (P3.Z-P1.Z) - (P3.Y-P1.Y) * (P2.Z-P1.Z);
     B := (P3.X-P1.X) * (P2.Z-P1.Z) - (P2.X-P1.X) * (P3.Z-P1.Z);
     C := (P2.X-P1.X) * (P3.Y-P1.Y) - (P3.X-P1.X) * (P2.Y-P1.Y);
   end;

 var
   A1, B1, C1: Double;
   A2, B2, C2: Double;
 begin
   Determinant(A, B, C, A1, B1, C1);
   Determinant(D, E, F, A2, B2, C2);
   Result := ArcCos((A1*A2 + B1*B2 + C1*C2) /
     (Sqrt(A1*A1 + B1*B1 + C1*C1) * Sqrt(A2*A2 + B2*B2 + C2*C2)));
 end;

 function CoordXYZ(Lambda, Phi: Double): TCoordPoint;
 begin
   Result.X := R * cos(Lambda) * cos(Phi);
   Result.Y := R * sin(Lambda) * cos(Phi);
   Result.Z := R * sin(Phi);
 end;

 function Hemisphere(Index: Integer): Integer;
 begin
  {case Index of
     0: Result := 1;  // северная широта или восточная долгота
     1: Result :=-1;  // южная широта или западная долгота
   end;}
   Result := Index * (-2) + 1;
 end;

const
 Eps = 3.14/180/60/60/10;
var
 Phi1, Phi2: Real;
 Lambda1, Lambda2: Real;
 ro: Real;
 Theta: Double;
 Pt1, Pt2, Pt1E, PtZ: TCoordPoint;
begin
 if FLockCalc then Exit;
 try
   Phi1 := Hemisphere(rgPhi1.ItemIndex) *
     (edPhiG1.Number + edPhiM1.Number / 60 + edPhiS1.Number / 60 / 60);
   Phi2 := Hemisphere(rgPhi2.ItemIndex) *
     (edPhiG2.Number + edPhiM2.Number / 60 + edPhiS2.Number / 60 / 60);
   Lambda1 := Hemisphere(rgLam1.ItemIndex) *
     (edLamG1.Number + edLamM1.Number / 60 + edLamS1.Number / 60 / 60);
   Lambda2 := Hemisphere(rgLam2.ItemIndex) *
     (edLamG2.Number + edLamM2.Number / 60 + edLamS2.Number / 60 / 60);

   // Перевод в радианы
   Lambda1 := Pi / 180 * Lambda1;
   Lambda2 := Pi / 180 * Lambda2;
   Phi1 := Pi / 180 * Phi1;
   Phi2 := Pi / 180 * Phi2;

   ro := arccos(cos(Phi1) * cos(Phi2) * (cos(Lambda1) * cos(Lambda2) +
     sin(Lambda1) * sin(Lambda2)) + sin(Phi1) * sin(Phi2));
   Distance := R * ro;

   if (Abs(Lambda1 - Lambda2) < 1e-6) and (Abs(Phi1 - Phi2) < 1e-6) then
     edAzimuth.Text := "";

     Pt1 := CoordXYZ(Lambda1, Phi1);
     Pt2 := CoordXYZ(Lambda2, Phi2);
     Pt1E := CoordXYZ(Lambda1, Phi1 + Eps);
     FillChar(PtZ, SizeOf(PtZ), 0);

     Theta := PlaneDegrees(Pt1, Pt1E, PtZ, Pt1, Pt2, PtZ);

     if Lambda2 < 0 then Theta := 2 * Pi - Theta;
     edAzimuth.Text := GetGMSStr(gmsSecond, Theta * 180 / Pi, 0);
   end;

 except end;
end;


---

← →
demon_god ©   (2005-11-19 02:48) [2]

Благодарю за ответ, PKT. Попробую разобраться. Получится - выложу.

---

← →
palva ©   (2005-11-19 08:29) [3]

Советую сначала разобраться, что такое теорема косинусов в сферической геометрии.

---

Страницы: 1 вся ветка

Текущий архив: 2005.12.04;
Скачать: CL | DM;

Наверх

Память: 0.48 MB
Время: 0.031 c