新增函数&修复bug

This commit is contained in:
2020-07-14 15:38:51 +08:00
parent 08fcbf9f89
commit 1c4397c9dc
12 changed files with 526 additions and 551 deletions
+47 -4
View File
@@ -70,12 +70,12 @@ func dt_cal(y float64) float64 { //传入年, 返回世界时UT与原子时(
if y >= y0 {
jsd := float64(31) // sjd是y1年之后的加速度估计
// 瑞士星历表jsd=31, NASA网站jsd=32, skmap的jsd=29
if y > y0+100 {
if y > y0+100.00 {
return dt_ext(y, jsd)
}
v := dt_ext(y, jsd) //二次曲线外推
dv := dt_ext(y0, jsd) - t0 // ye年的二次外推与te的差
return (v - dv*(y0+100-y)/100)
return (v - dv*(y0+100.00-y)/100.00)
}
d := dt_at
var i int
@@ -85,7 +85,7 @@ func dt_cal(y float64) float64 { //传入年, 返回世界时UT与原子时(
// 判断年所在的区间
}
}
t1 := (y - d[i]) / (d[i+5] - d[i]) * 10 //////// 三次插值, 保证精确性
t1 := (y - d[i]) / (d[i+5] - d[i]) * 10.00 //////// 三次插值, 保证精确性
t2 := t1 * t1
t3 := t2 * t1
res := d[i+1] + d[i+2]*t1 + d[i+3]*t2 + d[i+4]*t3
@@ -94,7 +94,8 @@ func dt_cal(y float64) float64 { //传入年, 返回世界时UT与原子时(
func DeltaT(Date float64, IsJDE bool) (Result float64) { //传入年或儒略日,传出为秒
var Year float64
if IsJDE {
Year = (Date-2451545.0)/365.25 + 0.1 + 2000
dates := JDE2Date(Date)
Year = float64(dates.Year()) + float64(dates.YearDay())/365.0
} else {
Year = Date
}
@@ -103,6 +104,7 @@ func DeltaT(Date float64, IsJDE bool) (Result float64) { //传入年或儒略日
return
}
if Year < 2100 && Year >= 2010 {
//fmt.Println(Year)
Result = dt_cal(Year) //-3.2-(Year-2017)*0.029915;
return
}
@@ -167,6 +169,41 @@ func JDE2Date(JD float64) time.Time {
return dates
}
func JDE2DateByZone(JD float64, tz *time.Location) time.Time {
JD = JD + 0.5
Z := float64(int(JD))
F := JD - Z
var A, B, Years, Months, Days float64
if Z < 2299161.0 {
A = Z
} else {
alpha := math.Floor((Z - 1867216.25) / 36524.25)
A = Z + 1 + alpha - math.Floor(alpha/4)
}
B = A + 1524
C := math.Floor((B - 122.1) / 365.25)
D := math.Floor(365.25 * C)
E := math.Floor((B - D) / 30.6001)
Days = B - D - math.Floor(30.6001*E) + F
if E < 14 {
Months = E - 1
}
if E == 14 || E == 15 {
Months = E - 13
}
if Months > 2 {
Years = C - 4716
}
if Months == 1 || Months == 2 {
Years = C - 4715
}
tms := (Days - math.Floor(Days)) * 24 * 3600
Days = math.Floor(Days)
dates := time.Date(int(Years), time.Month(int(Months)), int(Days), 0, 0, 0, 0, tz)
dates = time.Unix(dates.Unix()+int64(tms), int64((tms-math.Floor(tms))*1000000000))
return dates
}
func GetLunar(year, month, day int) (lmonth, lday int, leap bool, result string) {
jde := JDECalc(year, month, float64(day)) //计算当前JDE时间
if month == 11 || month == 12 { //判断当前日期属于前一年周期还是后一年周期
@@ -321,3 +358,9 @@ func GetSolar(year, month, day int, leap bool) float64 {
jde := moon[min-1+month] + float64(day) - 1
return jde
}
// Date2JDE 日期转儒略日
func Date2JDE(date time.Time) float64 {
day := float64(date.Day()) + float64(date.Hour())/24.0 + float64(date.Minute())/24.0/60.0 + float64(date.Second())/24.0/3600.0 + float64(date.Nanosecond())/1000000000.0/3600.0/24.0
return JDECalc(date.Year(), int(date.Month()), day)
}
+1 -1
View File
@@ -6,5 +6,5 @@ import (
)
func Test_LoBo(t *testing.T) {
fmt.Sprintf("%.14f", LoToRa(22, 33, 2451545.0))
fmt.Printf("%.9f", dt_cal(2020.5))
}
+2 -2
View File
@@ -1420,9 +1420,9 @@ func GetMoonRiseTime(JD, Lon, Lat, TZ, ZS float64) float64 {
}
if moonheight < 0 && moonang > 180 {
tms = (180 - moonang) / 15
JD1 += (tms/24 + (tms/24*12.0)/15.0/24.0)
JD1 = JD1 - (tms/24 + (tms/24*12.0)/15.0/24.0)
} else if moonheight < 0 && moonang < 180 {
tms = (-180 + moonang) / 15
tms = (180 - moonang) / 15
JD1 += (tms/24 + (tms/24*12.0)/15.0/24.0)
}
now := MoonTimeAngle(JD1, Lon, Lat, TZ)
+22 -9
View File
@@ -903,21 +903,34 @@ func HSunTrueLo(JD float64) float64 {
return L
}
func HSunSeeLo(JD float64) float64 {
t := (JD - 2451545) / 365250.0
L := HSunTrueLo(JD)
R := planet.WherePlanet(-1, 2, JD)
t2 := t * t
t3 := t2 * t //千年数的各次方
R += (-0.0020 + 0.0044*t + 0.0213*t2 - 0.0250*t3)
L = L + HJZD(JD) - 20.4898/R/3600
func HSunTrueBo(JD float64) float64 {
L := planet.WherePlanet(0, 1, JD)
return L
}
func HSunSeeLo(JD float64) float64 {
L := HSunTrueLo(JD)
/*
t := (JD - 2451545) / 365250.0
R := planet.WherePlanet(-1, 2, JD)
t2 := t * t
t3 := t2 * t //千年数的各次方
R += (-0.0020 + 0.0044*t + 0.0213*t2 - 0.0250*t3)
L = L + HJZD(JD) - 20.4898/R/3600.00
*/
L = L + HJZD(JD) + SunLoGXC(JD)
return L
}
func SunLoGXC(JD float64) float64 {
R := planet.WherePlanet(0, 2, JD)
return -20.49552 / R / 3600
}
func EarthAway(JD float64) float64 {
//t=(JD - 2451545) / 365250;
//R=Earth_R5(t)+Earth_R4(t)+Earth_R3(t)+Earth_R2(t)+Earth_R1(t)+Earth_R0(t);
return planet.WherePlanet(0, 2, 2555555)
return planet.WherePlanet(0, 2, JD)
}
func HSunSeeRaDec(JD float64) (float64, float64) {
+8 -4
View File
@@ -7,10 +7,14 @@ import (
)
func Test_Jq(t *testing.T) {
//fmt.Println(GetOneYearJQ(2019))
fmt.Println(JDE2Date(GetWHTime(2019, 10)))
fmt.Println(JDE2Date(GetJQTime(2019, 15)))
fmt.Println(JDE2Date(GetJQTime(2019, 0)))
data := GetOneYearJQ(2019)
for i := 1; i < 25; i++ {
fmt.Println(JDE2Date(data[i]))
}
//fmt.Println(JDE2Date(GetWHTime(2019, 10)))
//fmt.Println(JDE2Date(GetJQTime(2020, 0)))
//date := TD2UT(GetJQTime(2020, 0), true)
//fmt.Println(HSunSeeLo(date))
}
func Test_SunLo(t *testing.T) {