1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
|
/*
* Copyright (C) 2014 ~ 2018 Deepin Technology Co., Ltd.
*
* Author: jouyouyun <jouyouwen717@gmail.com>
*
* This program is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* any later version.
*
* This program is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program. If not, see <http://www.gnu.org/licenses/>.
*/
package lunar
import (
"math"
"github.com/linuxdeepin/go-lib/calendar/util"
)
/**
* 按儒略日计算地球的日心黄经
*
* 参数： jd 儒略日
* 返回 地球的日心黄经，单位是弧度(rad)
*/
func GetSunEclipticLongitudeForEarth(jd float64) float64 {
t := util.GetJulianThousandYears(jd)
L0 := GetEarthL0(t)
L1 := GetEarthL1(t)
L2 := GetEarthL2(t)
L3 := GetEarthL3(t)
L4 := GetEarthL4(t)
L5 := GetEarthL5(t)
L := ((((L5*t+L4)*t+L3)*t+L2)*t+L1)*t + L0
return Mod2Pi(L)
}
/**
* 按儒略日计算地球的日心黄纬
*
* 参数 jd 儒略日
* 返回 地球的日心黄纬，单位是弧度(rad)
*/
func GetSunEclipticLatitudeForEarth(jd float64) float64 {
t := util.GetJulianThousandYears(jd)
B0 := GetEarthB0(t)
B1 := GetEarthB1(t)
B2 := GetEarthB2(t)
B3 := GetEarthB3(t)
B4 := GetEarthB4(t)
B := ((((B4*t)+B3)*t+B2)*t+B1)*t + B0
return B
}
/**
* 按照儒略日计算地球和太阳的距离
*
* 参数 jd 儒略日
* 返回 地球和太阳的距离，单位是天文单位(au)
*/
func GetSunRadiusForEarth(jd float64) float64 {
t := util.GetJulianThousandYears(jd)
R0 := GetEarthR0(t)
R1 := GetEarthR1(t)
R2 := GetEarthR2(t)
R3 := GetEarthR3(t)
R4 := GetEarthR4(t)
R5 := GetEarthR5(t)
R := ((((R5*t+R4)*t+R3)*t+R2)*t+R1)*t + R0
return R
}
/**
* 用于把vsop87理论算出来的经度转换成fk5目视系统的经度的修正值，参考 Jean Meeus 的 Astronomical
* Algorithms 第二版(1998)第32章219页(32.3)式
*
* 参数 l
* vsop87经度(rad)
* 参数 b
* vsop87纬度(rad)
* 参数 jd
* 儒略日
* 返回 修正量(rad)
*/
func Vsop2Fk5LongitudeCorrection(l float64, b float64, jd float64) float64 {
t := util.GetJulianCentury(jd)
lp := l - ToRadians(1.397)*t - ToRadians(0.00031)*t*t
return SecondsToRadians(-0.09033 + 0.03916*(math.Cos(lp)+math.Sin(lp))*math.Tan(b))
}
/**
* 计算修正后的太阳的地心视黄经
*
* 参数 jd
* 儒略日
* 返回 修正后的地心黄经(rad)
*/
// 常量
var lightAberration = SecondsToRadians(20.4898)
func GetEarthEclipticLongitudeForSun(jd float64) float64 {
// 计算地球的日心黄经
l := GetSunEclipticLongitudeForEarth(jd)
// 计算地球的日心黄纬
b := GetSunEclipticLatitudeForEarth(jd)
// 修正章动
l += CalcEarthLongitudeNutation(util.GetJulianCentury(jd))
// 转换到fk5
l += Vsop2Fk5LongitudeCorrection(l, b, jd)
// 转换成太阳的地心黄经
l = Mod2Pi(l + math.Pi)
// 计算光行差
// 计算日地距离
r := GetSunRadiusForEarth(jd)
// 太阳到地球的光行差参数
l -= lightAberration / r
return l
}
|
