천문 프로그래밍과 관련된 이야기를 나눌 수 있는 곳입니다.
총 게시물 144건, 최근 0 건
   
지구의 세차운동(歲差運動, Precession) 계산하기
글쓴이 : 지용호 별님  (121.♡.140.207) 날짜 : 2009-04-14 (화) 13:13 조회 : 14417

세차운동(歲差運動, Precession)은 일반적으로 강체에 돌림힘(Torque)이 작용할 때, 회전하는 물체의 축이 어떤 부동축의 둘레를 회전하는 현상을 말한다. 천문학에서 말하는 세차운동은 지구의 자전축이 불변인 황도면의 축의 둘레를 2만 5800년 주기로 회전하는 운동을 말한다. 이 현상은 지구에 작용하는 달, 태양의 중력과 지구가 타원체라는 것에서부터 기인한다.

 

천문 계산을 하는데 있어서 세차를 고려해야 과거,미래의 천체위치를 제대로 예측할 수 있다. 왜냐하면 세차운동으로 인해 지구 자전축이 변한다는 것은 천체 위치의 기준점인 춘분점(Vernal Equinox)의 위치도 변한다는 것을 의미하기 때문이다.

 

여기서는 적도좌표계(Equatorial Coordinate System)황도좌표계(Ecliptic Coordinate System)에서 세차운동 계산 방법을 학습하는 것을 주 목표로 한다. 세차운동의 공식 유도 및 사용되는 각종 다항식이 어떻게 만들어졌는지 알아내는 과정은 생략한다.

 

이 글의 참고 서적은 Jean Meeus의 Astronomical Algorithms 2nd임을 밝혀둔다. 하지만 이 책에서 소개한 대로만 하지 않으며 더 효과적으로 계산하는 방법도 알아본다.

 

이 글은 세차운동 계산시 정밀한 방법(Rigorous Method)을 이용한다. 낮은 정밀도를 가지는 방법은 책을 참고하길 바란다.

 

Julian Century

1 Julian Century는 100 Julian Year와 같다.  Julian Year는 86000 SI 초를 기반으로하는 365.25일을 1로 측정하는 단위이다. Julian Century는 Julian Year의 100배 이므로 36525일을 1로 측정하는 단위이다.  보통 J2000.0 = JDE 2451545.0 = 2000년 1월 1일 12h DT를 0으로 시작하게 된다. 여기서 JDE는 Julian Ephemeris Date, DT는 Dynamic Time을 의미한다.

T = (JD0 – J2000) / 36525

세차운동의 경우 어떤 시작원기(the initial epoch)부터 마지막 원기(the final epoch)까지의 세차운동값을 계산해야한다. 그러므로 시작원기와 마지막원기의 Julian Century차이는 다음과 같다.

t = (JD – JD0) / 36525

T, t에서 JD0는 시작원기, JD는 마지막 원기를 의미한다.

 

적도좌표계에서 세차운동적용 (일반 공식 이용)

Meeus의 책에 보면 다음과 같은 적도좌표계에서 세차운동을 계산하기 위한 수치적 표현공식이 적혀있다.

ζ = (2306”.2181 + 1”.39656T – 0”.000139 T2 ) t
      + (0”.30188 – 0”.000344T)t2 + 0”.017998 t3

z = ζ + (0”.79280+0”.000411*T)t2 + 0”.000205t3

ϑ =(2004”.3109 – 0”.85330T + 0”.000217T2)t
       - ( 0”.42665+0”.000217T )t2 + 0.041833t3

이 식은 IAU(International Astronomical Union, 1976)에서 지정한 값으로 관측을 통해 얻어지는 데이타로 만들어지는 식이다.

시작원기의 적경,적위를 (α0, δ0)라고 하고 세차운동이 적용된 마지막원기의 적경,적위를 (α, δ)로 한다면 위 수치적 표현공식을 이용해 다음과 같이 (α, δ)를 계산할 수 있다. 이 식을 유도하는 과정은 구면천문학에 관련된 책을 참고하기 바란다.

A = cos(δ0)  sin( α0 + ζ )

B = cos( ϑ ) cos( δ0 ) cos( α0 + ζ ) – sin( ϑ ) sin( δ0 )

C = sin( ϑ ) cos ( δ0 ) cos( α0 + ζ ) + cos( ϑ ) sin( δ0 )

tan(  α – z ) = A/B

sin( δ ) = C

 

이 식만으로도 정확한 계산을 할 수 있다. 그러나 보통 컴퓨터 프로그램을 이용하여 계산을 한다면 이 식은 매우 비효율적인 식이 된다. 가령, 수천, 수만개의 별이 있다고 하자. 그 별들은 어떤 시작원기를 가지는 각각 다른 α0, δ0값을 가지고 있을 것이다. 이들 α0, δ0을 세차운동이 적용된 새로운 α, δ를 계산하기 위해 일일히 cos, sin을 계산한다. 컴퓨터 프로그램을 한 사람은 알겠지만 cos, sin은 컴퓨터의 많은 자원을 소비하는 급수형태의 계산식이다. cos, sin을 너무 많이 쓰면 계산의 효율이 급격하게 떨어져서 1분이면 도출할 수 있는 계산을 1시간~2시간 이상식 소비할 수도 있게 된다.

 

컴퓨터 프로그램을 이용한 계산을 하려면 되도록이면 cos, sin과 같이 매우 비싼 자원을 소비하는 함수 사용을 지양하도록 해야한다. 그래야 빠른 계산을 할 수 있기 때문이다.

 

적도좌표계에서 세차운동적용 ( 변환행렬이용 )

앞서 세차운동이 적용된 새로운 α, δ를 계산하는 식은 한개의 천체 위치를 계산할 때마다 α0, δ0이 포함된 cos, sin을 계산해야한다. 그래서 새로운 α0, δ0가 주어져도 cos, sin을 최소한으로 이용하여 계산하는 방법을 여기서 소개한다. 그것은 변환행렬(Transformation Matrix)를 이용하는 방법이다.

 

변환행렬을 이용하면 다음과 같은 계산이 가능해진다. 세차운동을 적용하는 3×3 변환행렬을 P라고 하자. 그리고 α0, δ0를 직교좌표계로 바꾼 백터(vector)를 r0=[x0,y0,z0]라고 하고 세차운동이 적용된 α, δ를 직교좌표계를 바꾼 백터를 r=[x,y,z]라고 한다. 그럼 다음과 같이 계산하면 되겠다.

DRW000051a879c6

변환행렬 P는 세차운동이 적용되어 앞서 설명한 공식과 다르게 α0, δ0이 분리가 되었다. 행렬 P만 구하면 cos, sin과 같이 큰 자원을 소비하는 함수를 계산할 때마다 이용할 필요가 없어지며 대신, 단순하게 행렬과 백터 곱으로 수천,수만개의 별을 빠른 시간안에 세차운동을 적용할 수 있게 된다.

 

그럼 P는 어떻게 구할까? 그것은 회전변환행렬을 이용하면 된다. 앞서 구한 ζ, z, ϑ 값은 모두 세차운동에 의해 계산된 각도이다. 이 각도를 이용해 시작원기때의 적도좌표를 3번 회전변환하여 마지막원기때의 적도좌표를 구하는 것이다.

회전변환행렬에서 X축 회전 행렬을 Rx, Y축 회전 행렬을 Ry, Z축 회전행렬을 Rz라고 하자. 이들 행렬은 모두 3×3 행렬이다.

DRW000051a879be

 

DRW000051a879bc

 

DRW000051a879ba

다른 참고를 보면 행렬안에 sin의 부호가 바뀐것을 볼 수 있는데 회전의 대상을 어떻게 정하느냐에 따라서 다른 결과이다. 천문계산에서 회전행렬은 위 정의대로 하는 것이 좋다.

 

적도좌표계에서 세차운동을 적용한 변환행렬 P는 다음과 같이 주어진다.

 

DRW000051a879b8

 

위 변환행렬은 ζ, z, ϑ 이 적용된 행렬곱이다. 이 행렬곱의 결과는 다음과 같다.

DRW000051a879c4

 

컴퓨터 프로그램으로 계산할 때는 이렇게 변환행렬을 만들어 계산하면 되겠다.

 

황도좌표계에서 세차운동적용 ( 변환행렬이용 )

 

황도좌표계에서도 적도좌표계와 동일한 방식이 적용된다.

Π = 174˚.876383889 + 3289”.4789 + 0”.60622 T2  
      - (869”.8089 + 0”.50491T)t + 0”.03536 t2

π = (5029”.0966 + 2”.22226T – 0”.000042T2)t 
     + (1”.11113 – 0”.000042T)t2 – 0”.000006t3

η = (47”.0029 – 0”.06603T – 0”.000598T2)t
     + (-0”.03302 + 0”.000598T)t2 + 0.000060t3

 

시작원기 황경(ecliptic longitude),황위(ecliptic latitude)가 λ0,β0 이고 마지막원기 λ,β이라고 하자. 이때 각각의 직교좌표값을 r0, r이라고 한다면 변환행렬 P는 다음과 같다.

 

DRW000051a879b6

DRW000051a879c6

 

결과적으로 P는 다음과 같이 결정된다.

DRW000051a879c2

 

응용하기

지금까지 세차운동이 적용된 2가지 변환행렬을 알아보았다. 이 행렬은 필요에 따라서 선택해서 사용한다. 한가지 예를 들어보겠다. 적도좌표계에서 세차운동을 적용한 변환행렬을 Pequ, 황도좌표계에서 세차운동을 적용은 변환행렬을 Pecl이라고 하자. 황도좌표계에서 적도좌표계로 변환하는 행렬을 EclToEqu라고 가정한다. 만약 어떤 행성의 황경과 황위를 계산했다고 하면 다음과 같은 2가지 방법으로 새로운 적경,적위를 계산할 수 있을 것이다.

r = EclToEqu Pecl  r0

r = Pequ Ecl2Equ  r0

여기서 r0는 황경,황위를 담은 직각좌표계에서의 벡터값, r은 적경,적위를 담은 직각좌표계에서의 벡터값

 

세차운동 계산기(Precession Calculator)

이 프로그램은 지금까지 설명한 내용을 토대로 만든 프로그램이다.  

 

[세차운동 계산기 실행]

 

참고

글쓴이 : 지돌스타(http://jidolstar.com/blog/archives/1198)

 


송용준별님 (163.♡.175.176) 2009-04-14 (화) 16:16
참고 입니다.
2006년부터 IAU에서 지정한 표준 Precession Formla는 IAU2000A 모델입니다. 실제 Astronomical Almanac에서도 2006년부터 이 모델을 채택하고 있는데요. 위의 Precession Formula의 다항식이 다음과 같이 바뀌었습니다.
zetaA = +2''.5976176 + 2306''.0809506T + 0''.3019015T^2 + 0''.0179663T^3 -32''.7*10^-6T^4 - 0''.2 * 10^-6 T^5

zA = -2''.5976176 + 2306''.0803226T +1''.0947790T^2 + 0''.0182273T^3 + 47''.0 * 10^-6 T^4 -0''.3 * 10^-6 *T^5

thetaA = +2004''.1917476T - 0''.4269353T^2 - 0''.0418251T^3 - 60''. * 10^-6 T^4 - 0''.1 * 10^-6 T^5

T = (JD - 245545.0)/365.25
''는 arcsec 단위입니다.
댓글주소
송용준별님 (163.♡.175.176) 2009-04-14 (화) 16:16
이것은 적도 좌표계에서의 Precession Formulae이고, 황도 좌표계에서의 Precession Formulae 도 역시 변경되었으니 자세한것은 Astronomical Almanac을 참고하시기 바랍니다.
댓글주소
지용호별님 (121.♡.140.207) 2009-04-14 (화) 17:17
송용준//
좋은 정보 감사합니다.
IAU2000A 에 대한 것은 참고자료로 링크가 있습니다.
댓글주소
이형철별님 (203.♡.137.5) 2009-04-14 (화) 23:23
좋은 글이네요. 도움이 많이 되었습니다.
Astroview에서도 IAU2000A를 사용하고 있는 것으로 기억합니다.
세계 표준이 정해지기 전에는 다양한 공식들이 존재하다가
표준이 정해지면 game over 인 듯한 느낌이 듭니다.
댓글주소
송용준별님 (163.♡.175.176) 2009-04-15 (수) 11:11
헉.. T =(JD-2451545.0)/365.25 입니다. 중간에 1이 빠졌네요.
댓글주소
   

총 게시물 144건, 최근 0 건
번호 제목 글쓴이 날짜 조회 추천
144  음양력 변환 프로그램 CalTime 3.4 +1 김동빈별님 11-24 7091 0
143  현재 시간에 따른 태양, 달 위치 계산. +1 조재훈별님 08-30 6340 0
142  [DreamSpark] Microsoft의 정품 프로그램을 무료로 … +6 백승우별님 01-21 8224 0
141  일식 예측을 위한 태양과 달의 위치 계산 +1 이형철별님 11-21 10750 0
140  사이토구니치의 古天文學 번역판을 공개하며 +1 이형철별님 08-22 11151 2
139  [AstroDev의 C언어 강좌] 나왔다 Hello World!!! +3 유환용별님 07-01 8574 0
138  [AstroDev의 C언어 강좌] 컴퓨터 데이터 표현 방… 유환용별님 06-28 9844 0
137  IAU2000 장동 모델 +2 김창환별님 06-25 11766 0
136  [AstroDev의 C언어 강좌]메모리 안에서는 무슨일… +10 유환용별님 06-10 11921 0
135  [AstroDev의 C언어 강좌]우리는 무엇을 준비해야… +7 유환용별님 06-08 9816 0
134  [AstroDev의 C언어 강좌]우리는 무엇을 준비해야… +3 유환용별님 06-07 9392 0
133  [AstroDev의 C언어 강좌]C언어를 공부하기 전에 … +4 유환용별님 06-06 9981 0
132  IAU2006 세차 모델 +1 김창환별님 05-18 10559 0
131  각도변환문제 : 도(degree), 시(hour), 라디안(radia 지용호별님 04-15 17930 0
130  지구의 세차운동(歲差運動, Precession) 계산하… +5 지용호별님 04-14 14418 0
 1  2  3  4  5  6  7  8  9  10  맨끝
 
Since 2001.2.7 미래창조과학부 등록 비영리민간단체 천문노트. Copyright All rights reserved.
단체명 : 천문노트  |    고유번호 : 101-82-15888  |    대표자명 : 김태욱, 조우성  |    주소 : 138-804 서울특별시 송파구 가락동 93 금강빌딩 7층 710호  |    전화 : 02-543-3295  |    Fax : 02-6918-6888  |    통신판매신고번호 : 종로 제01-5696호  |    개인정보관리책임자 및 사이트관리자 : 지용호