Stack Overflow на русском Asked on August 30, 2021
Знаю, что вопрос не новый, но у меня большие проблемы с реализацией, даже при успешном гуглении интересующих формул из-за того, что необходимо сопоставлять долготу и широту с радиусом заданным в метрах.
Есть два типа точек:
Первый тип точек являются центрами окружностей, у которых заданы интересующие сектора (на схеме показана одна точка (X1,Y1) у которой интересующий сектор задан синим цветом), второй тип точек это положение на местности некоторого объекта (X2,Y2), для каждой отдельно взятой задачи этот объект будет всего 1, но каждый раз разный.
Известен азимут сектора в градусах (на схеме это Альфа значения 0-360),
Радиус в метрах (на схеме R),
Ширина сектора в градусах (на схеме это Бета 0-360))
Необходимо определить все точки, которые относятся к первому типу, в сектора которых попадает точка второго типа
На первом этапе я отсек точки, которые не попадают в окружности описываемые координатами точек и радиусами секторов код на MS SQL, но решить задачу с секторами у меня не получается
DECLARE float @point_latitude = 57.2; DECLARE float @point_logitude = 65.5539; DECLARE int @max_radius = 1500
SELECT * FROM Test_table
WHERE 2*ASIN(SQRT(POWER(SIN(RADIANS((@point_latitude-latitude)/2)),2)+COS(RADIANS(@point_latitude))*COS(RADIANS(latitude))*POWER(SIN(RADIANS((@point_logitude-longitude)/2)),2)))*6378245) <= @max_radius
С sql-реализацией, я к сожалению, вам помочь не могу, но могу предложить реализацию на псевдокоде.
Пусть Вектор(X1, Y1, X2, Y2)
- это вектор из точки 1 в точку 2. Он может быть определен как точка с координатами (X2 - X1, Y2 - Y1)
.
Пусть Угол(X1, Y1)
- это азимут вектора. Азимут у вас считается как угол по часовой стрелке от северного направления. Если север - это луч OY+, то на тригонометрической окружности он совпадает с точкой pi/2
. Тогда азимут произвольного вектора можно вычислить как (pi/2 - ATAN2(X1, Y1) + 2*pi) % (2*pi)
. Далее его нужно перевести в градусы по формуле Угол_в_градусах = (Угол / 2pi) * 360
.
Теперь условием для нахождения сектора будет:
Выбрать все сектора такие, что
(X2 - X1)^2 + (Y2 - Y1)^2 <= R1^2
И
Alpha1 - Beta1 / 2 <= Угол(Вектор(X1, Y1, X2, Y2)) <= Alpha1 + Beta1 / 2
Где X1
и Y1
- центр окружности; X2
и Y2
- точка; R1
- радиус окружности; Alpha1
- азимут сектора; Beta1
- ширина сектора.
Upd: Спасибо MBo за дополнение. Может случиться такое, что Alpha1 - Beta1 / 2 < 0
или Alpha1 + Beta1 / 2 > 360
. Кажется это можно решить двумя дополнительными неравенствами, записав их к основному через ИЛИ: Alpha1 - Beta1 / 2 + 360 <= Угол <= 360
, 0 <= Угол <= Alpha1 - Beta1 / 2 - 360
. Но это нужно протестировать, я не полностью уверен, что данные условия покрывают все крайние случаи.
Correct answer by EzikBro on August 30, 2021
Подход, который обрабатывает ситуации с переходом сектора через ноль автоматически, и не боится больших секторов:
Квадрат расстояния:
dx = X2-X1
dy = Y2-Y1
r2 = dx^2 + dy^2
Если он меньше R^2
, то продолжаем
r = sqrt(r2)
ux = dx / r
uy = dy / r
Теперь смотрим, что угол, образуемый вектором на точку, отклоняется от среднего азимута alpha
менее, чем на beta/2
. Для этого достаточно найти скалярное произведение единичного вектора в направлении alpha
и единичного вектора на точку ux, uy
, и сравнить его с косинусом beta/2
if ux * cos(alpha) + uy * sin(alpha) >= cos(beta/2)
то точка внутри сектора
Answered by MBo on August 30, 2021
Наваял на Компонентном Паскале прототипчик. Два радара - (0,0) и (1,0), первый смотрит на Pi/4 с шириной охвата +/- pi/8, второй - на 3*Pi/4 с такой же широтой охвата. Радиус действия - 1.0. Мишени две - (0.5, 0.5) и (0.5, 0.5 + 0.1), - обе попадают в сектора обоих радаров.
Результат работы программы:
Target point 0 in sector of radar 0
Target point 1 in sector of radar 0
Target point 0 in sector of radar 1
Target point 1 in sector of radar 1
Исходник для BlackBox Component Builder:
MODULE DESLibSectors;
IMPORT M := Math, Log := StdLog;
TYPE
Point = RECORD
x, y : REAL;
END;
PropertyRadar = RECORD
radius, directionAngle, widenessAngle : REAL;
END;
PointVector = POINTER TO ARRAY OF Point;
RadarPropertiesVector = POINTER TO ARRAY OF PropertyRadar;
PROCEDURE Atan2 (x, y : REAL) : REAL;
VAR
result, pi, pi_2 : REAL;
BEGIN
pi := M.Pi();
pi_2 := pi/2;
IF x = 0 THEN
IF y > 0.0 THEN
RETURN pi_2;
ELSIF y < 0.0 THEN
RETURN -pi_2;
END;
END;
result := M.ArcTan (y/x);
IF x > 0 THEN
RETURN result;
ELSE
IF y >= 0 THEN
RETURN result + pi;
ELSE
RETURN result - pi;
END;
END;
END Atan2;
PROCEDURE IsPointInRadarSector (pointRadar, pointTarget : Point; radarProperty : PropertyRadar) : BOOLEAN;
VAR
result : BOOLEAN;
r, rx, ry, angle, angle_min, angle_max : REAL;
normalized : Point;
BEGIN
rx := pointTarget.x - pointRadar.x;
ry := pointTarget.y - pointRadar.y;
r := M.Sqrt(rx*rx + ry*ry);
IF r = 0.0 THEN
RETURN TRUE; (* если объект прямо на радаре - видимо он таки в секторе :) *)
END;
result := (r <= radarProperty.radius); (* дистанция меньше радиуса охвата радара *)
angle_min := radarProperty.directionAngle - 0.5*radarProperty.widenessAngle;
angle_max := radarProperty.directionAngle + 0.5*radarProperty.widenessAngle;
normalized.x := rx/r; (* после нормализации - равен косинусу угла на мишень *)
normalized.y := ry/r; (* после нормализации - равен синусу угла на мишень *)
angle := Atan2 (normalized.x, normalized.y);
result := result & (angle >= angle_min) & (angle <= angle_max); (* угол на мишень больше минимального угла радара и меньше максимального *)
RETURN result;
END IsPointInRadarSector;
PROCEDURE TestPoints (pointsRadar, pointsTarget : PointVector; radarProperties : RadarPropertiesVector);
VAR
indRadar, indTarget : INTEGER;
BEGIN
FOR indRadar := 0 TO LEN (pointsRadar) - 1 DO
FOR indTarget := 0 TO LEN (pointsTarget) - 1 DO
IF IsPointInRadarSector (pointsRadar[indRadar], pointsTarget[indTarget], radarProperties[indRadar]) THEN
Log.String ("Target point "); Log.Int (indTarget); Log.String (" in sector of radar "); Log.Int (indRadar); Log.Ln;
END;
END;
END;
END TestPoints;
PROCEDURE Run*;
VAR
numOfRadars, numOfTargets : INTEGER;
pointsRadar, pointsTarget : PointVector;
propertiesRadar : RadarPropertiesVector;
BEGIN
numOfRadars := 2;
numOfTargets := 2;
NEW (pointsRadar, 2);
pointsRadar[0].x := 0.0;
pointsRadar[0].y := 0.0;
pointsRadar[1].x := 1.0;
pointsRadar[1].y := 0.0;
NEW (propertiesRadar, 2);
propertiesRadar[0].radius := 1.0;
propertiesRadar[0].directionAngle := M.Pi()/4;
propertiesRadar[0].widenessAngle := M.Pi()/4;
propertiesRadar[1].radius := 1.0;
propertiesRadar[1].directionAngle := 3*M.Pi()/4;
propertiesRadar[1].widenessAngle := M.Pi()/4;
NEW (pointsTarget, 2);
pointsTarget[0].x := 0.5;
pointsTarget[0].y := 0.5;
pointsTarget[1].x := 0.5;
pointsTarget[1].y := 0.5 + 0.1;
TestPoints (pointsRadar, pointsTarget, propertiesRadar);
END Run;
END DESLibSectors.
DESLibSectors.Run
Answered by serg.tortilliani on August 30, 2021
Get help from others!
Recent Answers
Recent Questions
© 2024 TransWikia.com. All rights reserved. Sites we Love: PCI Database, UKBizDB, Menu Kuliner, Sharing RPP