package basic import ( "math" "b612.me/astro/planet" . "b612.me/astro/tools" ) const ( MERCURY_S_PERIOD = 1 / ((1 / 87.9691) - (1 / 365.256363004)) mercuryConjunctionDerivativeStepDay = 2e-5 * 36525.0 mercuryLightTimeDaysPerAU = 0.0057755183 mercuryEventSearchN = 16 mercuryStationWindowDays = 30.0 mercuryStationDerivativeStepDay = 0.01 mercuryStationCoarseStepDay = 2.0 mercuryStationHalfWindowDay = 2.0 mercuryStationMotionTolerance = 1e-3 // mercuryStationAnchorAttempts 类型化「留」在锚点不满足侧向不变量时,最多推进/回退几个会合周期。 mercuryStationAnchorAttempts = 4 // mercuryConjunctionSameInstantDegrees 合搜索「同刻」快速路径的视黄经差触发阈值(度)。 // 0.05 度约合 30 分钟的时间跨度,远大于截断级数在根附近的误差, // 从而保证「查询落在合附近」时一定走精确解 + 侧向判定,而不是被方向扫描跨过去。 mercuryConjunctionSameInstantDegrees = 5.0e-2 // mercuryConjunctionScanStepDay / mercuryConjunctionScanMaxSteps 方向性括号扫描参数: // 8 天一步 × 80 步 = 640 天,覆盖一个水星会合周期(115.88 天)且有大量余量。 mercuryConjunctionScanStepDay = 8.0 mercuryConjunctionScanMaxSteps = 80 // mercuryConjunctionPolishToleranceDay 全项抛光的时间容差(0.01 秒): // 割线法从 8 天括号收敛到该量级只需多 1~2 次求值,但把合时刻精度保持在毫秒级 // (与旧实现牛顿法的 1e-5 天口径一致)。 mercuryConjunctionPolishToleranceDay = 0.01 / 86400.0 ) type mercuryConjunctionLBR struct { lo float64 bo float64 r float64 } type mercuryConjunctionGeo struct { lo float64 bo float64 dist float64 } type mercuryConjunctionResult struct { diff float64 sunLightDays float64 geoLightDays float64 } func mercuryHelioN(planetIndex int, jd float64, n int) mercuryConjunctionLBR { return mercuryConjunctionLBR{ lo: planet.WherePlanetN(planetIndex, 0, jd, n), bo: planet.WherePlanetN(planetIndex, 1, jd, n), r: planet.WherePlanetN(planetIndex, 2, jd, n), } } func mercuryGeocentric(planetPos, earthPos mercuryConjunctionLBR) mercuryConjunctionGeo { x := planetPos.r*Cos(planetPos.bo)*Cos(planetPos.lo) - earthPos.r*Cos(earthPos.bo)*Cos(earthPos.lo) y := planetPos.r*Cos(planetPos.bo)*Sin(planetPos.lo) - earthPos.r*Cos(earthPos.bo)*Sin(earthPos.lo) z := planetPos.r*Sin(planetPos.bo) - earthPos.r*Sin(earthPos.bo) dist := math.Sqrt(x*x + y*y + z*z) return mercuryConjunctionGeo{ lo: Limit360(math.Atan2(y, x) * 180 / math.Pi), bo: math.Atan2(z, math.Sqrt(x*x+y*y)) * 180 / math.Pi, dist: dist, } } func mercuryConjunctionAngleDelta(diff float64) float64 { diff = Limit360(diff) if diff > 180 { diff -= 360 } if diff < -180 { diff += 360 } return diff } func mercuryConjunctionHeliocentricDelta(jd, targetDeg float64, n int) float64 { planetLo := planet.WherePlanetN(1, 0, jd, n) earthLo := planet.WherePlanetN(-1, 0, jd, n) return mercuryConjunctionAngleDelta(planetLo - earthLo - targetDeg) } func mercuryConjunctionDifference(jd float64, n int, targetDeg, sunLightDays, geoLightDays float64) mercuryConjunctionResult { earthForSun := mercuryHelioN(-1, jd-sunLightDays, n) sunLo := Limit360(earthForSun.lo + 180) earth := mercuryHelioN(-1, jd-geoLightDays, n) planetPos := mercuryHelioN(1, jd-geoLightDays, n) geo := mercuryGeocentric(planetPos, earth) return mercuryConjunctionResult{ diff: mercuryConjunctionAngleDelta(geo.lo - sunLo - targetDeg), sunLightDays: earthForSun.r * mercuryLightTimeDaysPerAU, geoLightDays: geo.dist * mercuryLightTimeDaysPerAU, } } func mercuryConjunctionExactDelta(jd float64) float64 { return mercuryConjunctionAngleDelta(MercuryApparentLo(jd) - HSunApparentLo(jd)) } func mercuryConjunctionApproxTT(seed float64, inferior bool) float64 { heliocentricTarget := 180.0 if inferior { heliocentricTarget = 0 } jd := seed for i := 0; i < 6; i++ { jd -= mercuryConjunctionHeliocentricDelta(jd, heliocentricTarget, 8) / (360.0 / MERCURY_S_PERIOD) } startSample := mercuryConjunctionDifference(jd, 8, 0, 0, 0) nextSample := mercuryConjunctionDifference(jd+mercuryConjunctionDerivativeStepDay, 8, 0, 0, 0) diffSlope := mercuryConjunctionAngleDelta(nextSample.diff-startSample.diff) / mercuryConjunctionDerivativeStepDay refined := mercuryConjunctionDifference(jd, 40, 0, startSample.sunLightDays, startSample.geoLightDays) jd -= refined.diff / diffSlope final := mercuryConjunctionDifference(jd, -1, 0, refined.sunLightDays, refined.geoLightDays) jd -= final.diff / diffSlope return jd } func mercuryConjunctionExactTT(seed float64, inferior bool) float64 { estimateJD := mercuryConjunctionApproxTT(seed, inferior) converged := false for i := 0; i < eventNewtonMaxIterations; i++ { prevJD := estimateJD longitudeDelta := mercuryConjunctionExactDelta(prevJD) longitudeSlope := (mercuryConjunctionExactDelta(prevJD+0.000005) - mercuryConjunctionExactDelta(prevJD-0.000005)) / 0.00001 nextJD := prevJD - longitudeDelta/longitudeSlope estimateJD = nextJD if math.Abs(nextJD-prevJD) <= 0.00001 { converged = true break } } if !converged { return math.NaN() } return estimateJD } // mercuryConjunction 在 jde 的指定方向上求最近一次水星合(内合/外合)。 // // 旧实现用「启发式跳 + 2 天一步走」:跳过头落到合之后时,|Δ| 会持续变大, // 于是走满一个会合周期、跳过一次合(例如 2008-01-22 查到 2008-06-07 而不是 2008-02-06), // 并且下游的类型化「留」会因此返回错事件、甚至让 Last 返回未来。 // 新实现改成方向性括号扫描:截断级数逐段找异号区间(顺序扫描 ⇒ 一定取最近的一个合), // 再用全项级数在括号内抛光;扫满 640 天仍无括号时有界返回 NaN。 func mercuryConjunction(jde float64, next uint8) float64 { //0=last 1=next if !isFiniteFloat(jde) { return math.NaN() } if math.Abs(mercuryConjunctionExactDelta(jde)) <= mercuryConjunctionSameInstantDegrees { best := math.NaN() consider := func(inferior bool) { eventUT := TD2UT(mercuryConjunctionExactTT(jde, inferior), false) if !isFiniteFloat(eventUT) { return } if next == 0 && !eventUTQueryBeforeOrEqual(eventUT, jde) { return } if next == 1 && !eventUTQueryAfterOrEqual(eventUT, jde) { return } if math.IsNaN(best) || math.Abs(eventUTQueryTTDelta(eventUT, jde)) < math.Abs(eventUTQueryTTDelta(best, jde)) { best = eventUT } } consider(true) consider(false) if !math.IsNaN(best) { return best } } direction := 1.0 if next == 0 { direction = -1 } leftJD := jde leftValue := mercuryConjunctionDeltaN(leftJD, mercuryEventSearchN) if !isFiniteFloat(leftValue) { return math.NaN() } for i := 0; i < mercuryConjunctionScanMaxSteps; i++ { rightJD := jde + direction*mercuryConjunctionScanStepDay*float64(i+1) rightValue := mercuryConjunctionDeltaN(rightJD, mercuryEventSearchN) if !isFiniteFloat(rightValue) { return math.NaN() } if leftValue == 0 || rightValue == 0 || leftValue*rightValue < 0 { return mercuryConjunctionPolish(jde, leftJD, rightJD, direction) } leftJD, leftValue = rightJD, rightValue } return math.NaN() } // mercuryConjunctionDeltaN 截断级数下的水星-太阳视黄经差(度,[-180,180]),用于方向性括号扫描。 func mercuryConjunctionDeltaN(jd float64, n int) float64 { return mercuryConjunctionAngleDelta(MercuryApparentLoN(jd, n) - HSunApparentLoN(jd, n)) } // mercuryConjunctionPolish 用全项级数在截断级数给出的括号内抛光。 // 截断误差可能让括号两端在全项函数上同号(罕见),此时沿扫描方向再扩一两个扫描步; // 若仍未被确认(典型情形:查询几乎正好落在合上,截断级数在根两侧的符号与全项不一致), // 退回全项级数的方向扫描,保证有界且不返回 NaN。 func mercuryConjunctionPolish(jde, leftJD, rightJD, direction float64) float64 { for attempt := 0; attempt < 3; attempt++ { leftValue := mercuryConjunctionExactDelta(leftJD) rightValue := mercuryConjunctionExactDelta(rightJD) if !isFiniteFloat(leftValue) || !isFiniteFloat(rightValue) { return math.NaN() } if leftValue == 0 || rightValue == 0 || leftValue*rightValue < 0 { root, ok := eventBracketSecantRoot(leftJD, rightJD, leftValue, rightValue, mercuryConjunctionPolishToleranceDay, mercuryConjunctionExactDelta) if !ok { return math.NaN() } return TD2UT(root, false) } if direction > 0 { rightJD += mercuryConjunctionScanStepDay continue } leftJD -= mercuryConjunctionScanStepDay } return mercuryConjunctionFullDirectionalScan(jde, direction) } // mercuryConjunctionFullDirectionalScan 全项级数的方向扫描(截断括号未被确认时的兜底)。 func mercuryConjunctionFullDirectionalScan(jde, direction float64) float64 { leftJD := jde leftValue := mercuryConjunctionExactDelta(leftJD) if !isFiniteFloat(leftValue) { return math.NaN() } for i := 0; i < mercuryConjunctionScanMaxSteps; i++ { rightJD := jde + direction*mercuryConjunctionScanStepDay*float64(i+1) rightValue := mercuryConjunctionExactDelta(rightJD) if !isFiniteFloat(rightValue) { return math.NaN() } if leftValue == 0 || rightValue == 0 || leftValue*rightValue < 0 { root, ok := eventBracketSecantRoot(leftJD, rightJD, leftValue, rightValue, mercuryConjunctionPolishToleranceDay, mercuryConjunctionExactDelta) if !ok { return math.NaN() } return TD2UT(root, false) } leftJD, leftValue = rightJD, rightValue } return math.NaN() } func LastMercuryConjunction(jde float64) float64 { return inclusiveLastSimpleEvent(jde, LastMercuryConjunctionStrict, NextMercuryConjunctionStrict) } func NextMercuryConjunction(jde float64) float64 { return inclusiveNextSimpleEvent(jde, LastMercuryConjunctionStrict, NextMercuryConjunctionStrict) } func LastMercuryConjunctionStrict(jde float64) float64 { return mercuryConjunction(jde, 0) } func NextMercuryConjunctionStrict(jde float64) float64 { return mercuryConjunction(jde, 1) } func NextMercuryInferiorConjunction(jde float64) float64 { date := NextMercuryConjunctionStrict(jde) if EarthMercuryAway(date) > EarthAway(date) { return NextMercuryConjunctionStrict(date + 2) } return date } func NextMercurySuperiorConjunction(jde float64) float64 { date := NextMercuryConjunctionStrict(jde) if EarthMercuryAway(date) < EarthAway(date) { return NextMercuryConjunctionStrict(date + 2) } return date } func LastMercuryInferiorConjunction(jde float64) float64 { date := LastMercuryConjunctionStrict(jde) if EarthMercuryAway(date) > EarthAway(date) { return LastMercuryConjunctionStrict(date - 2) } return date } func LastMercurySuperiorConjunction(jde float64) float64 { date := LastMercuryConjunctionStrict(jde) if EarthMercuryAway(date) < EarthAway(date) { return LastMercuryConjunctionStrict(date - 2) } return date } func mercuryRADerivative(jde, delta float64) float64 { sub := MercuryApparentRa(jde+delta) - MercuryApparentRa(jde-delta) if sub > 180 { sub -= 360 } if sub < -180 { sub += 360 } return sub / (2 * delta) } func mercuryRADerivativeN(jde, delta float64, n int) float64 { sub := MercuryApparentRaN(jde+delta, n) - MercuryApparentRaN(jde-delta, n) if sub > 180 { sub -= 360 } if sub < -180 { sub += 360 } return sub / (2 * delta) } func mercuryStationInWindow(startTT, endTT float64) float64 { bestJD := zeroEventInWindow(startTT, endTT, mercuryStationCoarseStepDay, mercuryStationHalfWindowDay, 30.0/86400.0, func(jd float64) float64 { return mercuryRADerivativeN(jd, mercuryStationDerivativeStepDay, mercuryEventSearchN) }, func(jd float64) float64 { return mercuryRADerivative(jd, mercuryStationDerivativeStepDay) }) return TD2UT(bestJD, false) } func mercuryStationBetween(startTT, endTT float64) bool { if endTT < startTT { startTT, endTT = endTT, startTT } if endTT-startTT <= 0 { return false } // 跨度超过一个回留窗口时截断扫描已不可靠:返回 false 不再阻断快速路径。 // 此时候选事件由 mercuryConjunction 的完整方向扫描保证,快速路径无需再复核。 // A span beyond the station window makes the truncated scan unreliable; report false so the // fast path is not blocked. The typed candidate is already guaranteed by mercuryConjunction's // exhaustive directional scan. if endTT-startTT > mercuryStationWindowDays { return false } // 截断扫描足以判断单候选快速路径是否安全 / A truncated scan is enough to decide whether the one-candidate fast path is safe. left := startTT leftValue := mercuryRADerivativeN(left, mercuryStationDerivativeStepDay, mercuryEventSearchN) for left < endTT { right := left + mercuryStationCoarseStepDay if right > endTT { right = endTT } rightValue := mercuryRADerivativeN(right, mercuryStationDerivativeStepDay, mercuryEventSearchN) if leftValue == 0 || leftValue*rightValue < 0 || rightValue == 0 { return true } left = right leftValue = rightValue } return false } func mercuryProgradeToRetrogradeAroundInferior(inferiorUT float64) float64 { inferiorTT := TD2UT(inferiorUT, true) return mercuryStationInWindow(inferiorTT-mercuryStationWindowDays, inferiorTT) } func mercuryRetrogradeToProgradeAroundInferior(inferiorUT float64) float64 { inferiorTT := TD2UT(inferiorUT, true) return mercuryStationInWindow(inferiorTT, inferiorTT+mercuryStationWindowDays) } // 类型化「留」统一结构:以下合为锚点求候选,再用侧向不变量(Next 不得早于查询、 // Last 不得晚于查询)校验;不满足则把锚点推进/回退一个会合周期重试,最多 // mercuryStationAnchorAttempts 次,最终仍有界失败为 NaN。 // 这保证 Last* 永远不会返回未来、Next* 永远不会返回过去的事件。 func NextMercuryProgradeToRetrograde(jde float64) float64 { if !isFiniteFloat(jde) { return math.NaN() } inferior := NextMercuryInferiorConjunction(jde) for i := 0; i < mercuryStationAnchorAttempts && isFiniteFloat(inferior); i++ { date := mercuryProgradeToRetrogradeAroundInferior(inferior) if isFiniteFloat(date) && stationUTQueryAfterOrEqual(date, jde) { return date } inferior = NextMercuryInferiorConjunction(eventUTNextQueryTT(inferior)) } return math.NaN() } func NextMercuryRetrogradeToPrograde(jde float64) float64 { if !isFiniteFloat(jde) { return math.NaN() } inferior := LastMercuryInferiorConjunction(jde) for i := 0; i < mercuryStationAnchorAttempts && isFiniteFloat(inferior); i++ { date := mercuryRetrogradeToProgradeAroundInferior(inferior) if isFiniteFloat(date) && stationUTQueryAfterOrEqual(date, jde) { return date } inferior = NextMercuryInferiorConjunction(eventUTNextQueryTT(inferior)) } return math.NaN() } func LastMercuryProgradeToRetrograde(jde float64) float64 { if !isFiniteFloat(jde) { return math.NaN() } inferior := NextMercuryInferiorConjunction(jde) for i := 0; i < mercuryStationAnchorAttempts && isFiniteFloat(inferior); i++ { date := mercuryProgradeToRetrogradeAroundInferior(inferior) if isFiniteFloat(date) && stationUTQueryBeforeOrEqual(date, jde) { return date } inferior = LastMercuryInferiorConjunction(eventUTLastQueryTT(inferior)) } return math.NaN() } func LastMercuryRetrogradeToPrograde(jde float64) float64 { if !isFiniteFloat(jde) { return math.NaN() } inferior := LastMercuryInferiorConjunction(jde) for i := 0; i < mercuryStationAnchorAttempts && isFiniteFloat(inferior); i++ { date := mercuryRetrogradeToProgradeAroundInferior(inferior) if isFiniteFloat(date) && stationUTQueryBeforeOrEqual(date, jde) { return date } inferior = LastMercuryInferiorConjunction(eventUTLastQueryTT(inferior)) } return math.NaN() } // earliestFiniteEventUT / latestFiniteEventUT 在候选中取最早/最晚的有限值(全部非有限时返回 NaN)。 func earliestFiniteEventUT(candidates ...float64) float64 { best := math.NaN() for _, candidate := range candidates { if !isFiniteFloat(candidate) { continue } if math.IsNaN(best) || candidate < best { best = candidate } } return best } func latestFiniteEventUT(candidates ...float64) float64 { best := math.NaN() for _, candidate := range candidates { if !isFiniteFloat(candidate) { continue } if math.IsNaN(best) || candidate > best { best = candidate } } return best } func nextMercuryRetrogradeFromTyped(jde float64) float64 { return earliestFiniteEventUT(NextMercuryProgradeToRetrograde(jde), NextMercuryRetrogradeToPrograde(jde)) } func NextMercuryRetrograde(jde float64) float64 { if !isFiniteFloat(jde) { return math.NaN() } motion := mercuryRADerivative(jde, mercuryStationDerivativeStepDay) if motion > mercuryStationMotionTolerance { p2r := NextMercuryProgradeToRetrograde(jde) if isFiniteFloat(p2r) && !mercuryStationBetween(jde, TD2UT(p2r, true)) { return p2r } best := earliestFiniteEventUT(p2r, NextMercuryRetrogradeToPrograde(jde)) if !stationUTQueryAfterOrEqual(best, jde) { return math.NaN() } return best } if motion < -mercuryStationMotionTolerance { r2p := NextMercuryRetrogradeToPrograde(jde) if isFiniteFloat(r2p) && !mercuryStationBetween(jde, TD2UT(r2p, true)) { return r2p } best := earliestFiniteEventUT(NextMercuryProgradeToRetrograde(jde), r2p) if !stationUTQueryAfterOrEqual(best, jde) { return math.NaN() } return best } return nextMercuryRetrogradeFromTyped(jde) } func lastMercuryRetrogradeFromTyped(jde float64) float64 { return latestFiniteEventUT(LastMercuryProgradeToRetrograde(jde), LastMercuryRetrogradeToPrograde(jde)) } func LastMercuryRetrograde(jde float64) float64 { if !isFiniteFloat(jde) { return math.NaN() } motion := mercuryRADerivative(jde, mercuryStationDerivativeStepDay) if motion > mercuryStationMotionTolerance { r2p := LastMercuryRetrogradeToPrograde(jde) if isFiniteFloat(r2p) && !mercuryStationBetween(TD2UT(r2p, true), jde) { return r2p } best := latestFiniteEventUT(LastMercuryProgradeToRetrograde(jde), r2p) if !stationUTQueryBeforeOrEqual(best, jde) { return math.NaN() } return best } if motion < -mercuryStationMotionTolerance { p2r := LastMercuryProgradeToRetrograde(jde) if isFiniteFloat(p2r) && !mercuryStationBetween(TD2UT(p2r, true), jde) { return p2r } best := latestFiniteEventUT(p2r, LastMercuryRetrogradeToPrograde(jde)) if !stationUTQueryBeforeOrEqual(best, jde) { return math.NaN() } return best } return lastMercuryRetrogradeFromTyped(jde) } func LastMercuryRetrogradeStrict(jde float64) float64 { return LastMercuryRetrograde(jde) } func NextMercuryRetrogradeStrict(jde float64) float64 { return NextMercuryRetrograde(jde) } func MercurySunElongation(jde float64) float64 { lo1, bo1 := MercuryApparentLoBo(jde) lo2 := HSunApparentLo(jde) bo2 := HSunTrueBo(jde) return StarAngularSeparation(lo1, bo1, lo2, bo2) } func mercurySunElongationN(jde float64, n int) float64 { lo1, bo1 := MercuryApparentLoBoN(jde, n) lo2 := HSunApparentLoN(jde, n) bo2 := HSunTrueBoN(jde, n) return StarAngularSeparation(lo1, bo1, lo2, bo2) } // mercuryGreatestElongationInWindow 求窗口内大距:目标是公开 MercurySunElongation 的极大(视距角), // 而不是忽略光行差/视位置修正的真距角,否则返回的时刻不是调用方能量到的那个极值。 // 窗口两端是世界时,目标函数收力学时,因此逐次换算。 func mercuryGreatestElongationInWindow(start, end float64) float64 { return maximizeInWindow(start, end, 2.0, func(utJD float64) float64 { return mercurySunElongationN(TD2UT(utJD, true), mercuryEventSearchN) }, func(utJD float64) float64 { return MercurySunElongation(TD2UT(utJD, true)) }) } func mercuryEastElongationWindowEndingAt(inferior float64) (float64, float64) { lastSuperior := LastMercurySuperiorConjunction(eventUTLastQueryTT(inferior)) return lastSuperior + innerEventWindowPadding, inferior - innerEventWindowPadding } func mercuryWestElongationWindowEndingAt(superior float64) (float64, float64) { lastInferior := LastMercuryInferiorConjunction(eventUTLastQueryTT(superior)) return lastInferior + innerEventWindowPadding, superior - innerEventWindowPadding } func mercuryEastElongationWindowContaining(jde float64) (float64, float64) { nextInferior := NextMercuryInferiorConjunction(jde) start, end := mercuryEastElongationWindowEndingAt(nextInferior) if eventUTQueryBeforeOrEqual(start, jde) { return start, end } currentInferior := LastMercuryInferiorConjunction(jde) return mercuryEastElongationWindowEndingAt(currentInferior) } func mercuryWestElongationWindowContaining(jde float64) (float64, float64) { nextSuperior := NextMercurySuperiorConjunction(jde) start, end := mercuryWestElongationWindowEndingAt(nextSuperior) if eventUTQueryBeforeOrEqual(start, jde) { return start, end } currentSuperior := LastMercurySuperiorConjunction(jde) return mercuryWestElongationWindowEndingAt(currentSuperior) } func nextMercuryGreatestElongationTyped(jde float64, east bool) float64 { if !isFiniteFloat(jde) { return math.NaN() } if east { start, windowEnd := mercuryEastElongationWindowContaining(jde) for i := 0; i < eventDirectionalSearchIterations; i++ { date := mercuryGreatestElongationInWindow(start, windowEnd) if !isFiniteFloat(date) { return math.NaN() } if eventUTQueryAfterOrEqual(date, jde) { return date } nextInferior := NextMercuryInferiorConjunction(eventUTNextQueryTT(windowEnd)) start, windowEnd = mercuryEastElongationWindowEndingAt(nextInferior) } return math.NaN() } start, windowEnd := mercuryWestElongationWindowContaining(jde) for i := 0; i < eventDirectionalSearchIterations; i++ { date := mercuryGreatestElongationInWindow(start, windowEnd) if !isFiniteFloat(date) { return math.NaN() } if eventUTQueryAfterOrEqual(date, jde) { return date } nextSuperior := NextMercurySuperiorConjunction(eventUTNextQueryTT(windowEnd)) start, windowEnd = mercuryWestElongationWindowEndingAt(nextSuperior) } return math.NaN() } func lastMercuryGreatestElongationTyped(jde float64, east bool) float64 { if !isFiniteFloat(jde) { return math.NaN() } if east { start, windowEnd := mercuryEastElongationWindowContaining(jde) for i := 0; i < eventDirectionalSearchIterations; i++ { date := mercuryGreatestElongationInWindow(start, windowEnd) if !isFiniteFloat(date) { return math.NaN() } if eventUTQueryBeforeOrEqual(date, jde) { return date } prevInferior := LastMercuryInferiorConjunction(eventUTLastQueryTT(start)) start, windowEnd = mercuryEastElongationWindowEndingAt(prevInferior) } return math.NaN() } start, windowEnd := mercuryWestElongationWindowContaining(jde) for i := 0; i < eventDirectionalSearchIterations; i++ { date := mercuryGreatestElongationInWindow(start, windowEnd) if !isFiniteFloat(date) { return math.NaN() } if eventUTQueryBeforeOrEqual(date, jde) { return date } prevSuperior := LastMercurySuperiorConjunction(eventUTLastQueryTT(start)) start, windowEnd = mercuryWestElongationWindowEndingAt(prevSuperior) } return math.NaN() } // mercuryElongationWindowAt 返回包含 jde 的该侧大距窗口;查询落在合的 4 秒内边距里时两侧都不包含。 func mercuryElongationWindowAt(jde float64, east bool) (float64, float64, bool) { if east { start, end := mercuryEastElongationWindowEndingAt(NextMercuryInferiorConjunction(jde)) return start, end, eventUTQueryBeforeOrEqual(start, jde) && eventUTQueryAfterOrEqual(end, jde) } start, end := mercuryWestElongationWindowEndingAt(NextMercurySuperiorConjunction(jde)) return start, end, eventUTQueryBeforeOrEqual(start, jde) && eventUTQueryAfterOrEqual(end, jde) } // 无东西侧参数的 Next/Last 先只看查询所在窗口那一侧:该侧大距未过就是答案,已过则另一侧紧接着的 // 下一个才是答案,因此只有在查询贴住合(两侧窗口都不含)时才退化为两侧都算。赤经差在合附近会提前 // 变号,不能用它定窗口。 func NextMercuryGreatestElongation(jde float64) float64 { if !isFiniteFloat(jde) { return math.NaN() } for _, east := range [2]bool{true, false} { start, end, ok := mercuryElongationWindowAt(jde, east) if !ok { continue } if date := mercuryGreatestElongationInWindow(start, end); isFiniteFloat(date) && eventUTQueryAfterOrEqual(date, jde) { return date } if date := nextMercuryGreatestElongationTyped(jde, !east); isFiniteFloat(date) && eventUTQueryAfterOrEqual(date, jde) { return date } break } return earliestFiniteEventUT(nextMercuryGreatestElongationTyped(jde, true), nextMercuryGreatestElongationTyped(jde, false)) } func LastMercuryGreatestElongation(jde float64) float64 { if !isFiniteFloat(jde) { return math.NaN() } for _, east := range [2]bool{true, false} { start, end, ok := mercuryElongationWindowAt(jde, east) if !ok { continue } if date := mercuryGreatestElongationInWindow(start, end); isFiniteFloat(date) && eventUTQueryBeforeOrEqual(date, jde) { return date } if date := lastMercuryGreatestElongationTyped(jde, !east); isFiniteFloat(date) && eventUTQueryBeforeOrEqual(date, jde) { return date } break } return latestFiniteEventUT(lastMercuryGreatestElongationTyped(jde, true), lastMercuryGreatestElongationTyped(jde, false)) } func LastMercuryInferiorConjunctionInclusive(jde float64) float64 { return inclusiveLastSimpleEvent(jde, LastMercuryInferiorConjunction, NextMercuryInferiorConjunction) } func NextMercuryInferiorConjunctionInclusive(jde float64) float64 { return inclusiveNextSimpleEvent(jde, LastMercuryInferiorConjunction, NextMercuryInferiorConjunction) } func LastMercurySuperiorConjunctionInclusive(jde float64) float64 { return inclusiveLastSimpleEvent(jde, LastMercurySuperiorConjunction, NextMercurySuperiorConjunction) } func NextMercurySuperiorConjunctionInclusive(jde float64) float64 { return inclusiveNextSimpleEvent(jde, LastMercurySuperiorConjunction, NextMercurySuperiorConjunction) } func LastMercuryRetrogradeInclusive(jde float64) float64 { return inclusiveLastSimpleEvent(jde, LastMercuryRetrograde, NextMercuryRetrograde) } func NextMercuryRetrogradeInclusive(jde float64) float64 { return inclusiveNextSimpleEvent(jde, LastMercuryRetrograde, NextMercuryRetrograde) } func LastMercuryProgradeToRetrogradeInclusive(jde float64) float64 { return inclusiveLastSimpleEvent(jde, LastMercuryProgradeToRetrograde, NextMercuryProgradeToRetrograde) } func NextMercuryProgradeToRetrogradeInclusive(jde float64) float64 { return inclusiveNextSimpleEvent(jde, LastMercuryProgradeToRetrograde, NextMercuryProgradeToRetrograde) } func LastMercuryRetrogradeToProgradeInclusive(jde float64) float64 { return inclusiveLastSimpleEvent(jde, LastMercuryRetrogradeToPrograde, NextMercuryRetrogradeToPrograde) } func NextMercuryRetrogradeToProgradeInclusive(jde float64) float64 { return inclusiveNextSimpleEvent(jde, LastMercuryRetrogradeToPrograde, NextMercuryRetrogradeToPrograde) } func LastMercuryGreatestElongationInclusive(jde float64) float64 { return inclusiveLastSimpleEvent(jde, LastMercuryGreatestElongation, NextMercuryGreatestElongation) } func NextMercuryGreatestElongationInclusive(jde float64) float64 { return inclusiveNextSimpleEvent(jde, LastMercuryGreatestElongation, NextMercuryGreatestElongation) } func LastMercuryGreatestElongationEastInclusive(jde float64) float64 { return inclusiveLastSimpleEvent(jde, LastMercuryGreatestElongationEast, NextMercuryGreatestElongationEast) } func NextMercuryGreatestElongationEastInclusive(jde float64) float64 { return inclusiveNextSimpleEvent(jde, LastMercuryGreatestElongationEast, NextMercuryGreatestElongationEast) } func LastMercuryGreatestElongationWestInclusive(jde float64) float64 { return inclusiveLastSimpleEvent(jde, LastMercuryGreatestElongationWest, NextMercuryGreatestElongationWest) } func NextMercuryGreatestElongationWestInclusive(jde float64) float64 { return inclusiveNextSimpleEvent(jde, LastMercuryGreatestElongationWest, NextMercuryGreatestElongationWest) } func NextMercuryGreatestElongationEast(jde float64) float64 { return nextMercuryGreatestElongationTyped(jde, true) } func NextMercuryGreatestElongationWest(jde float64) float64 { return nextMercuryGreatestElongationTyped(jde, false) } func LastMercuryGreatestElongationEast(jde float64) float64 { return lastMercuryGreatestElongationTyped(jde, true) } func LastMercuryGreatestElongationWest(jde float64) float64 { return lastMercuryGreatestElongationTyped(jde, false) }