Проверить, принадлежит ли какая-либо точка отрезка заданной (заданным) окружностям

#cpp #геометрия

Привет. Прошу помощи. В общем задача: на координатной плоскости задается отрезок
с началом в точке А и с концом в точке В. Также задаются круги (можно и 1 круг) с центрами
в точках M1,M2,M3,...M100 и с радиусами R1,R2,R3,...,R100. Нужно проверить, принадлежит
ли какая-либо точка отрезка заданной (заданным) окружностям.
В общем-то это задача из олимпиады, код компилируется на сервере. Мой код: http://ideone.com/J5MJEH 
double xA, yA, xB, yB;
int count = 0;
int n;

cin >> xA >> yA >> xB >> yB; // Координаты A и B отрезка: 
cin >> n; // Количество окружностей

double * * Arr = new double * [3]; // Задаем массив для кругов, в каждой строчке
которого будут хранится x,y,r круга
for (int i = 0; i < n; i++)
    Arr[i] = new double[3];

for (int i = 0; i < n; i++) {
    for (int j = 0; j < 3; j++) {
        cin >> Arr[i][j];
    }
}
double r12 = sqrt((xB - xA) * (xB - xA) + (yB - yA) * (yB - yA)); // Длина отрезка
от A до B
for (int i = 0; i < n; i++) {
    double r1 = sqrt((Arr[i][0] - xA) * (Arr[i][0] - xA) + (Arr[i][1] - yA) * (Arr[i][1]
- yA)); // Длина отрезка от т. M до A
    double r2 = sqrt((Arr[i][0] - xB) * (Arr[i][0] - xB) + (Arr[i][1] - yB) * (Arr[i][1]
- yB)); // Длина отрезка от т. M до B

    if (((r12 * r12 + r1 * r1) - r2 * r2) >= 0 && ((r12 * r12 + r2 * r2) - r1 * r1)
>= 0) { // Если перпендикуляр опущенный на прямую, проведенную через A u B лежит на
отрезке AB, то проверяем, равен или больше ли радиус окружности
        double d = (abs((yB - yA) * Arr[i][0] + (xB - xA) * Arr[i][1] + (xA * yB
- yA * xB))) / sqrt((yB - yA) * (yB - yA) + (xB - xA) * (xB - xA));
        if ((Arr[i][2] - d) >= 0) {
            count++;
        }
    } else { // Иначе, если перпендикуляр опущенный на прямую лежит вне отрезка,
то проверяем, равен или больше ли радиус одному из отрезков от точки М до А или В
        if ((Arr[i][2] - r1) >= 0 || (Arr[i][2] - r2) >= 0) {
            count++;
        }
    }
}
cout << count;

В комментариях кода написал, что,где, как нахожу.
То, как определить, падает ли перпендикуляр на отрезок, я взял отсюда: http://algolist.manual.ru/maths/geom/misc/perp.php,
а именно неравенства: B2 <= A2 + C2, C2 <= A2 + B2.
В общем-то программа работает, но, вероятно, не всегда правильно. Уже достаточно
много перепробовал тестовых примером, всё правильно определяет. Но, когда отправляю
решение на сервер олимпиады, сервер говорит: "Wrong Answer". Я все грешу на типы. Сможете
помочь?     

Ответы

Ответ 1

Вы не проверяете случай, когда часть отрезка лежит в одной окружности, а часть -
в другой. Например, A(10, 10), B(10, 20), M1(10, 10), R1(5), M2(10, 20), R2(5). тут
половина отрезка лежит в первой окружности, а половина - во второй.
http://ideone.com/TdMCuz
Я бы предложил следующий алгоритм. Для каждой окружности проверяем, лежит ли отрезок
в ней полностью (т. е. лежат ли оба конца отрезка внутри окружности). Если да, то ответ
- "да", выходим. Если оба конца отрезка вне окружности, переходим к следующей окружности.
Если один конец внутри, а другой - снаружи, то обезаем отрезок, оставляя только ту
его часть, что снаружи, и продолжаем проверки. Если после проверки всех окружностей
остался ненулевой участок отрезка, ответ - "нет".


Ответ 2

Проверку того, что какая-либо точка отрезка принадлежит окружности, можно заменить
на проверку всего лишь трёх точек: концов отрезка и проекции центра окружности на прямую,
содержащую отрезок (третью точку берём только если она принадлежит отрезку).
Вы же проверяете третью точку, но не проверяете края. То есть, теряете случай, когда
очень короткий отрезок пересекает границу окружности. (Плюс наверное где-то ошибки
реализации из-за очень «низкоуровневого» кода.)
Реализация не очень сложна (но придётся вспомнить геометрию за 9-ый класс):
bool GetProjection(const Segment& segment, const Point& source,
                   Point& projection)
{
    Vector direction = segment.end - segment.begin;
    double segmentLength = direction.Length();
    Vector unitDirection = direction / segmentLength;
    Vector beginToSource = source - segment.begin;
    double projectionLength = scalarProduct(unitDirection, beginToSource);
    if (projectionLength < 0 || projectionLength > segmentLength)
        return false;
    projection = segment.begin + unitDirection * projectionLength;
    return true;
}

bool IsPointInCircle(const Point& point, const Circle& circle)
{
    var radiusVector = point - circle.center;
    return radiusVector.Length() <= circle.radius;
}

bool AnyPointInCircle(const Segment& segment, const Circle& circle)
{
    Point projection;
    return
        IsPointInCircle(segment.begin, circle) ||
        IsPointInCircle(segment.end, circle) ||
        (GetProjection(segment, circle.center, projection) &&
            IsPointInCircle(projection, circle));
}

Классы Point, Vector и прочие определяются самоочевидным образом.


Ответ 3

На рисунке изображены два случая непересечения отрезком окружности:


Когда прямая не пересекает окружность.
Когда прямая пересекает окружность, а отрезок не пересекает.




Для каждой окружности (R,M) и отрезка AB можно применить следующий алгоритм: 


Вычисляем коэффициенты уравнения
(A+qB)/(1+q) on (R,M), или
(xa + q* xb - (1+q)* xm)2 + (ya + q* yb - (1+q)* ym)2 = (1+q)2*R2:
С2q2+2C1q+C0 = 0, где
C0 = (xa - xm)2 + (ya - ym)2 - R2,
C1 = (xa - xm) * (xb - xm) + (ya - ym) * (yb - ym) - R2,
C2 = (xb - xm)2 + (yb - ym)2 - R2.  
Если C0 <= 0, точка A попала в круг.
Если C2 <= 0, точка B попала в круг.
Если det = C12 - C0C2 <= 0, прямая АB с кругом не пересекается.
Вычисляем корни уравнения q1 и q2 (если есть).
Если q1 >= 0 ИЛИ q2 >= 0, то попали.
Иначе не попали.  


Программа:

function are_crossed($a, $b, $r, $m){
    $rr = $r*$r;        
    $am = array($a[0] - $m[0], $a[1] - $m[1]);      // вектор AM 
    $bm = array($b[0] - $m[0], $b[1] - $m[1]);      // вектор BM 
    $eq = array(
        $am[0]*$am[0] + $am[1]*$am[1] - $rr,
        $am[0]*$bm[0] + $am[1]*$bm[1] - $rr,
        $bm[0]*$bm[0] + $bm[1]*$bm[1] - $rr
    );
    $det = $eq[1]*$eq[1] - $eq[0]*$eq[2];

    if($eq[0] <= 0) return 1;                       // начало отрезка попало в круг
    if($eq[2] <= 0) return 2;                       // конец отрезка попал в круг
    if($det < 0) return -1;                         // линия не пересекает круг

    $d = sqrt($det); 
    $root1 = ( -$eq[1] - $d) / $eq[2];              // первая точка пересечения линии
и круга 
    $root2 = (-$eq[1]+$d)/$eq[2];                   // вторая точка пересечения линии
и круга

    if(($root1>=0) return 3;                        // первая точка принадлежит отрезку
    if(($root2>=0) return 4;                        // вторая точка принадлежит отрезку
    return -2;                                      // обе точки пересечения - вне
отрезка 
} 

$a = array(5, 20);
$b = array(25, 5);
$m = array(5, 5);
$r = 12;

$result = are_crossed($a, $b, $r, $m);
print "A = "; print_r($a);
print "
B = "; print_r($b);
print "
M = "; print_r($m);
print "
R = $r";
print("
result = $result
");
print($result>0 ? "пересекаются" : "не пересекаются");


Результаты:



A = Array ( [0] => 5 [1] => 20 ) 
B = Array ( [0] => 25 [1] => 5 ) 
M = Array ( [0] => 5 [1] => 5 ) 
R = 12
result = 3
пересекаются