2013년 7월 19일 금요일

정확한 24절기 시각을 계산하는 방법


정확한 24절기 입기시각을 계산하는 방법


정확한 24절기를 계산하기 위해서는 다음 방식을 사용한다.

우선 계산하는데 필요한 관련식은 다음에 소개하는 (1),(2),(4)를 사용한다.

태양의 평균경도(평균황경)   L= 279.827287 + 0.98564736*d      (1)
태양의 근점이각                  M = 356.666444 + 0.98560028*d     (2)
d는 2013년 1월 0일 0시(한국시각 1월 0일 9시) 이후 경과 일수

위의 내용은 한국천문연구원에서 편찬한 2013년 역서의 page 113에 나온 식이다.

먼저 위 두 식의 의미를 이해할 수 있어야 이 식을 사용할수 있다.


이 식들의 의미는 다음과 같다.

먼저, d=0(UT 2012년 12월 31일 0시 0분, KST 09:00)일때
태양의 평균경도는 279.927287도이고 태양의 근점이각은 356.666444도라는 뜻이다.

그리고 하루가 지날때마다(=d가 1증가할때마다)
평균경도는 0.98564736도 증가하고, 근점이각은 0.98560028도 증가한다.
평균경도의 1일 증가분 0.98564736도는 360도를  1회귀년(365.242190일)의 값으로 나눈 값이다. 따라서 365.242190일간격으로 같은 평균경도를 갖게되고 평균경도의 범위는 0~360도이다.
근점이각의 1일 증가분 0.98560028도는 360도를 1근점년(365.259636일)으로 나눈 값이다.
이 역시 365.259636일 간격으로 같은 근점이각을 갖게되고 범위는 0~360도이다.

(여기서 사용한 1회귀년과 1근점년의 길이는 2013년 역서의 page 113 에 있는 값이다.)

평균경도 식(1)에서 상수항 279.827287은 2013년 3월의 춘분시각에 태양의 진황경을 0도로 맞추기 위해서 미리 계산된 상수값이다.
근점이각 식(2)에서 상수항 356.666444의 값 역시 2013년의 1월의 근일점(태양과 지구가 가장 가까워지는 시점)에 0도의 값을 갖도록 미리 계산된 상수값이다.

이제 평균경도와 근점이각이 0이 되는 d의 값을 구하면 다음과 같다.

(1)에서 279.827287 + 0.98564736*d=0
 ==>  0.98564736*d=-279.827287 ==> d=-279.827287/0.98564736 = -283.9020306

(2)에서  356.666444 + 0.98560028*d=0
==> 0.98560028*d=-356.666444 ==>  d=-356.666444/0.98560028 = -361.8773769

이 값들의 의미는 UT 2012년 12월 31일 0시에서
 -283.9020306일 떨어진 시각(2012-03-22 02:21분, KST 11:21)이 태양의 평균경도(평균황경)이 0도이고,
 -361.8773769일 떨어진 시각(2012-01-04 02:56분, KST 11:56)이 태양의 근점이각이 0도라는 뜻이다.

이제 시험삼아 평균경도와 근점이각이 360이 되는 d의 값을 구하면 각각 다음과 같다.

(1)에서 279.827287 + 0.98564736*d=360
 ==>  d=(360-279.827287)/0.98564736 = 81.34015902

(2)에서  356.666444 + 0.98560028*d=360
==>   d=(360-356.666444)/0.98560028 = 3.382259591

2012년 12월 31일 0시에 위의 각 값들을 더하면
2013-03-22 08:09분(KST 17:09)에 태양이 있는 위치의  평균황경 0도(=360도)이고, 2013-01-03 9:10분(KST 18:10)에 태양이 있는 위치의 근점이각이 0도(=360도)이다.

(1),(2)식의 결과값은 d의 일차식이기 때문에  d의 값에 따라 평균황경과 근점이각이 -무한대에서 +무한대까지의 값을 갖는다. 하지만 평균경도나 근점이각의 범위는 0~360도까지이므로, 결과값을 360으로 나눈 나머지만 사용한다.

이제 특정시각이 지정되면 해당 시각에서부터 2012년 12월 31일 0시까지의 떨어진 시간 d를 구한뒤 이 d를 (1),(2)식에 대입하면 지정된 특정시각의 평균황경과 근점이각을 알수 있다.

이상이 (1),(2)식의 의미이다.


이제 진황경을 구하는 방법을 알아보자.

만일 지구가 태양의 주위를 등속 원궤도로 공전한다면 평균황경이 진황경이다.
하지만 지구가 불규칙한 속도의 타원궤도로 공전하기 때문에 평균황경과 진황경은 차이가 있다.

지구의 공전궤도에서 평균황경과 진황경이 일치하는지점은 지구가 태양에 가장 가까운 지점(근점에 있을때)과 태양에서 가장 먼지점(원점에 있을때)의 두시점뿐이다.

지구가 타원궤도를 그리며 태양 주변을 공전할때 근점에 있을때를 0도라 하고, 원점에 있을때를 180도라고 하면 지구는 근점(근점이각 0도,1월 3~4일경)에서 시작해서 봄(근점이각 대략 90도)을 지나 여름에 원점(근점이각 180도, 7월 초)을 지나고 가을(근점이각 대략 270도)을 지나서 다시 근점와 와서 타원궤도를 한바퀴 돌고 1년이 지나간다.

식(1)에서 지구는 1일에 평균 0.98564736도씩 태양 주변을 이동한다고 했는데, 실제로는 근점에 가까울때 가장 빨리 이동하고, 근점에서 가장 멀때 가장 느리게 이동한다.
(이 사실은 로마시대부터 알려져 있고, 북반구의 겨울은 지구가 공전궤도를 빨리 지나가서 짧고, 여름에는 지구가 느리게 지나가서 길다. 그래서 율리우스력에서 2월이 짧고, 7,8월이 길게 만들어졌다. 괜히 로마황제들이 마음대로 지정해서 날짜가 길어지고 짧아진 것이 아니다.)

이제 태양의 황경 변화량을  근점에서의 하루의 변화량은 조금 늘리고, 원점에서는 조금 줄이면 진태양황경을 구할 수 있다. 그런데 근점에서 얼마나 늘리고, 원점에서 얼마나 줄여야 될까? 이 문제는 케플러가 해결했다.

( 우선 참고삼아 진태양황경은 이런 형식으로 구할수 있다는것을 생각해 두자.
진태양황경=평균황경+(특정계수)*sin(근점이각)으로 표현할수 있다.
근점이각이 0도일때와 180도일때, sin(0도)와 sin(180도)의 값은 0이므로 진황경과 평균황경이 같아진다. 그리고 각각의 근점이각에 따라 특정계수만큼의 보정값을 평균황경에 더하면 진황경이 된다. )

케플러는 

 E = M + e*sin(E)              식 (3)

라는 케플러 방정식을 제시했다.
( E:이심율이각, M:근점이각, e:타원궤도의 이심률)

이 식은 특정시점에 태양의 평균황경과 진황경의 차이를 E라고 하고 우리말로 이심율이각(혹은 중심차)라고 표시하면, 이 E를 구할수 있다는 의미이다. 하지만 왜 좌변과 우변에 같은 E(이심률이각)가 있을까?

우선 특정 시점이 지정되면 식(1),(2)를 이용하여 해당 시점의 평균황경(L1)과 근점이각(M1)을  구할수 있고, 근점이각에 의해서 계산되는 보정값(이심률이각)을 평균황경(L1)과 근점이각(M1)에 더하면 새로운, 진황경에 좀더 가까운 평균황경(L2)과 근점이각(M2)을 구할수 있고,  또 다시 L2, M2에 의해서 계산되는 보정값을 다시 L2,M2에 더하면 L3, M3가 나오고.. 이를 계속 반복해서 새로운 값L,M과 바로 직전의 L,M의 값의 차이가 극히 미미할때, L값이 진황경이고, M값이 진근점이각이다.
따라서 식(3)의 이심율이각은

 E0=M
 E1=M+e*sin(E0)
 E2=M+e*sin(E1)
 E3=M+e*sin(E2)

의 형태로 계속 진행되어 원하는 정밀도에 이를때까지 계산한다.
이렇게 해서 마지막으로 계산된 이심율이각 E값을 평균황경에 더하면 진황경을 구할수 있다.
케플러 방정식의 풀이방법은 다양한 방법이 있다. 이와 관련한 사항은 여기서는 더이상 언급하지 않는다. 더 자세한 사항은 인터넷을 검색하거나 관련서적을 찾아보자.

 특정 이심률e가 지정되면 이렇게 무한급수로 계산하지 않고, 특정 공식을 이용해서 이심율이각을 계산할수 있다.

지구공전궤도의 이심률이 대략 0.016703일때,
E=1.919도*sin(M) + 0.020도*sin(2*M)-0.005도*sin(M*T)    식(4)
(E는 도, 2000년대에 T=1)
로 계산할수 있다. 이 식은 역법의 원리분석 page 288의 (13)식이다.

이제 특정시점이 주어지면 식(1),(2)에서 평균황경 L, 근점이각 M을 구할수 있고, 근점이각 M을 이용해서 식(4)에서 이심율이각 E를 구하고, 진황경은 L+E의 값이 된다.

이제 2013역서에 나온 춘분점 KST 3월 20일 20시 2분(UT 11:02)의 진황경을 구해본다.
2012-12-31 0:0부터 2013-03-20 11:02분까지의 차이 d=79.45972222이다.

식(1)에서  평균황경 L=358.1465524도
식(2)에서  근점이각 M=74.98196847도=1.3086822라디안
식(4)에서 이심률이각
   E=1.919도*sin(1.3086822라디안)+0.020도*sin(2*1.3086822라디안)
      -0.005*sin(1.3086822라디안)
     =1.8586369도

이제 진태양황경=L+E=358.1465524도+1.8586369도=360.0051894도=0.0051894도

결국 KST 2013년 3월 20일 20시 02분(UT 11:02분)의 태양 진황경은=0.0051894도이다.
각도로 0.0051894 차이는 시간으로 약 7분 30초 정도이다.

여기서 설명한 방식으로 2013년 3월의 춘분점을 계산하면 대략 3월 20일 19:55분무렵이라는 의미이다. ( 하고 싶은 말은 많지만 정밀도와 오차범위는 이 글에서 논의하지 않는다.)

지금까지 특정시점 d가 주어졌을때, 진태양황경을 구할수 있는 방법에 대해서 설명했다.

이와 반대로 진태양황경이 주어졌을때 d를 구하는 특별한 방법은 존재하지 않는다.
따라서 24절기의 입기 시각을 구하려면 임의의 시각을 선택해서 이를  d로 변환하고, 이시각의  진태양황경이 15도의 배수에 해당하는지를 파악하는 수 밖에 없다.

만일 2013년 2월 입춘 입기시각을 구하려면 다음과 같은 과정으로 구한다.
우선 2월 초의 특정한 하루의 0시를 선택해서 d0라 하고, 그 다음날의 0시를 d1으로 놓는다.
입춘의 진태양황경은 315도이다.
이제 위에서 설명한대로 계산해서 d0의 진태양황경이 315도보다 작고, d1의 진태양황경이 315도보다 큰 d0를 선택한다. 만일 d0와 d1이 둘다 315도보다 작다면 d0와 d1을 1일씩 증가시키고, 둘다 315도보다 크다면 1일씩 감소시키면 된다.(d0를 2월1일로 선택했다면 진황경 계산회수는 대략 5~6회이내이다),
이제 진태양황경 315도를 포함하는 d0,d1이 계산되면 d0날의 특정시간을 dh0로 놓고 한시간 뒤를 dh1로 놓고, dh0, dh1시각의 진태양황경을 구하고, 진태양황경 315도가 dh0와 dh1사이에 존재하는 dh0를 찾는다.(최대 계산회수 23회), 같은 방식으로 분단위로 dhm0와 dhm1을 계산한다.(최대계산회수 59회).  진태양황경 315도를 포함하는 dhm0와 dhm1을 찾아내면 해당 dhm이 입춘의 입기시각이다. 대략 88회 정도 계산하면 태양이 진태양황경 315도를 지나는 시각을 찾아낼 수 있다.
이 88여회의 계산은 컴퓨터가 하도록 프로그램을 짜면 된다. 능력있는 프로그래머라면 2진검색방법으로 15~16회 이하의 진태양황경 구하기 계산으로 하나의 절기시각을 찾아낼 수 있으리라 생각한다.

이상이 정확한(?) 24절기 구하기 방법이다.

이글은 제가 혼자서 독자적으로 작성한 글입니다.
이 글을 퍼가실때에는 작성자의 이름을 꼭 표시하고, 글의 출처 링크를 꼭 밝혀주시길 바랍니다.

2013년 7월 19일
고영창 작성
Email : kohyc@chol.com
홈페이지 http://user.chol.com/~kohyc