계산 기하 Part1

도입

• 계산 기하(computational geometry) 알고리즘: 점, 선, 다각형과 원 등 각종 기하학적 도형을 다루는 알고리즘
  
• 계산 기하는 3d 그래픽이나 캐드, 로보틱스 등 다양한 분야에서 사용을 하기 때문에 전산학에서는 중요한 역할
• 많은 주제를 포함하고 있지만 학부 선형 대수나 고등학교 수준의 기하학을 요구
• 2차원 평면 도형과 3차원 입체 도형을 모두 다룸.(주로 2차원을 많이 사용)

계산 기하의 도구들

• 벡터
• 정의 및 특징
• 방향 + 크기 = 단위 벡터 + 스칼라(스칼라: 크기만 가지는 양)
• 크기와 방향이 같으면 동일한 벡터
• 계산기하의 기본적인 도구
• 2차원 평면 위에 있는 서로 다른 두 점의 상대적 위치

• 벡터의 종류

• 위치 벡터(Position Vector)

• 원점에서 임의의 점까지 향하는 벡터

• 그 점은 위치 좌표들에 의해 정의되어짐

• 일반적으로 측정점을 지칭하며, 보통 r로 표기

• 

  
  

• 위치 벡터는 주로 기준점으로부터 변위를 나타내는 벡터를 말함

• 변위: 방향을 고려한 거리

• 거리 벡터(Distance Vector)

• 두 점간의 거리를 나타내는 벡터 표기 : R = r-r'

• 

  
  

• 변위 벡터(Displacement Vector)

• 중간 경로에 관계없이 시작점과 끝점을 가장 짧게 연결하는 벡터

• 

  
  

• 법선 벡터(Normal Vector)

• 벡터 x에 수직(직교)하는 벡터 n

• nx = 0 (법선의 벡터 방벙식)

• 

  
  

• 방향 벡터(Direction Vector)
  

• x=tx(직선의 벡터 방정식)

• 

  
  

• 단위 벡터(Unit Vector)

• 크기가 1인 벡터

• ux = x/ ||x||

• ||ux|| = 1

• 벡터의 덧셈

• 삼각형법

• 평행사변형법

• 벡터의 뺄셈
  

• 벡터의 스칼라곱

• 벡터에 실수배를 곱하면 된다.

• 점, 직선, 선분의 표현

• 점: (x1, y1), (x2, y2)

• 선: 점과 점 사이를 잇는 선분

• 직선: 점들을 넘어서까지 잇는 선분

//c++ code #include <iostream> #include <cstdlib> #include <cstdio> #include <algorithm> #include <math.h> #include <vector> const double PI = 2.0 * acos(0.0); //2차원 벡터를 표현 struct vector2 { double x, y; //생성자를 explicit으로 지정하면 vector2를 젛을 곳에 잘못해 //실수가 들어가는 일을 방지 explicit vector2(double x_ = 0, double y_ = 0) : x(x_), y(y_) {} //두 벡터의 비교 bool operator == (const vector2& rhs) const { return x == rhs.x && y == rhs.y; } bool operator < (const vector2& rhs) const { return x != rhs.x ? x < rhs.x : y < rhs.y; } //벡터의 덧셈과 뺄셈 vector2 operator + (const vector2& rhs) const { return vector2(x + rhs.x, y + rhs.y); } vector2 operator -(const vector2& rhs) const { return vector2(x - rhs.x, y - rhs.y); } //실수로 곱셈 vector2 operator * (double rhs) const { return vector2(x*rhs, y*rhs); } //벡터의 길이를 반환 double norm() const { return hypot(x, y); } //방향이 같은 단위 벡터(unit vector)를 반환 //영벡터에 대해 호출한 경우 반환 값은 정의되지 않음 vector2 normalize() const { return vector2(x / norm(), y / norm()); } //x축의 양의 방향으로부터 이 벡터까지 반시계 방향으로 잰 각도 double polar() const { return fmod(atan2(y, x) + 2 * PI, 2 * PI); } //내적/외적의 구현 double dot(const vector2& rhs) const { return x* rhs.x + y*rhs.y; } double cross(const vector2& rhs) const { return x* rhs.y - rhs.x*y; } //이 벡터를 rhs에 사영한 결과 vector2 project(const vector2& rhs) const { vector2 r = rhs.normalize(); return r * r.dot(*this); } //점 p, a, b가 주어졌을 때, a가 b보다 p에 얼마나 더 가까운지 변환하는 함수 double howMuchCloser(vector2 p, vector2 a, vector2 b) { return (b - p).norm() - (a - p).norm(); } //점 p, a, b를 두개의 실수로 표현했을 때의 함수 double howMuchCloser(double px, double py, double ax, double ay, double bx, double by) { return sqrt((bx - px) * (bx - px) + (by - py) * (by - py)) - sqrt((ax - px) * (ax - px) + (ay - py) * (ay - py)); } };
//java code Vector addVector(Vector v){ this.x = v.x + x; this.y = v.y + y; return returnVec; //벡터의 합 } Vector subtractVector(Vector v){ returnVec.x = v.getX() - this.getX(); returnVec.y = v.getY() - this.getY(); return returnVec; //벡터의 차 } Vector scalaVector(double v){ //스칼라 : 방향은 고정되고 속도만 증가하는것. Vector returnVec= new Vector(0,0); this.x = this.getX()*v; this.y = this.getY()*v; return returnVec; //원점으로 시작점을 옮긴다 -(?) //백터의 뺄샘은 빼려는 것의 방향을 바꾸고 더하는것과같음. } void normalVector( ){ //x y z를 벡터의 크기로 나누는것. 정규화라고함. //두 벡터의 길이를 정규화해서 1미만으로 만들면 사이의 길이 값이 1이된다. double size = this.lengthVector(); this.x = x/size; this.y = y/size; } public double lengthVector(){ return Math.sqrt(Math.pow(x,2)+ Math.pow(y,2)); //벡터 사이의 길이를 구함 } public double innterProduct(Vector v){ //벡터의 내적 //계산 하고나면 스칼라값 즉 숫자 하나로 나온다. double result = this.x*v.getX() + this.y*v.getY(); double cos =0.0; cos = result/ (v.lengthVector()* this.lengthVector()); return cos; }

• 벡터의 내적
• 두 벡터 사이의 각도를 구하기 위해 사용
• http://j1w2k3.tistory.com/627

• 벡터의 사이각 구하기
  

• 벡터의 직각 여부 확인
• 벡터가 직각인 경우 가 항상 PI/2, 3PI/2일 수 밖에 없음.
• 내적이 0이라면 두 벡터는 항상 직각
  
• 벡터의 사영
• 벡터 b에 수직으로 빛이 내리쬘때 벡터 a가 b에 드리우는 그림자
• |a|cos
• 벡터의 외적
• 두 개의 3차원 벡터에 정의되는 연산, 3차원 벡터 a, b가 주어졌을 때 이 두 벡터에 모두 수직인 다른 벡터를 반환
• z축이 0인 경우 2차원 벡터에서도 정의가 가능
• http://blog.naver.com/PostView.nhn?blogId=mindo1103&logNo=90103361104&parentCategoryNo=&categoryNo=&viewDate=&isShowPopularPosts=false&from=postView

• 면적 계산하기

• 외적의 절대값은 a,b 두변으로 하는 평행 사변형의 넓이

• 원점과 a, b를 꼭지점으로 하는 삼각형 넓이의 정확히 두 배, 외적의 절대값을 구하고 이를 반으로 나누어 해당 삼각형의 넓이를 구하는 것이 가능

• 두 벡터의 방향 판별

• 각도가 단순히 두 벡터 사이각이 아니라 a에서 b까지 반시계 방향으로 얼마나 회전해야 하는 가를 나타냄
  
• sin- = -sin
• 외적이 양수인 경우

• 외적이 음수인 경우

//원점에서 벡터 b가 벡터ㅓ a의 반시계 방향이면 양수, 시계 방향이면 음수 //평행이면 0을 반환한다. double ccw(vector2 a, vector2 b) { return a.cross(b); } //점 p를 기준으로 벡터 b가 벡터 a의 반시계 방향이면 양수, 시계방향이면 음수, //평행이면 0을 반환 double ccw(vector2 p, vector2 a, vector2 b) { return ccw(a - p, b - p); }

교차와 거리, 면적

• 직선의 교차를 작성할 수 있는 간단한 방법
• a+pb 형태로 표현
• a+pb = c+qd 방정식 필요

(b와 d가 평행이면 교점의 정의되지 않음으로)

#define EPSILON 0.0001 // Define your own tolerance //(a,b)를 포함하는 선과 (c,d)를 포함하는 선의 교점을 x에 반환한다. //두 선이 평행이면 (겹치는 경우를 포함) 거짓을, 아니면 참을 반환한다. bool lineIntersection(vector2 a, vector2 b, vector2 c, vector2 d, vector2& x) { double det = (b - a).cross(d - c); //두 선이 평행인 경우 if (fabs(det) < EPSILON) return false; x = a + (b - a)*((c - a).cross(d - c) / det); return true; }

• 선분과 선분의 교차
• 직선과 직선의 교차를 계산 후 이 교차점이 두 선분 위에 오는지 확인
• 한 직선위의 두 선분의 관계
• 두 선분이 서로 겹치지 않음
• 두 선분이 한 점에서 닿음
  
• 두 선분이 겹쳐짐
• 한 선분이 다른 선분 안에 포함됨

/*선분과 선분의 교차점 계산 */ //(a,b)와 (c, d)가 평한 두 선분일 때 이들이 한 점에서 겹치는지 확인한다 bool parallelSegments(vector2 a, vector2 b, vector2 c, vector2 d, vector2& p) { if (b < a) std::swap(a, b); if (d < c) std::swap(c, d); //한 직선 위에 없거나 두 선분이 겹치지 않는 경우를 우선 걸러낸다. if (ccw(a, b, c) != 0 || b < c || d < a) return false; //두 선분은 확실히 겹친다. 교차점을 하나 찾자. if (a < c) p = c; else p = a; return true; } //p가 (a, b)를 감싸면서 각 변이 x, y축에 평행한 최소 사각형 내부에 있는지 확인한다. //a, b, p는 일직선 상에 있다고 가정한다. bool inBoundingRectangle(vector2 p, vector2 a, vector2 b) { if (b < a) std::swap(a, b); return p == a || p == b || (a < p && p < b); } //(a, b) 선분과 (c, d)선분의 교점을 p에 반환한다. //교점이 여러 개일 경우 아무 점이나 반환한다. //두 선분이 교차하지 않을 경우 false를 반환한다. bool segmentIntersection(vector2 a, vector2 b, vector2 c, vector2 d, vector2& p) { //두 직선이 평행인 경우를 우선 예외로 처리한다. if (!lineIntersection(a, b, c, d, p)) return parallelSegments(a, b, c, d, p); return inBoundingRectangle(p, a, b) && inBoundingRectangle(p, c, d); }

• 선분과 선분의 교차: 교차점이 필요 없을 때

• ccw(a, b, c) * ccw(a, b, d) < 0
• ccw(c, d, a) * ccw(c, d, a) < 0
• a<b, c<d인 네 점 a,b,c,d가 일직선 상에 있을 경우 두 선분이 겹치지 않게 하기 위해서는 b<c or d<a이어야 함.
• 교차하지 않고 접촉하는 경우도 인식하기 위해 ab or cd중 하나가 0인 경우에도 참을 반환

//두 선분이 서로 접촉하는지 여부를 반환한다. bool sementIntersects(vector2 a, vector2 b, vector2 c, vector2 d) { double ab = ccw(a, b, c)*ccw(a, b, d); double cd = ccw(c, d, a)*ccw(c, d, b); //두 선분이 한 직선에 위에 있거나 끝점이 겹치는 경우 if (ab == 0 && cd == 0) { if (b < a) std::swap(a, b); if (d < c) std::swap(c, d); return !(b< c || d <a); } return ab <= 0 && cd <= 0; }

• 점과 선 사이의 거리
• 한 직선에 포함된 두 점의 벡터 a,  b가 주어졌을 때 이 직선과 다른 점 p 사이의 거리
• 벡터의 내적을 이용하면 손쉽게 계산 가능
• 벡터 p-a가 벡터 b-a에 내린 벡터 사영을 이용해 계산 가능
• 수선의 발을 먼저 계산하고 나서 이와의 거리를 계산해 점과 선 사이의 거리를 계산

/*점과 선 사이의 거리를 계산하는 함수 pointToLine()*/ //점 p에서 (a,b) 직선에 내린 수선의 발을 구한다. vector2 perpendicularFoot(vector2 p, vector2 a, vector2 b) { return a + (p - a).project(b - a); } //점 p와 (a,b) 직선 사이의 거리를 구한다. double pointToLine(vector2 p, vector2 a, vector2 b) { return (p - perpendicularFoot(p, a, b)).norm(); }

기하 개념 code 종합

#include <iostream> #include <cstdlib> #include <cstdio> #include <algorithm> #include <math.h> #include <vector> #include <utility> #define EPSILON 0.0001 // Define your own tolerance const double PI = 2.0 * acos(0.0); //2차원 벡터를 표현 struct vector2 { double x, y; //생성자를 explicit으로 지정하면 vector2를 젛을 곳에 잘못해 //실수가 들어가는 일을 방지 explicit vector2(double x_ = 0, double y_ = 0) : x(x_), y(y_) {} //두 벡터의 비교 bool operator == (const vector2& rhs) const { return x == rhs.x && y == rhs.y; } bool operator < (const vector2& rhs) const { return x != rhs.x ? x < rhs.x : y < rhs.y; } //벡터의 덧셈과 뺄셈 vector2 operator + (const vector2& rhs) const { return vector2(x + rhs.x, y + rhs.y); } vector2 operator -(const vector2& rhs) const { return vector2(x - rhs.x, y - rhs.y); } //실수로 곱셈 vector2 operator * (double rhs) const { return vector2(x*rhs, y*rhs); } //벡터의 길이를 반환 double norm() const { return hypot(x, y); } //방향이 같은 단위 벡터(unit vector)를 반환 //영벡터에 대해 호출한 경우 반환 값은 정의되지 않음 vector2 normalize() const { return vector2(x / norm(), y / norm()); } //x축의 양의 방향으로부터 이 벡터까지 반시계 방향으로 잰 각도 double polar() const { return fmod(atan2(y, x) + 2 * PI, 2 * PI); } //내적/외적의 구현 double dot(const vector2& rhs) const { return x* rhs.x + y*rhs.y; } double cross(const vector2& rhs) const { return x* rhs.y - rhs.x*y; } //이 벡터를 rhs에 사영한 결과 vector2 project(const vector2& rhs) const { vector2 r = rhs.normalize(); return r * r.dot(*this); } /*점과 직선, 선분의 표현*/ //점 p, a, b가 주어졌을 때, a가 b보다 p에 얼마나 더 가까운지 변환하는 함수 double howMuchCloser(vector2 p, vector2 a, vector2 b) { return (b - p).norm() - (a - p).norm(); } //점 p, a, b를 두개의 실수로 표현했을 때의 함수 double howMuchCloser(double px, double py, double ax, double ay, double bx, double by) { return sqrt((bx - px) * (bx - px) + (by - py) * (by - py)) - sqrt((ax - px) * (ax - px) + (ay - py) * (ay - py)); } /* 두 벡터의 방향성을 판단하는 ccw 함수*/ //원점에서 벡터 b가 벡터 a의 반시계 방향이면 양수, 시계 방향이면 음수 //평행이면 0을 반환한다. double ccw(vector2 a, vector2 b) { return a.cross(b); } //점 p를 기준으로 벡터 b가 벡터 a의 반시계 방향이면 양수, 시계방향이면 음수, //평행이면 0을 반환 double ccw(vector2 p, vector2 a, vector2 b) { return ccw(a - p, b - p); } /*두 직선의 교차점을 계산*/ //(a,b)를 포함하는 선과 (c,d)를 포함하는 선의 교점을 x에 반환한다. //두 선이 평행이면 (겹치는 경우를 포함) 거짓을, 아니면 참을 반환한다. bool lineIntersection(vector2 a, vector2 b, vector2 c, vector2 d, vector2& x) { double det = (b - a).cross(d - c); //두 선이 평행인 경우 if (fabs(det) < EPSILON) return false; x = a + (b - a)*((c - a).cross(d - c) / det); return true; } /*선분과 선분의 교차점 계산 */ //(a,b)와 (c, d)가 평한 두 선분일 때 이들이 한 점에서 겹치는지 확인한다 bool parallelSegments(vector2 a, vector2 b, vector2 c, vector2 d, vector2& p) { if (b < a) std::swap(a, b); if (d < c) std::swap(c, d); //한 직선 위에 없거나 두 선분이 겹치지 않는 경우를 우선 걸러낸다. if (ccw(a, b, c) != 0 || b < c || d < a) return false; //두 선분은 확실히 겹친다. 교차점을 하나 찾자. if (a < c) p = c; else p = a; return true; } //p가 (a, b)를 감싸면서 각 변이 x, y축에 평행한 최소 사각형 내부에 있는지 확인한다. //a, b, p는 일직선 상에 있다고 가정한다. bool inBoundingRectangle(vector2 p, vector2 a, vector2 b) { if (b < a) std::swap(a, b); return p == a || p == b || (a < p && p < b); } //(a, b) 선분과 (c, d)선분의 교점을 p에 반환한다. //교점이 여러 개일 경우 아무 점이나 반환한다. //두 선분이 교차하지 않을 경우 false를 반환한다. bool segmentIntersection(vector2 a, vector2 b, vector2 c, vector2 d, vector2& p) { //두 직선이 평행인 경우를 우선 예외로 처리한다. if (!lineIntersection(a, b, c, d, p)) return parallelSegments(a, b, c, d, p); return inBoundingRectangle(p, a, b) && inBoundingRectangle(p, c, d); } //두 선분이 서로 접촉하는지 여부를 반환한다. bool sementIntersects(vector2 a, vector2 b, vector2 c, vector2 d) { double ab = ccw(a, b, c)*ccw(a, b, d); double cd = ccw(c, d, a)*ccw(c, d, b); //두 선분이 한 직선에 위에 있거나 끝점이 겹치는 경우 if (ab == 0 && cd == 0) { if (b < a) std::swap(a, b); if (d < c) std::swap(c, d); return !(b< c || d <a); } return ab <= 0 && cd <= 0; } /*점과 선 사이의 거리를 계산하는 함수 pointToLine()*/ //점 p에서 (a,b) 직선에 내린 수선의 발을 구한다. vector2 perpendicularFoot(vector2 p, vector2 a, vector2 b) { return a + (p - a).project(b - a); } //점 p와 (a,b) 직선 사이의 거리를 구한다. double pointToLine(vector2 p, vector2 a, vector2 b) { return (p - perpendicularFoot(p, a, b)).norm(); } };

문제: 핀볼 시뮬레이션(난이도 상)

• 문제: https://algospot.com/judge/problem/read/PINBALL

동혁이는 새로 나온 핀볼 게임을 하고 있습니다. 이 핀볼 게임은 아주 큰 게임판에 여러 개의 장애물을 놓고, 밖에서 공을 던져 가장 많은 장애물을 맞추는 것을 목표로 합니다. 공과 장애물은 완전한 원형이며, 공은 장애물에 부딪히면 항상 입사각과 반사각이 같도록 정반사됩니다. 장애물은 게임판에 고정되어 있기 때문에 움직이지 않습니다. 위 그림은 예제 게임판에서 동혁이가 던진 공이 장애물들에 순서대로 부딪히는 과정을 보여줍니다.

장애물들의 위치와 크기가 주어지고, 공의 위치와 방향이 주어질 때 공이 부딪히는 장애물들의 번호를 계산하는 프로그램을 작성하세요. 공은 반지름 1인 원이며, 항상 직선으로 움직인다고 가정합니다.

입력

입력의 첫 줄에는 테스트 케이스의 수 C(C <= 50)가 주어집니다. 각 테스트 케이스의 첫 줄에는 5개의 정수 x, y, dx, dy, N 이 주어집니다. 이 때 x, y (0 <= x, y <= 100)는 현재 공 중심의 위치를 나타내며, dx, dy (-10 <= dx, dy <= 10)는 공의 방향을 나타냅니다. 이 때 공은 현재 (x,y) 에서 (x+dx, y+dy)를 향해 일직선으로 굴러가고 있습니다. N(0 <= N <= 50)은 장애물의 개수를 나타내며, 각 장애물에는 0번부터 N-1번까지의 번호가 매겨져 있습니다. 그 다음 N줄에는 0번부터 각 장애물의 정보가 3개의 정수 x_i, y_i, r_i 로 주어집니다. (0 <= x_i, y_i, r_i <= 100) 이 때 i번 장애물은 (x_i, y_i)을 중심으로 한 반지름 r_i인 원 형태를 하고 있습니다.

장애물은 게임판 밖에 걸쳐 있을 수도 있습니다. 두 장애물 사이의 거리가 2 이하인 경우는 없습니다. 시작 위치에서 공은 장애물과 겹쳐 있거나 닿아 있지 않습니다.

공은 게임판의 벽에 닿는 순간 게임판에서 떨어진다고 가정합니다. 입력에 공이 게임판의 벽과 장애물에 동시에 닿는 경우는 주어지지 않습니다.

출력

각 테스트 케이스마다 한 줄에 공이 닿는 순서대로 각 장애물의 번호를 출력합니다. 만약 공이 장애물에 100번 넘게 부딪힐 경우 첫 100번 부딪히는 장애물의 번호만을 출력합니다.

주어진 입력에서 장애물의 크기나 위치가 10^-7 만큼 변하더라도 답은 변하지 않는다고 가정해도 좋습니다.

• 조건 확인

• x,y : 공의 중심 위치

• dx, dy : 공의 방향

• N : 장애물 개수

• 장애물은 게임판 밖에 존재할 수도 있음

• 두 장애물 사이의 거리는 2보다 작지 않음

• 부딪히는 장애물의 번호 출력

• 공이 100번 넘게 부딪히는 경우 장애물만 출력

• 풀이

• 공은 반지름이 1인 원 => 점으로 바꾸고 장애물들의 반지름을 1씩 증가

• 공의 위치와 방향을 알아야 함

• 공의 위치를 나타내는 위치 벡터

• 공의 방향을 나타내는 단위 벡터

• 시뮬레이션

• 출력후 장애물 번호 출력

• 방향과 위치 변경을 매 부딪히는 순간마다 해주어야 함

• 반복 실행

• 공과 장애물의 충돌 시점

• 위치벡터 : p, 방향 벡터: q

• 1초에 d만큼 굴러간다고 가정

• f(t) = p+t·d(t>=0)

• 중심이 c, 반지름이 r인 장애물

• |c-f(t)| = r

• |a|^2 = a·a

• (c-f(t))·(c-f(t)) = r^2

• d·dt^2 +2(p-c)·d·t + c·c + p·p -2(c·p) -r^2 =0

• 코드

#include<cstdio> #include<cmath> #include<iostream> #include<string> #include<vector> using namespace std; int n; // 각 원 중심의 위치와 반지름 vector2 center[50]; int radius[50]; const double EPSILON = 1e-8; const double INFTY = 1e200; struct vector2 { double x, y; //생성자를 explicit으로 지정하면 vector2를 젛을 곳에 잘못해 //실수가 들어가는 일을 방지 explicit vector2(double x_ = 0, double y_ = 0) : x(x_), y(y_) {} //두 벡터의 비교 bool operator == (const vector2& rhs) const { return x == rhs.x && y == rhs.y; } bool operator < (const vector2& rhs) const { return x != rhs.x ? x < rhs.x : y < rhs.y; } //벡터의 덧셈과 뺄셈 vector2 operator + (const vector2& rhs) const { return vector2(x + rhs.x, y + rhs.y); } vector2 operator -(const vector2& rhs) const { return vector2(x - rhs.x, y - rhs.y); } //실수로 곱셈 vector2 operator * (double rhs) const { return vector2(x*rhs, y*rhs); } //벡터의 길이를 반환 double norm() const { return hypot(x, y); } //방향이 같은 단위 벡터(unit vector)를 반환 //영벡터에 대해 호출한 경우 반환 값은 정의되지 않음 vector2 normalize() const { return vector2(x / norm(), y / norm()); } //x축의 양의 방향으로부터 이 벡터까지 반시계 방향으로 잰 각도 double polar() const { return fmod(atan2(y, x) + 2 * PI, 2 * PI); } //내적/외적의 구현 double dot(const vector2& rhs) const { return x* rhs.x + y*rhs.y; } double cross(const vector2& rhs) const { return x* rhs.y - rhs.x*y; } //이 벡터를 rhs에 사영한 결과 vector2 project(const vector2& rhs) const { vector2 r = rhs.normalize(); return r * r.dot(*this); } /*점과 직선, 선분의 표현*/ //점 p, a, b가 주어졌을 때, a가 b보다 p에 얼마나 더 가까운지 변환하는 함수 double howMuchCloser(vector2 p, vector2 a, vector2 b) { return (b - p).norm() - (a - p).norm(); } //점 p, a, b를 두개의 실수로 표현했을 때의 함수 double howMuchCloser(double px, double py, double ax, double ay, double bx, double by) { return sqrt((bx - px) * (bx - px) + (by - py) * (by - py)) - sqrt((ax - px) * (ax - px) + (ay - py) * (ay - py)); } /* 두 벡터의 방향성을 판단하는 ccw 함수*/ //원점에서 벡터 b가 벡터 a의 반시계 방향이면 양수, 시계 방향이면 음수 //평행이면 0을 반환한다. double ccw(vector2 a, vector2 b) { return a.cross(b); } //점 p를 기준으로 벡터 b가 벡터 a의 반시계 방향이면 양수, 시계방향이면 음수, //평행이면 0을 반환 double ccw(vector2 p, vector2 a, vector2 b) { return ccw(a - p, b - p); } /*두 직선의 교차점을 계산*/ //(a,b)를 포함하는 선과 (c,d)를 포함하는 선의 교점을 x에 반환한다. //두 선이 평행이면 (겹치는 경우를 포함) 거짓을, 아니면 참을 반환한다. bool lineIntersection(vector2 a, vector2 b, vector2 c, vector2 d, vector2& x) { double det = (b - a).cross(d - c); //두 선이 평행인 경우 if (fabs(det) < EPSILON) return false; x = a + (b - a)*((c - a).cross(d - c) / det); return true; } /*선분과 선분의 교차점 계산 */ //(a,b)와 (c, d)가 평한 두 선분일 때 이들이 한 점에서 겹치는지 확인한다 bool parallelSegments(vector2 a, vector2 b, vector2 c, vector2 d, vector2& p) { if (b < a) std::swap(a, b); if (d < c) std::swap(c, d); //한 직선 위에 없거나 두 선분이 겹치지 않는 경우를 우선 걸러낸다. if (ccw(a, b, c) != 0 || b < c || d < a) return false; //두 선분은 확실히 겹친다. 교차점을 하나 찾자. if (a < c) p = c; else p = a; return true; } //p가 (a, b)를 감싸면서 각 변이 x, y축에 평행한 최소 사각형 내부에 있는지 확인한다. //a, b, p는 일직선 상에 있다고 가정한다. bool inBoundingRectangle(vector2 p, vector2 a, vector2 b) { if (b < a) std::swap(a, b); return p == a || p == b || (a < p && p < b); } //(a, b) 선분과 (c, d)선분의 교점을 p에 반환한다. //교점이 여러 개일 경우 아무 점이나 반환한다. //두 선분이 교차하지 않을 경우 false를 반환한다. bool segmentIntersection(vector2 a, vector2 b, vector2 c, vector2 d, vector2& p) { //두 직선이 평행인 경우를 우선 예외로 처리한다. if (!lineIntersection(a, b, c, d, p)) return parallelSegments(a, b, c, d, p); return inBoundingRectangle(p, a, b) && inBoundingRectangle(p, c, d); } //두 선분이 서로 접촉하는지 여부를 반환한다. bool sementIntersects(vector2 a, vector2 b, vector2 c, vector2 d) { double ab = ccw(a, b, c)*ccw(a, b, d); double cd = ccw(c, d, a)*ccw(c, d, b); //두 선분이 한 직선에 위에 있거나 끝점이 겹치는 경우 if (ab == 0 && cd == 0) { if (b < a) std::swap(a, b); if (d < c) std::swap(c, d); return !(b< c || d <a); } return ab <= 0 && cd <= 0; } /*점과 선 사이의 거리를 계산하는 함수 pointToLine()*/ //점 p에서 (a,b) 직선에 내린 수선의 발을 구한다. vector2 perpendicularFoot(vector2 p, vector2 a, vector2 b) { return a + (p - a).project(b - a); } //점 p와 (a,b) 직선 사이의 거리를 구한다. double pointToLine(vector2 p, vector2 a, vector2 b) { return (p - perpendicularFoot(p, a, b)).norm(); } }; /* 장애물과의 충돌 시점을 찾는 함수 */ // 2차방정식 ax^2 + bx + c = 0 의 모든 실근을 크기 순서대로 반환한다 vector<double> solve2(double a, double b, double c) { double d = b*b - 4 * a*c; //해가 없는 경우 if (d < -EPSILON) return vector<double>(); //중근 if (d < EPSILON) return vector<double>(1, -b / (2 * a)); vector<double> ret; ret.push_back((-b - sqrt(d)) / (2 * a)); ret.push_back((-b + sqrt(d)) / (2 * a)); return ret; } // here 에 있던 공이 1초마다 dir 만큼 굴러갈 때, center 를 중심으로 하고 // 반지름 radius 인 장애물과 몇 초 후에 충돌하는지 반환한다 // 충돌하지 않을 경우 '아주 큰 값' INFTY를 반환한다. double hitCircle(vector2 here, vector2 dir, vector2 center, int radius) { double a = dir.dot(dir); double b = 2 * dir.dot(here - center); double c = center.dot(center) + here.dot(here) - 2 * here.dot(center) - radius * radius; vector<double> sols = solve2(a, b, c); if (sols.empty() || sols[0] < EPSILON) return INFTY; return sols[0]; } // here 에 있던 공이 dir 방향으로 굴러와 center 를 중심으로 하는 장애물에서 // contact 위치에서 충돌했을 때 공의 새로운 방향을 반환한다 // 시작 벡터 -d, 중심점에 대해 사영, 사영까지 변화량 확인, 다시 반대방향으로 추가 d 구하기 vector2 reflect(vector2 here, vector2 dir, vector2 center, vector2 contact) { return (dir - dir.project(contact - center) * 2).normalize(); } // 공의 현재 위치와 방향이 주어질 때 최대 100번의 충돌을 출력한다 // 시간 복잡도 O(nk) void simulate(vector2 here, vector2 dir) { // 방향은 항상 단위 벡터로 저장하자 dir = dir.normalize(); int hitCount = 0; while (hitCount < 100) { // fprintf(stderr, "here (%.2lf, %.2lf)\n", here.x, here.y); // 이번에 충돌할 장애물을 찾는다. int circle = -1; double time = INFTY*0.5; // 각 장애물을 순회하며 가장 먼저 만나는 장애물을 찾는다(반지름 1 더해주는 것 잊으면 안된다) for (int i = 0; i < n; ++i) { double cand = hitCircle(here, dir, center[i], radius[i] + 1); if (cand < time) { time = cand; circle = i; } } // 더 이상 장애물에 충돌하지 않고 게임판을 벗어난 경우 if (circle == -1) break; // 부딪히는 장애물의 번호를 출력한다 if (hitCount++) cout << " "; cout << circle; // 공의 새 위치를 계산한다 vector2 contact = here + dir * time; // 부딪힌 위치와 새 방향으로 here 와 dir 을 변경한다 dir = reflect(here, dir, center[circle], contact); here = contact; } cout << endl; } int main() { int cases; cin >> cases; for (int cc = 0; cc < cases; ++cc) { int x, y, dx, dy; cin >> x >> y >> dx >> dy >> n; vector2 here, dir; here = vector2(x, y); dir = vector2(dx, dy).normalize(); for (int i = 0; i < n; i++) { int x, y, r; cin >> x >> y >> radius[i]; center[i] = vector2(x, y); } simulate(here, dir); } }