원뿔을 평면으로 잘랐을 때 나타나는 곡선인 conic section은 직교 좌표계에서 (x,y)(x,y)에 대한 2차 형식으로 쓰인다.

F(x,y)=ax2+bxy+cy2+dx+ey+f=0F(x,y)=ax2+bxy+cy2+dx+ey+f=0

이 conic section이 타원을 기술할 때 parameter {a,b,c,d,e,fa,b,c,d,e,f}를 이용해서 타원의 중심, 장축과 단축의 길이, 그리고 회전각에 대한 공식을 구해보자. 2차 항을 행렬을 써서 표현하면

(x,y)(ab/2b/2c)(xy)+=0

따라서, 적절한 회전 변환을 하면 두 좌표의 곱으로 주어지는 xy-항을 없앨 수 있다. 회전 변환을 시행하면 행렬은 eigenvalue를 성분으로 하는 대각 행렬이 된다. 회전 변환이 determinant를 보존하므로 determinant는 행렬의 두 eigenvalue의 곱으로 주어짐을 알 수 있다.

(ab/2b/2c)(λ100λ2)det=acb2/4

회전 후의 방정식이 타원의 방정식(원점이 이동된)을 기술하기 위해서는 det>0 이어야 한다 (회전시킨 후 식에서 x2y2의 계수는 두 eigenvalue로 주어지므로 같은 부호를 가져야 한다.)

F=ellipseb24ac<0

conic section F가 타원을 기술한다면, 다음과 같이 평행이동을 시켜서 타원의 중심 (x0,y0)이 원점에 놓이게 하고, 회전변환을 시켜서 xy 항을 없애도록 하자.

xx0+xcosθysinθ,yy0+xsinθ+ycosθ

여기서 회전각 θx을 기준으로 측정된다. F에 적용해서 xy-항이 없어지는 조건과 1차 항이 사라지도록 하는 조건을 찾으면 타원의 중심 (x0,y0)θ는 

tan2θ=bacellipse center:  (x0,y0)=(2cdbeb24ac,2aebdb24ac)

로 주어짐을 확인할 수 있다. Eigenvalue의 제곱근의 역수가 두 축의 반지름을 결정하므로 음수가 되어서는 안된다 (둘 중 하나가 0이면 직선이고, 둘 모두 0이면 한 점에 해당). 이는 대칭행렬의 고유값이 항상 0보다 작지 않다는 사실에서 기인한다. 위에서 구한 회전각 θ를 이용해서 두 고유값을 표현하면 ,

λ1=acos2θ+bcosθsinθ+csin2θ

λ2=asin2θbcosθsinθ+ccos2θ

위의 회전변환식을 대입했을 떄 나머지 상수항은 (첫 번째 등호는 계산해서 확인할 수 있음)

ax20+bx0y20+cy20+dx0+ey0+f=f(ax20+bx0y0+cy20)scale1

로 주어짐을 알 수 있다. 따라서 회전시킨 타원은 표준형 꼴

λ1x2+λ2y2=scale1

로 표현된다. 이 표준형 타원의 두 축의 반지름은 각각

rx=scale1λ1,ry=scale1λ2 로 주어진다.

// conic_params(a, b, c, d, e, f): ax^2 + bxy + cy^2 + dx + ey + f = 0;
// ellipse_params(radius_x, radius_y, center_x, center_y, tilt_angle w.r.t x-axis);
bool conic_to_ellipse(double conic_params[6], double ellipse_params[5]) {
    const double a = conic_params[0];
    const double b = conic_params[1];
    const double c = conic_params[2];
    const double d = conic_params[3];
    const double e = conic_params[4];
    const double f = conic_params[5];
    // get ellipse orientation w.r.t x-axis;
    const double theta = 0.5 * atan2(b, a - c);
    // get scaled x/y radius;
    const double ct = cos(theta);
    const double st = sin(theta);
    const double ap = a * ct * ct + b * ct * st + c * st * st;
    const double cp = a * st * st - b * ct * st + c * ct * ct;
    // get center of ellipse;
    const double cx = (2 * c * d - b * e) / (b * b - 4 * a * c);
    const double cy = (2 * a * e - b * d) / (b * b - 4 * a * c);
    // get scale factor
    const double val = a * cx * cx + b * cx * cy + c * cy * cy;
    const double scale_inv = val - f;
    if (scale_inv / ap <= 0 || scale_inv / cp <= 0) {
        TRACE("Error! ellipse parameters are imaginary\n");
        return 0;
    }
    ellipse_params[0] = sqrt(scale_inv / ap);
    ellipse_params[1] = sqrt(scale_inv / cp);
    ellipse_params[2] = cx;
    ellipse_params[3] = cy;
    ellipse_params[4] = theta;
    return 1;
};
728x90
,