diff --git a/basic/api_default.go b/basic/api_default.go index 8048031..83edbe5 100644 --- a/basic/api_default.go +++ b/basic/api_default.go @@ -14,7 +14,7 @@ func ClosestJupiterGalileanPhenomenonEvent(jd float64, satellite int, phenomenon return ClosestJupiterGalileanPhenomenonEventWithInstance(jd, satellite, phenomenonType, nil) } -// ClosestMercuryTransit 最近一次地心水星凌日 / closest geocentric transit of Mercury. +// ClosestMercuryTransit 最近一次全球水星凌日 / closest global transit of Mercury. func ClosestMercuryTransit(jde float64) PlanetTransitResult { return ClosestMercuryTransitWithInstance(jde, nil) } @@ -34,7 +34,7 @@ func ClosestMoonPlanetConjunction(jde float64, planet MoonPlanetConjunctionPlane return ClosestMoonPlanetConjunctionWithInstance(jde, planet, nil) } -// ClosestVenusTransit 最近一次地心金星凌日 / closest geocentric transit of Venus. +// ClosestVenusTransit 最近一次全球金星凌日 / closest global transit of Venus. func ClosestVenusTransit(jde float64) PlanetTransitResult { return ClosestVenusTransitWithInstance(jde, nil) } @@ -236,7 +236,7 @@ func LastJupiterGalileanPhenomenonEvent(jd float64, satellite int, phenomenonTyp return LastJupiterGalileanPhenomenonEventWithInstance(jd, satellite, phenomenonType, nil) } -// LastMercuryTransit 上一次地心水星凌日 / previous geocentric transit of Mercury. +// LastMercuryTransit 上一次全球水星凌日 / previous global transit of Mercury. func LastMercuryTransit(jde float64) PlanetTransitResult { return LastMercuryTransitWithInstance(jde, nil) } @@ -256,7 +256,7 @@ func LastMoonPlanetConjunction(jde float64, planet MoonPlanetConjunctionPlanet) return LastMoonPlanetConjunctionWithInstance(jde, planet, nil) } -// LastVenusTransit 上一次地心金星凌日 / previous geocentric transit of Venus. +// LastVenusTransit 上一次全球金星凌日 / previous global transit of Venus. func LastVenusTransit(jde float64) PlanetTransitResult { return LastVenusTransitWithInstance(jde, nil) } // LocalSolarEclipseDiagramIAUSingleK 站心日食视圆图几何（IAU Single-K） / local solar eclipse disk diagram, IAU Single-K. @@ -560,7 +560,7 @@ func NextJupiterGalileanPhenomenonEvent(jd float64, satellite int, phenomenonTyp return NextJupiterGalileanPhenomenonEventWithInstance(jd, satellite, phenomenonType, nil) } -// NextMercuryTransit 下一次地心水星凌日 / next geocentric transit of Mercury. +// NextMercuryTransit 下一次全球水星凌日 / next global transit of Mercury. func NextMercuryTransit(jde float64) PlanetTransitResult { return NextMercuryTransitWithInstance(jde, nil) } @@ -580,7 +580,7 @@ func NextMoonPlanetConjunction(jde float64, planet MoonPlanetConjunctionPlanet) return NextMoonPlanetConjunctionWithInstance(jde, planet, nil) } -// NextVenusTransit 下一次地心金星凌日 / next geocentric transit of Venus. +// NextVenusTransit 下一次全球金星凌日 / next global transit of Venus. func NextVenusTransit(jde float64) PlanetTransitResult { return NextVenusTransitWithInstance(jde, nil) } // PlanetOccultationDiagram 已求解固定地点行星月掩的月心切平面轨迹，事件无效时不含帧 / tangent-plane track of a solved fixed-site planetary occultation. @@ -596,6 +596,16 @@ func PlanetOccultationFootprintsAt(at time.Time, planet OccultationPlanet) (Plan return PlanetOccultationFootprintsAtWithInstance(at, planet, nil) } +// PlanetTransitVisibilityColumns 按 2° 经度、1° 纬度探测站心可见区并细化边界 / sample topocentric visibility at 2-degree longitude and 1-degree latitude spacing and refine boundaries. +func PlanetTransitVisibilityColumns(transit PlanetTransitResult) []PlanetTransitVisibilityColumn { + return PlanetTransitVisibilityColumnsWithInstance(transit, nil) +} + +// PlanetTransitVisibilityPolygons 可见区的经度条带多边形 / longitude-strip polygons of the visible region. +func PlanetTransitVisibilityPolygons(transit PlanetTransitResult) [][]PlanetTransitIsochronePoint { + return PlanetTransitVisibilityPolygonsWithInstance(transit, nil) +} + // SaturnAscendingNodeN 土星升交点黄经（n 变体） / n-parameterised ascending node longitude of Saturn. func SaturnAscendingNodeN(jd float64, n int) float64 { return SaturnAscendingNodeNWithInstance(jd, n, nil) diff --git a/basic/planet_transit.go b/basic/planet_transit.go index 00040ef..9385351 100644 --- a/basic/planet_transit.go +++ b/basic/planet_transit.go @@ -7,20 +7,17 @@ import ( ) const ( - planetTransitMeanSolarMotionDegPerDay = 360.0 / 365.2422 - planetTransitTropicalYearDays = 365.2422 - planetTransitSeasonProbeStepDays = 20.0 - planetTransitSearchLimit = 2400 - planetTransitSearchEpsilonDays = 1.0 / 86400.0 - planetTransitGreatestWindowDays = 1.2 - planetTransitGreatestToleranceDays = 0.25 / 86400.0 - planetTransitContactStepDays = 0.02 - planetTransitContactSpanDays = 1.0 - planetTransitContactToleranceDays = 0.25 / 86400.0 - planetTransitCoarseN = 16 + planetTransitSearchLimit = 2400 + planetTransitSearchEpsilonDays = 1.0 / 86400.0 + planetTransitGreatestWindowDays = 1.2 + planetTransitGreatestToleranceDays = 0.25 / 86400.0 + planetTransitContactStepDays = 0.02 + planetTransitContactSpanDays = 1.0 + planetTransitContactToleranceDays = 0.25 / 86400.0 + planetTransitCoarseN = 16 ) -// PlanetTransitResult 表示一次地心行星凌日结果。 +// PlanetTransitResult 全球行星凌日与地心接触 / global planetary transit and geocentric contacts. // // Valid 为 false 时表示没有找到有效凌日。所有时刻字段均为 UT 儒略日。 // MinimumSeparationArcsec、SunSemidiameterArcsec、PlanetSemidiameterArcsec 的单位均为角秒。 @@ -45,6 +42,15 @@ type PlanetTransitResult struct { HasExternal bool HasInternal bool + + // WindowStart、WindowEnd 是站心求解窗口，不是接触时刻。 + WindowStart float64 + WindowEnd float64 + // Gamma 是凌甚时影轴到地心的有符号距离，以地球赤道半径为单位，北侧为正。 + Gamma float64 + // PenumbralRadius 是地心基本面上的半影半径，以地球赤道半径为单位。 + PenumbralRadius float64 + Series PlanetTransitSeries } type planetTransitConfig struct { @@ -53,7 +59,6 @@ type planetTransitConfig struct { synodicPeriodDays float64 anchorInferiorTT float64 - seasonWindowDays float64 latitudePrefilter float64 conjunctionStepDay float64 @@ -62,7 +67,6 @@ type planetTransitConfig struct { apparentRaDecN func(float64, int) (float64, float64) semidiameterN func(float64, int) float64 earthDistanceN func(float64, int) float64 - nodeN func(float64, int) float64 } type planetTransitState struct { @@ -82,7 +86,6 @@ func mercuryTransitConfig(instances ...*Instance) planetTransitConfig { planetIndex: 1, synodicPeriodDays: MERCURY_S_PERIOD, anchorInferiorTT: UTC2TT(JDCalc(2019, 11, 11+(15+21.0/60+40.0/3600)/24)), - seasonWindowDays: 12, latitudePrefilter: 1.0, conjunctionStepDay: 0.00001, apparentLoN: func(jd float64, n int) float64 { return MercuryApparentLoNWithInstance(jd, n, inst) }, @@ -90,7 +93,6 @@ func mercuryTransitConfig(instances ...*Instance) planetTransitConfig { apparentRaDecN: func(jde float64, n int) (float64, float64) { return MercuryApparentRaDecNWithInstance(jde, n, inst) }, semidiameterN: func(jde float64, n int) float64 { return MercurySemidiameterNWithInstance(jde, n, inst) }, earthDistanceN: func(jd float64, n int) float64 { return EarthMercuryAwayNWithInstance(jd, n, inst) }, - nodeN: func(jd float64, n int) float64 { return MercuryAscendingNodeNWithInstance(jd, n, inst) }, } } @@ -101,7 +103,6 @@ func venusTransitConfig(instances ...*Instance) planetTransitConfig { planetIndex: 2, synodicPeriodDays: VENUS_S_PERIOD, anchorInferiorTT: UTC2TT(JDCalc(2012, 6, 6+(1+29.0/60)/24)), - seasonWindowDays: 8, latitudePrefilter: 0.8, conjunctionStepDay: 0.00001, apparentLoN: func(jd float64, n int) float64 { return VenusApparentLoNWithInstance(jd, n, inst) }, @@ -109,7 +110,6 @@ func venusTransitConfig(instances ...*Instance) planetTransitConfig { apparentRaDecN: func(jd float64, n int) (float64, float64) { return VenusApparentRaDecNWithInstance(jd, n, inst) }, semidiameterN: func(jde float64, n int) float64 { return VenusSemidiameterNWithInstance(jde, n, inst) }, earthDistanceN: func(jd float64, n int) float64 { return EarthVenusAwayNWithInstance(jd, n, inst) }, - nodeN: func(jd float64, n int) float64 { return VenusAscendingNodeNWithInstance(jd, n, inst) }, } } @@ -174,18 +174,24 @@ func searchPlanetTransit(jd float64, cfg planetTransitConfig, direction int, inc } targetTT := UTC2TT(jd) - // 节点季可能略早于查询时刻而凌甚仍在其后，从一个季窗之前起探，方向由凌甚时刻再判。 - probeTT := targetTT - float64(direction)*cfg.seasonWindowDays + if !isFiniteFloat(targetTT) { + return PlanetTransitResult{}, false + } + cycle := math.Floor((targetTT-cfg.anchorInferiorTT)/cfg.synodicPeriodDays) - 1 + if direction < 0 { + cycle += 3 + } for i := 0; i < planetTransitSearchLimit; i++ { - seasonTT, ok := nextPlanetTransitSeasonTT(probeTT, cfg, direction, inst) - if !ok { - return PlanetTransitResult{}, false - } - - seedTT := nearestPlanetTransitInferiorSeedTT(seasonTT, cfg) - if math.Abs(seedTT-seasonTT) <= cfg.seasonWindowDays { - conjunctionTT := refinePlanetTransitInferiorConjunctionTT(seedTT, cfg, inst) - if math.Abs(conjunctionTT-seasonTT) <= cfg.seasonWindowDays+1 && isPotentialPlanetTransit(conjunctionTT, cfg, inst) { + seedTT := cfg.anchorInferiorTT + cycle*cfg.synodicPeriodDays + if planetTransitPredictedSeason(seedTT, cfg, inst) { + // 单个会合周算不出粗时刻（例如候选落在有限覆盖星历之外）时跳过该周，后面的周期仍可能落在覆盖内。 + coarseTT := planetTransitCoarseConjunction(seedTT, cfg, inst) + if !isFiniteFloat(coarseTT) { + cycle += float64(direction) + continue + } + if isPotentialPlanetTransit(coarseTT, cfg, inst) { + conjunctionTT := refinePlanetTransitInferiorConjunctionTT(coarseTT, cfg, inst) resultTT, ok := planetTransitAtInferiorConjunctionTT(conjunctionTT, cfg, inst) if ok && planetTransitMatchesDirection(resultTT.Greatest, targetTT, direction, includeCurrent) { return planetTransitResultTTToUT(resultTT), true @@ -193,102 +199,12 @@ func searchPlanetTransit(jd float64, cfg planetTransitConfig, direction int, inc } } - probeTT = seasonTT + float64(direction)*planetTransitSeasonProbeStepDays + cycle += float64(direction) } return PlanetTransitResult{}, false } -func nextPlanetTransitSeasonTT(jdTT float64, cfg planetTransitConfig, direction int, instances ...*Instance) (float64, bool) { - inst := instanceValue(instances) - - best := math.NaN() - for nodeOffset := 0; nodeOffset <= 1; nodeOffset++ { - candidate := estimatePlanetTransitSeasonTT(jdTT, cfg, nodeOffset, direction, inst) - candidate = refinePlanetTransitSeasonTT(candidate, cfg, nodeOffset, inst) - if !isFiniteFloat(candidate) { - continue - } - matchedDirection := false - for i := 0; i < eventDirectionalSearchIterations; i++ { - if planetTransitMatchesDirection(candidate, jdTT, direction, false) { - matchedDirection = true - break - } - candidate += float64(direction) * planetTransitTropicalYearDays - candidate = refinePlanetTransitSeasonTT(candidate, cfg, nodeOffset, inst) - if !isFiniteFloat(candidate) { - break - } - } - if !matchedDirection { - continue - } - if !isFiniteFloat(best) || math.Abs(candidate-jdTT) < math.Abs(best-jdTT) { - best = candidate - } - } - if !isFiniteFloat(best) { - return 0, false - } - return best, true -} - -func estimatePlanetTransitSeasonTT(jdTT float64, cfg planetTransitConfig, nodeOffset int, direction int, instances ...*Instance) float64 { - inst := instanceValue(instances) - - sunLongitude := HSunApparentLoNWithInstance(jdTT, planetTransitCoarseN, inst) - nodeLongitude := planetTransitNodeLongitude(jdTT, cfg, nodeOffset, planetTransitCoarseN) - if direction > 0 { - delta := Limit360(nodeLongitude - sunLongitude) - if delta <= planetTransitSearchEpsilonDays { - delta += 360 - } - return jdTT + delta/planetTransitMeanSolarMotionDegPerDay - } - delta := Limit360(sunLongitude - nodeLongitude) - if delta <= planetTransitSearchEpsilonDays { - delta += 360 - } - return jdTT - delta/planetTransitMeanSolarMotionDegPerDay -} - -func refinePlanetTransitSeasonTT(seedTT float64, cfg planetTransitConfig, nodeOffset int, instances ...*Instance) float64 { - inst := instanceValue(instances) - - current := seedTT - for i := 0; i < 8; i++ { - prev := current - value := planetTransitSunNodeLongitudeDelta(prev, cfg, nodeOffset, inst) - slope := (planetTransitSunNodeLongitudeDelta(prev+0.5, cfg, nodeOffset, inst) - - planetTransitSunNodeLongitudeDelta(prev-0.5, cfg, nodeOffset, inst)) / 1.0 - if slope == 0 || !isFiniteFloat(slope) { - break - } - current = prev - value/slope - if math.Abs(current-prev) <= 0.00001 { - break - } - } - return current -} - -func planetTransitSunNodeLongitudeDelta(jdTT float64, cfg planetTransitConfig, nodeOffset int, instances ...*Instance) float64 { - inst := instanceValue(instances) - - return planetTransitAngleDelta(HSunApparentLoNWithInstance(jdTT, planetTransitCoarseN, inst) - - planetTransitNodeLongitude(jdTT, cfg, nodeOffset, planetTransitCoarseN)) -} - -func planetTransitNodeLongitude(jdTT float64, cfg planetTransitConfig, nodeOffset int, n int) float64 { - return Limit360(cfg.nodeN(jdTT, n) + float64(nodeOffset)*180) -} - -func nearestPlanetTransitInferiorSeedTT(seasonTT float64, cfg planetTransitConfig) float64 { - k := math.Round((seasonTT - cfg.anchorInferiorTT) / cfg.synodicPeriodDays) - return cfg.anchorInferiorTT + k*cfg.synodicPeriodDays -} - func refinePlanetTransitInferiorConjunctionTT(seedTT float64, cfg planetTransitConfig, instances ...*Instance) float64 { inst := instanceValue(instances) @@ -342,7 +258,11 @@ func planetTransitAtInferiorConjunctionTT(conjunctionTT float64, cfg planetTrans greatestTT := greatestPlanetTransitTT(conjunctionTT, cfg, inst) greatestState := planetTransitStateAt(greatestTT, cfg, -1, inst) - if !isFiniteFloat(greatestState.externalContactMetric) || greatestState.externalContactMetric > 0 { + if !isFiniteFloat(greatestState.externalContactMetric) { + return PlanetTransitResult{}, false + } + gamma, penumbra, surfaceSeed, possible := planetTransitShadowGeometry(greatestTT, cfg, inst) + if !possible { return PlanetTransitResult{}, false } @@ -353,8 +273,19 @@ func planetTransitAtInferiorConjunctionTT(conjunctionTT float64, cfg planetTrans MinimumSeparationArcsec: greatestState.separationArcsec, SunSemidiameterArcsec: greatestState.sunSemidiameter, PlanetSemidiameterArcsec: greatestState.planetSemidiameter, - HasExternal: true, + HasExternal: greatestState.externalContactMetric <= 0, HasInternal: greatestState.internalContactMetric <= 0, + WindowStart: greatestTT - 0.25, + WindowEnd: greatestTT + 0.25, + Gamma: gamma, + PenumbralRadius: penumbra, + Series: planetTransitSeriesAt(greatestTT, cfg), + } + if !result.HasExternal { + if !planetTransitSurfaceOverlap(planetTransitResultTTToUT(result), surfaceSeed, inst) { + return PlanetTransitResult{}, false + } + return result, true } externalIngress, ok := refinePlanetTransitContactTT(greatestTT, cfg, -1, false, inst) @@ -367,6 +298,8 @@ func planetTransitAtInferiorConjunctionTT(conjunctionTT float64, cfg planetTrans } result.ExternalIngress = externalIngress result.ExternalEgress = externalEgress + result.WindowStart = externalIngress - planetTransitSnapshotMarginDays + result.WindowEnd = externalEgress + planetTransitSnapshotMarginDays if result.HasInternal { internalIngress, ok := refinePlanetTransitContactTT(greatestTT, cfg, -1, true, inst) @@ -527,8 +460,18 @@ func bisectPlanetTransitContactTT(leftJD, leftValue, rightJD, rightValue float64 func planetTransitResultTTToUT(result PlanetTransitResult) PlanetTransitResult { result.Greatest = TT2UTC(result.Greatest) - result.ExternalIngress = TT2UTC(result.ExternalIngress) - result.ExternalEgress = TT2UTC(result.ExternalEgress) + if result.ExternalIngress != 0 { + result.ExternalIngress = TT2UTC(result.ExternalIngress) + } + if result.ExternalEgress != 0 { + result.ExternalEgress = TT2UTC(result.ExternalEgress) + } + if result.WindowStart != 0 { + result.WindowStart = TT2UTC(result.WindowStart) + } + if result.WindowEnd != 0 { + result.WindowEnd = TT2UTC(result.WindowEnd) + } if result.InternalIngress != 0 { result.InternalIngress = TT2UTC(result.InternalIngress) } diff --git a/basic/planet_transit_geometry.go b/basic/planet_transit_geometry.go new file mode 100644 index 0000000..7fb160c --- /dev/null +++ b/basic/planet_transit_geometry.go @@ -0,0 +1,67 @@ +package basic + +import "math" + +const ( + planetTransitEarthRadiusKM = 6378.137 + planetTransitEarthEccentricitySquared = 0.00669437999014 +) + +func planetTransitShadowGeometry(tt float64, cfg planetTransitConfig, inst *Instance) (float64, float64, Observer, bool) { + sunRA, sunDec := HSunApparentRaDecNWithInstance(tt, -1, inst) + planetRA, planetDec := cfg.apparentRaDecN(tt, -1) + sunDistance := EarthAwayNWithInstance(tt, -1, inst) * angularDiameterAstronomicalUnitKM + planetDistance := cfg.earthDistanceN(tt, -1) * angularDiameterAstronomicalUnitKM + sun := eclipticCartesian(sunRA, sunDec, sunDistance) + body := eclipticCartesian(planetRA, planetDec, planetDistance) + axis, length := pn(Vector3{body[0] - sun[0], body[1] - sun[1], body[2] - sun[2]}) + projection := -(body[0]*axis[0] + body[1]*axis[1] + body[2]*axis[2]) + impact := Vector3{body[0] + projection*axis[0], body[1] + projection*axis[1], body[2] + projection*axis[2]} + _, distance := pn(impact) + planetRadius := venusEquatorialRadiusKM + if cfg.planetIndex == 1 { + planetRadius = mercuryEquatorialRadiusKM + } + sine := (sunEquatorialRadiusKM + planetRadius) / length + cosine := math.Sqrt(1 - sine*sine) + penumbra := (planetRadius + projection*sine) / cosine / planetTransitEarthRadiusKM + gamma := math.Copysign(distance/planetTransitEarthRadiusKM, impact[2]) + sidereal := ApparentSiderealTime(TT2UT1(tt)) * 15 + seed := Observer{ + Longitude: normalizeLongitude180(math.Atan2(impact[1], impact[0])*180/math.Pi - sidereal), + Latitude: math.Atan2(impact[2], math.Hypot(impact[0], impact[1])*(1-planetTransitEarthEccentricitySquared)) * 180 / math.Pi, + } + possible := isFiniteFloat(gamma) && isFiniteFloat(penumbra) && projection > 0 && planetDistance < sunDistance && math.Abs(gamma) <= penumbra+1/cosine + return gamma, penumbra, seed, possible +} + +func planetTransitSurfaceOverlap(transit PlanetTransitResult, seed Observer, inst *Instance) bool { + snapshot, ok := buildPlanetTransitSnapshot(transit, inst) + if !ok { + return false + } + best := math.Inf(1) + for step := 4.0; step > 0.0001; step /= 2 { + center := seed + for _, latitudeStep := range []float64{-step, 0, step} { + for _, longitudeStep := range []float64{-step, 0, step} { + observer := Observer{Longitude: normalizeLongitude180(center.Longitude + longitudeStep), Latitude: math.Max(-90, math.Min(90, center.Latitude+latitudeStep))} + tt, valid := snapshot.localGreatestTT(UTC2TT(transit.Greatest), observer) + state, stateOK := snapshot.stateAt(tt, observer) + if !valid || !stateOK { + continue + } + if state.external < 0 { + local := snapshot.local(transit, observer) + if local.Valid && local.Visible { + return true + } + } + if state.external < best { + best, seed = state.external, observer + } + } + } + } + return false +} diff --git a/basic/planet_transit_global_test.go b/basic/planet_transit_global_test.go new file mode 100644 index 0000000..d30eb09 --- /dev/null +++ b/basic/planet_transit_global_test.go @@ -0,0 +1,170 @@ +package basic + +import ( + "math" + "testing" + "time" + + "b612.me/astro/internal/ephem" +) + +func TestPlanetTransitGlobalGrazes(tester *testing.T) { + for _, model := range []ephem.Model{ephem.Semianalytic, ephem.Legacy} { + inst := ephem.NewModelInstance(model) + for _, year := range []int{-540, 2854} { + transit := NextVenusTransitWithInstance(JDCalc(year, 1, 1), inst) + if !transit.Valid || JD2DateByZone(transit.Greatest, time.UTC, false).Year() != year { + tester.Fatalf("model=%d year=%d: missed global event: %+v", model, year, transit) + } + if transit.HasExternal || transit.HasInternal || transit.ExternalIngress != 0 || transit.ExternalEgress != 0 { + tester.Fatalf("global graze fabricated geocentric contacts: %+v", transit) + } + if math.Abs(transit.Gamma) <= transit.PenumbralRadius || math.Abs(transit.Gamma) >= transit.PenumbralRadius+1.01 { + tester.Fatalf("unexpected shadow geometry: %+v", transit) + } + observer := Observer{Longitude: -135, Latitude: -60} + local := PlanetTransitLocalWithInstance(transit, observer, inst) + if !local.Valid || !local.Visible || local.HasInternal || !local.ExternalIngress.Valid || !local.ExternalEgress.Valid { + tester.Fatalf("model=%d year=%d: missing visible local graze: %+v", model, year, local) + } + if year == 2854 { + if other := PlanetTransitLocalWithInstance(transit, Observer{Latitude: 45}, inst); other.Valid { + tester.Fatal("northern site must miss the 2854 event") + } + } + } + } +} + +func TestPlanetTransitSeriesNumbers(tester *testing.T) { + for _, sample := range []struct{ planet, year, number, period int }{ + {2, 1631, 6, 243}, {2, 1639, 4, 243}, {2, 2004, 3, 243}, {2, 2012, 5, 243}, {2, 2611, 4, 243}, {2, 2854, 4, 243}, + {1, 1957, 9, 46}, {1, 2003, 9, 46}, {1, 2006, 8, 46}, {1, 2016, 7, 46}, {1, 2019, 6, 46}, {1, 2174, 11, 46}, + } { + query := JDCalc(sample.year, 1, 1) + result := NextVenusTransit(query) + if sample.planet == 1 { + result = NextMercuryTransit(query) + } + if result.Series.Number != sample.number || result.Series.PeriodYears != sample.period { + tester.Fatalf("planet=%d year=%d: series=%+v, want %d/%d", sample.planet, sample.year, result.Series, sample.number, sample.period) + } + } +} + +func TestPlanetTransitSeriesUnassigned(tester *testing.T) { + for _, sample := range []struct { + cfg planetTransitConfig + year int + }{ + {mercuryTransitConfig(), 1500}, {mercuryTransitConfig(), 2400}, + {venusTransitConfig(), -2100}, {venusTransitConfig(), 4100}, + } { + series := planetTransitSeriesAt(JDCalc(sample.year, 1, 1), sample.cfg) + if series.Number != 0 || series.PeriodYears == 0 { + tester.Fatalf("unassigned series: %+v", series) + } + } +} + +func TestPlanetTransitNearInternalEdge(tester *testing.T) { + for _, year := range []int{1631, 2611} { + transit := NextVenusTransit(JDCalc(year, 1, 1)) + if !transit.Valid || !transit.HasExternal || !transit.HasInternal { + tester.Fatalf("year=%d: missing geocentric contacts", year) + } + full, partial := false, false + for latitude := -75; latitude <= 75; latitude += 30 { + for longitude := -180; longitude < 180; longitude += 45 { + local := PlanetTransitLocal(transit, Observer{Longitude: float64(longitude), Latitude: float64(latitude)}) + if !local.Valid || !local.Visible { + continue + } + if local.HasInternal { + full = true + } else { + partial = true + } + } + } + if !full || !partial { + tester.Fatalf("year=%d: full=%v partial=%v", year, full, partial) + } + } +} + +func TestPlanetTransitInputBounds(tester *testing.T) { + for _, query := range []float64{math.NaN(), math.Inf(1), math.Inf(-1), math.MaxFloat64} { + if result := NextVenusTransit(query); result.Valid { + tester.Fatal("non-finite query produced a transit") + } + } +} + +func TestPlanetTransitResultRebuild(tester *testing.T) { + greatest := time.Date(2854, 12, 14, 13, 0, 0, 0, time.UTC) + if PlanetTransitResultFromTimes(2, time.Time{}, time.Time{}, greatest, time.Time{}, time.Time{}).Valid { + tester.Fatal("contact rebuild accepted missing C1/C4") + } + for _, window := range [][2]time.Time{ + {{}, greatest.Add(time.Hour)}, {greatest, greatest.Add(time.Hour)}, {greatest.Add(-time.Hour), greatest}, + } { + if PlanetTransitResultFromWindow(2, window[0], greatest, window[1]).Valid { + tester.Fatalf("window rebuild accepted %v", window) + } + } + transit := NextVenusTransit(JDCalc(2854, 1, 1)) + rebuilt := PlanetTransitResultFromWindow(2, JD2DateByZone(transit.WindowStart, time.UTC, false), + JD2DateByZone(transit.Greatest, time.UTC, false), JD2DateByZone(transit.WindowEnd, time.UTC, false)) + observer := Observer{Longitude: -135, Latitude: -60} + want, got := PlanetTransitLocal(transit, observer), PlanetTransitLocal(rebuilt, observer) + if !got.Valid || math.Abs(got.ExternalIngress.JD-want.ExternalIngress.JD)*86400 > 0.01 || math.Abs(got.ExternalEgress.JD-want.ExternalEgress.JD)*86400 > 0.01 { + tester.Fatalf("window rebuild changed local contacts: got %+v want %+v", got, want) + } +} + +// boundedCoverageEngine 把内置星历限制在一个时间窗内，窗外的请求按越界报错。 +type boundedCoverageEngine struct { + ephem.Engine + start, end float64 +} + +func (bounded boundedCoverageEngine) within(jd float64) bool { + return jd >= bounded.start && jd <= bounded.end +} + +func (bounded boundedCoverageEngine) State(target, center ephem.Body, jdTDB float64) (pos, vel [3]float64, err error) { + if !bounded.within(jdTDB) { + return [3]float64{}, [3]float64{}, ephem.ErrOutOfRange + } + return bounded.Engine.State(target, center, jdTDB) +} + +func (bounded boundedCoverageEngine) Apparent(target ephem.Body, jdTDB float64) (pos [3]float64, lightTime, distance float64, err error) { + if !bounded.within(jdTDB) { + return [3]float64{}, 0, 0, ephem.ErrOutOfRange + } + return bounded.Engine.(ephem.ApparentEngine).Apparent(target, jdTDB) +} + +func (bounded boundedCoverageEngine) Coverage(ephem.Body) (start, end float64, ok bool) { + return bounded.start, bounded.end, true +} + +// 有限覆盖星历：搜索起点前的候选落在覆盖之外时，不得因此丢掉覆盖内的凌日。 +func TestPlanetTransitBoundedEphemerisCoverage(tester *testing.T) { + base := ephem.NewModelInstance(ephem.Semianalytic).Engine() + start := JDCalc(2019, 1, 1) + bounded := ephem.NewInstance(boundedCoverageEngine{Engine: base, start: start, end: JDCalc(2021, 1, 1)}) + query := start + 60 + transit := NextMercuryTransitWithInstance(query, bounded) + if !transit.Valid { + tester.Fatal("bounded ephemeris lost the 2019-11-11 transit inside its coverage") + } + if year := JD2DateByZone(transit.Greatest, time.UTC, false).Year(); year != 2019 { + tester.Fatalf("bounded ephemeris returned year %d, want 2019", year) + } + if transit.HasExternal != true || transit.HasInternal != true { + tester.Fatalf("bounded ephemeris returned an incomplete transit: %+v", transit) + } +} diff --git a/basic/planet_transit_local.go b/basic/planet_transit_local.go index c25c7d5..7a0380b 100644 --- a/basic/planet_transit_local.go +++ b/basic/planet_transit_local.go @@ -7,7 +7,7 @@ import ( "time" ) -// 站心凌日：以地心凌日为种子，对太阳与行星的视赤经赤纬分别做观测者视差修正后重解接触与凌甚。 +// 站心凌日：以全球凌日窗口为种子，分别修正太阳与行星的视差后重解接触与凌甚。 // 太阳视半径按站心距离缩放（影响最多约 0.04″）；行星视半径沿用地心值，站心修正低于 0.01″。 // 结果里的太阳高度一律为几何中心高度；可见性、日出日落与等时线裁剪按太阳中心视高度 0° 判定， // 即计入标准大气（1010 hPa、10 °C）下的地平折射，几何高度阈值约 -34′。 @@ -91,8 +91,7 @@ type PlanetTransitLocalResult struct { Sunset float64 } -// PlanetTransitLocal 返回地心凌日 transit 在观测者 observer 处的站心接触与凌甚。 -// PlanetTransitLocal returns the topocentric contacts and greatest instant of the geocentric transit seen by observer. +// PlanetTransitLocal 观测地的站心接触与凌甚 / topocentric contacts and greatest transit at the site. func PlanetTransitLocal(transit PlanetTransitResult, observer Observer) PlanetTransitLocalResult { return PlanetTransitLocalWithInstance(transit, observer, nil) } @@ -118,7 +117,9 @@ type planetTransitSnapshot struct { step float64 // 按列存储：太阳赤经、赤纬、距离，行星赤经、赤纬、距离，太阳与行星视半径，格林尼治视恒星时（度）。 // 赤经与恒星时已展开为连续值，避免插值跨越 0°/360°。 - columns [planetTransitSnapshotColumns][]float64 + columns [planetTransitSnapshotColumns][]float64 + visibilityOnce sync.Once + visibility []PlanetTransitVisibilityColumn } const ( @@ -137,9 +138,9 @@ const ( const ( // planetTransitSnapshotStepDays 是缓存节点间距；四点拉格朗日插值在十分钟间距下误差远低于 1e-6″。 planetTransitSnapshotStepDays = 10.0 / 1440 - // planetTransitSnapshotMarginDays 是缓存相对地心一触、四触向外扩展的余量。 + // planetTransitSnapshotMarginDays 是缓存相对地心一触、四触或地心无接触时求解窗口向外扩展的余量。 planetTransitSnapshotMarginDays = 2.0 / 24 - // planetTransitMaxSpanDays 是地心一触到四触的跨度上限。实际凌日最长约 8 小时（水星）、 + // planetTransitMaxSpanDays 是地心一触到四触或求解窗口的跨度上限。实际凌日最长约 8 小时（水星）、 // 金星不到 9 小时；超出视为异常输入，避免按跨度分配巨大缓存。 planetTransitMaxSpanDays = 1.0 ) @@ -171,10 +172,16 @@ func newPlanetTransitSnapshot(transit PlanetTransitResult, inst *Instance) (*pla } func buildPlanetTransitSnapshot(transit PlanetTransitResult, inst *Instance) (*planetTransitSnapshot, bool) { - if !transit.Valid || !transit.HasExternal || !isFiniteFloat(transit.ExternalIngress) || !isFiniteFloat(transit.ExternalEgress) || - !isFiniteFloat(transit.Greatest) || transit.ExternalEgress <= transit.ExternalIngress || - transit.ExternalEgress-transit.ExternalIngress > planetTransitMaxSpanDays || - transit.Greatest < transit.ExternalIngress || transit.Greatest > transit.ExternalEgress { + start, end := transit.WindowStart, transit.WindowEnd + if transit.HasExternal { + start, end = transit.ExternalIngress, transit.ExternalEgress + } + if !transit.HasExternal && (transit.ExternalIngress != 0 || transit.ExternalEgress != 0 || transit.HasInternal) { + return nil, false + } + if !transit.Valid || !isFiniteFloat(start) || !isFiniteFloat(end) || start == 0 || end == 0 || + !isFiniteFloat(transit.Greatest) || end <= start || end-start > planetTransitMaxSpanDays || + transit.Greatest < start || transit.Greatest > end { return nil, false } var cfg planetTransitConfig @@ -186,8 +193,8 @@ func buildPlanetTransitSnapshot(transit PlanetTransitResult, inst *Instance) (*p default: return nil, false } - startTT := UTC2TT(transit.ExternalIngress) - planetTransitSnapshotMarginDays - endTT := UTC2TT(transit.ExternalEgress) + planetTransitSnapshotMarginDays + startTT := UTC2TT(start) - planetTransitSnapshotMarginDays + endTT := UTC2TT(end) + planetTransitSnapshotMarginDays count := int(math.Ceil((endTT-startTT)/planetTransitSnapshotStepDays)) + 1 snapshot := &planetTransitSnapshot{planetIndex: transit.PlanetIndex, startTT: startTT, step: planetTransitSnapshotStepDays} for column := range snapshot.columns { @@ -445,7 +452,11 @@ func (snapshot *planetTransitSnapshot) localContactTT( if seedValue > 0 { direction = -outward } - limit := int(math.Ceil(planetTransitLocalWindowDays / planetTransitLocalScanStepDays)) + span := planetTransitLocalWindowDays + if math.Abs(seedTT-greatestTT) < planetTransitLocalToleranceDays { + span = math.Max(greatestTT-snapshot.startTT, snapshot.startTT+float64(len(snapshot.columns[0])-1)*snapshot.step-greatestTT) + } + limit := int(math.Ceil(span / planetTransitLocalScanStepDays)) near, nearValue := seedTT, seedValue for step := 1; step <= limit; step++ { far := seedTT + direction*float64(step)*planetTransitLocalScanStepDays @@ -544,13 +555,13 @@ type PlanetTransitIsochrone struct { // PlanetTransitContactIsochrones 返回接触 kind 的站心等时线，时刻按 step 对齐到整刻度（UTC）； // 条数超过 64 时步长按倍数放大，保证覆盖整个接触时段； -// step 非正、kind 不是 C1–C4 或地心凌日没有该接触时返回 nil。 +// step 非正、kind 不是 C1–C4 或地表没有该接触时返回 nil。 // 时刻是站心接触时刻，受视差影响可落在地心 C1–C4 之外（金星约数分钟）； // 一触、四触与远离掠凌边界的二触、三触，线上地点的站心接触时刻与 PlanetTransitLocal 相差不到 1 秒； // 掠凌边界附近（站心内切余量 <0.5″）内切时刻梯度发散，偏差可达数十秒。 // PlanetTransitContactIsochrones returns topocentric isochrones of contact kind at UTC instants aligned to step. // When more than 64 lines would be needed, the step is doubled until the whole contact span fits. -// It returns nil for a non-positive step, a kind other than C1-C4, or a contact the geocentric transit lacks. +// It returns nil for a non-positive step, a kind other than C1-C4, or a contact absent everywhere on Earth. // Points agree with PlanetTransitLocal to within a second, except C2/C3 near the grazing limit (internal margin // below 0.5″), where the contact-time gradient diverges and errors reach tens of seconds. func PlanetTransitContactIsochrones(transit PlanetTransitResult, kind PlanetTransitContactKind, step time.Duration) []PlanetTransitIsochrone { @@ -567,9 +578,6 @@ func PlanetTransitContactIsochronesWithInstance( if step <= 0 || kind < PlanetTransitC1 || kind > PlanetTransitC4 { return nil } - if (kind == PlanetTransitC2 || kind == PlanetTransitC3) && !transit.HasInternal { - return nil - } snapshot, ok := newPlanetTransitSnapshot(transit, inst) if !ok { return nil @@ -971,9 +979,10 @@ func PlanetTransitLocalInfoFromResult(result PlanetTransitLocalResult, loc *time } } -// PlanetTransitResultFromTimes 由各接触时刻重建地心凌日结果，零值时刻表示该接触不存在。 -// 行星包的 TransitInfo 只保存 time.Time，站心计算用它还原 basic 结果。 -// PlanetTransitResultFromTimes rebuilds a geocentric transit result from contact times; a zero time means the contact is absent. +// PlanetTransitResultFromTimes 由各接触时刻重建地心凌日结果，一触、凌甚、四触必填，内切时刻为零表示不存在； +// 地心无接触事件改用 PlanetTransitResultFromWindow。 +// PlanetTransitResultFromTimes rebuilds a geocentric transit result from contact times; C1, greatest and C4 are required, +// and a zero internal contact means it is absent. Use PlanetTransitResultFromWindow for events without geocentric contacts. func PlanetTransitResultFromTimes(planetIndex int, start, internalStart, greatest, internalEnd, end time.Time) PlanetTransitResult { jd := func(value time.Time) float64 { if value.IsZero() { @@ -992,9 +1001,27 @@ func PlanetTransitResultFromTimes(planetIndex int, start, internalStart, greates HasExternal: true, } result.HasInternal = result.InternalIngress != 0 && result.InternalEgress != 0 + if result.Valid { + result.WindowStart = result.ExternalIngress - planetTransitSnapshotMarginDays + result.WindowEnd = result.ExternalEgress + planetTransitSnapshotMarginDays + } return result } +// PlanetTransitResultFromWindow 由凌甚与站心求解窗口重建地心无接触的凌日结果 / rebuild a transit without geocentric contacts from greatest and its local search window. +func PlanetTransitResultFromWindow(planetIndex int, windowStart, greatest, windowEnd time.Time) PlanetTransitResult { + if windowStart.IsZero() || greatest.IsZero() || windowEnd.IsZero() || !windowStart.Before(greatest) || !greatest.Before(windowEnd) { + return PlanetTransitResult{PlanetIndex: planetIndex} + } + return PlanetTransitResult{ + Valid: true, + PlanetIndex: planetIndex, + Greatest: Date2JD(greatest.UTC()), + WindowStart: Date2JD(windowStart.UTC()), + WindowEnd: Date2JD(windowEnd.UTC()), + } +} + // PlanetTransitLocalFrame 是站心日面上一个时刻的行星位置。 // PlanetTransitLocalFrame is the topocentric planet position on the solar disk at one instant. type PlanetTransitLocalFrame struct { diff --git a/basic/planet_transit_prediction.go b/basic/planet_transit_prediction.go new file mode 100644 index 0000000..c6cb89f --- /dev/null +++ b/basic/planet_transit_prediction.go @@ -0,0 +1,69 @@ +package basic + +import ( + "math" + + "b612.me/astro/internal/ephem" +) + +func planetTransitPredictedSeason(seedTT float64, cfg planetTransitConfig, inst *Instance) bool { + centuries := (seedTT - 2451545) / 36525 + if inst != nil && inst.Model() == ephem.Custom || centuries < -40 || centuries > 20 { + return true + } + node := 76.679920 + centuries*(0.9011206+centuries*(0.00040618-0.000000093*centuries)) + if cfg.planetIndex == 1 { + node = 48.330893 + centuries*(1.1861883+centuries*(0.00017542+0.000000215*centuries)) + } + solar := 280.46646 + centuries*(36000.76983+0.0003032*centuries) + phase := math.Abs(math.Remainder(solar-node, 180)) + return phase <= 24 +} + +func planetTransitCoarseConjunction(seedTT float64, cfg planetTransitConfig, inst *Instance) float64 { + current := seedTT + for iteration := 0; iteration < 5; iteration++ { + value := planetTransitLongitudeDeltaN(current, cfg, planetTransitCoarseN, inst) + slope := (planetTransitLongitudeDeltaN(current+0.02, cfg, planetTransitCoarseN, inst) - + planetTransitLongitudeDeltaN(current-0.02, cfg, planetTransitCoarseN, inst)) / 0.04 + if !isFiniteFloat(slope) || slope == 0 { + return math.NaN() + } + step := value / slope + current -= step + if math.Abs(step) < 0.0001 { + break + } + } + return current +} + +// PlanetTransitSeries 凌日周期族，Number 为零表示超出编号适用范围 / transit family; a zero Number means unassigned. +type PlanetTransitSeries struct { + Number int + PeriodYears int +} + +func planetTransitSeriesAt(greatestTT float64, cfg planetTransitConfig) PlanetTransitSeries { + period, synodicCycles := 243, 152 + firstYear, lastYear := -2000, 4000 + anchors := [][4]int{{-1892, 5, 21, 1}, {-1998, 11, 18, 2}, {-1884, 5, 19, 3}, {-1763, 11, 20, 4}, {554, 5, 22, 5}, {1631, 12, 7, 6}} + if cfg.planetIndex == 1 { + period, synodicCycles = 46, 145 + firstYear, lastYear = 1601, 2300 + anchors = [][4]int{{1644, 11, 9, 1}, {1631, 11, 7, 2}, {1628, 5, 5, 3}, {1618, 11, 4, 4}, {1615, 5, 3, 5}, {1605, 11, 1, 6}, {1740, 5, 2, 7}, {1776, 11, 2, 8}, {1957, 5, 6, 9}, {1993, 11, 6, 10}, {2174, 5, 8, 11}, {2210, 11, 9, 12}} + } + result := PlanetTransitSeries{PeriodYears: period} + if greatestTT < JDCalc(firstYear, 1, 1) || greatestTT >= JDCalc(lastYear+1, 1, 1) { + return result + } + cycle := int(math.Round((greatestTT - cfg.anchorInferiorTT) / cfg.synodicPeriodDays)) + for _, anchor := range anchors { + anchorCycle := int(math.Round((JDCalc(anchor[0], anchor[1], float64(anchor[2])) - cfg.anchorInferiorTT) / cfg.synodicPeriodDays)) + if (cycle-anchorCycle)%synodicCycles == 0 { + result.Number = anchor[3] + return result + } + } + return result +} diff --git a/basic/planet_transit_visibility.go b/basic/planet_transit_visibility.go new file mode 100644 index 0000000..2347618 --- /dev/null +++ b/basic/planet_transit_visibility.go @@ -0,0 +1,79 @@ +package basic + +import "math" + +// PlanetTransitVisibilityColumn 子午线上的可见纬度区间，端点单位为度 / visible latitude intervals on a meridian, in degrees. +type PlanetTransitVisibilityColumn struct { + Longitude float64 + Intervals [][2]float64 +} + +// PlanetTransitVisibilityColumnsWithInstance 按 2° 经度、1° 纬度探测站心可见区并细化边界 / sample topocentric visibility at 2-degree longitude and 1-degree latitude spacing and refine boundaries. +func PlanetTransitVisibilityColumnsWithInstance(transit PlanetTransitResult, inst *Instance) []PlanetTransitVisibilityColumn { + snapshot, ok := newPlanetTransitSnapshot(transit, inst) + if !ok { + return nil + } + snapshot.visibilityOnce.Do(func() { + for column := 0; column <= 180; column++ { + longitude := -180 + float64(column)*2 + visible := func(latitude float64) bool { + local := snapshot.local(transit, Observer{Longitude: longitude, Latitude: latitude}) + return local.Valid && local.Visible + } + value := PlanetTransitVisibilityColumn{Longitude: longitude} + previous, start := visible(-90), -90.0 + for latitude := -89.0; latitude <= 90; latitude++ { + current := visible(latitude) + if current != previous { + left, right := latitude-1, latitude + for right-left > 0.0001 { + middle := (left + right) / 2 + if visible(middle) == previous { + left = middle + } else { + right = middle + } + } + boundary := (left + right) / 2 + if current { + start = boundary + } else { + value.Intervals = append(value.Intervals, [2]float64{start, boundary}) + } + } + previous = current + } + if previous { + value.Intervals = append(value.Intervals, [2]float64{start, 90}) + } + snapshot.visibility = append(snapshot.visibility, value) + } + }) + result := make([]PlanetTransitVisibilityColumn, len(snapshot.visibility)) + for index, column := range snapshot.visibility { + result[index] = PlanetTransitVisibilityColumn{Longitude: column.Longitude, Intervals: append([][2]float64(nil), column.Intervals...)} + } + return result +} + +// PlanetTransitVisibilityPolygonsWithInstance 可见区的经度条带多边形 / longitude-strip polygons of the visible region. +func PlanetTransitVisibilityPolygonsWithInstance(transit PlanetTransitResult, inst *Instance) [][]PlanetTransitIsochronePoint { + columns := PlanetTransitVisibilityColumnsWithInstance(transit, inst) + var result [][]PlanetTransitIsochronePoint + for index := 1; index < len(columns); index++ { + left, right := columns[index-1], columns[index] + for _, interval := range left.Intervals { + for _, other := range right.Intervals { + if math.Min(interval[1], other[1]) <= math.Max(interval[0], other[0]) { + continue + } + result = append(result, []PlanetTransitIsochronePoint{ + {Longitude: left.Longitude, Latitude: interval[0]}, {Longitude: right.Longitude, Latitude: other[0]}, + {Longitude: right.Longitude, Latitude: other[1]}, {Longitude: left.Longitude, Latitude: interval[1]}, + }) + } + } + } + return result +} diff --git a/doc/img/transit-mercury-1937-partial-en.svg b/doc/img/transit-mercury-1937-partial-en.svg new file mode 100644 index 0000000..0418278 --- /dev/null +++ b/doc/img/transit-mercury-1937-partial-en.svg @@ -0,0 +1 @@ +1937-05-11 Transit of Mercury Global VisibilityGreatest 16:59:16 | Partial transit; no geocentric contacts (CST)γ -134.2117 | Series 5 (46 yr) | least separation 955.5″ | solar semidiameter 949.5″ | planet semidiameter 6.05″ (geocentric)Partial transit visibleOrthographic globe projection, centred on the topocentric visible footprint; one hemisphere only; shaded footprint uses topocentric diskoverlap and apparent-horizon visibility (standard refraction).All times are UTC (shown in CST, UTC+08:00) \ No newline at end of file diff --git a/doc/img/transit-mercury-1937-partial.svg b/doc/img/transit-mercury-1937-partial.svg new file mode 100644 index 0000000..e003f9f --- /dev/null +++ b/doc/img/transit-mercury-1937-partial.svg @@ -0,0 +1 @@ +1937-05-11 水星凌日全球可见图凌甚 16:59:16 | 仅部分凌日，地心无接触 (CST)γ -134.2117 | 周期族 5（46 年） | 最小角距 955.5″ | 太阳视半径 949.5″ | 行星视半径 6.05″（地心）部分凌日可见区正射球面投影，视点取站心可见区中心，只画朝向视点的半个地球；着色区域按站心圆面相交与视地平可见性计算（计入标准大气折射）。图中时刻为 UTC（显示时区 CST，UTC+08:00） \ No newline at end of file diff --git a/doc/img/transit-mercury-beijing-2032-en.svg b/doc/img/transit-mercury-beijing-2032-en.svg index 0df19f4..18c5230 100644 --- a/doc/img/transit-mercury-beijing-2032-en.svg +++ b/doc/img/transit-mercury-beijing-2032-en.svg @@ -1 +1 @@ -2032-11-13 Transit of Mercury, Local CircumstancesSite 116.4074°E 39.9042°N, height 0 mTopocentric full transit, visible | minimum separation 568.3″ | Sun 969.4″ | planet 4.97″NSEWC1 14:40:40C2 14:42:44Greatest 16:53:54C3 19:05:16C4 19:07:20ContactsTimeSun alt.Az.C114:40:40+20.8°221.3°C214:42:44+20.5°221.8°Greatest16:53:54+0.2°245.9°C319:05:16-24.1°265.9°C419:07:20-24.5°266.3°Duration C1–C4 4h26m41sSunset during transit 16:58:19Sky view, north up and east left. Topocentric contacts; Sun altitudes are geometric centre altitudes, grey parts havethe Sun's centre below the apparent horizon (standard refraction).All times are UTC (shown in CST, UTC+08:00) \ No newline at end of file +2032-11-13 Transit of Mercury, Local CircumstancesSite 116.4074°E 39.9042°N, height 0 mTopocentric full transit, visible | minimum separation 568.3″ | Sun 969.4″ | planet 4.97″NSEWC1 14:40:40C2 14:42:44Greatest 16:53:54C3 19:05:16C4 19:07:20ContactsTimeSun alt.Az.C114:40:40+20.8°221.3°C214:42:44+20.5°221.8°Greatest16:53:54+0.2°245.9°C319:05:16-24.1°265.9°C419:07:20-24.5°266.3°Duration C1–C4 4h26m41sSunset during transit 16:58:19Sky view, north up and east left. Topocentric contacts; Sun altitudes are geometric centre altitudes, grey parts havethe Sun's centre below the apparent horizon (standard refraction).All times are UTC (shown in CST, UTC+08:00) \ No newline at end of file diff --git a/doc/img/transit-mercury-beijing-2032-global-en.svg b/doc/img/transit-mercury-beijing-2032-global-en.svg index 8bbe5e9..58e9a8e 100644 --- a/doc/img/transit-mercury-beijing-2032-global-en.svg +++ b/doc/img/transit-mercury-beijing-2032-global-en.svg @@ -1 +1 @@ -2032-11-13 Transit of Mercury Global VisibilityC1 14:41:13 | C2 14:43:18 | Greatest 16:54:13 | C3 19:05:11 | C4 19:07:16 (CST)least separation 572.1″ | solar semidiameter 969.4″ | planet semidiameter 4.97″ (geocentric)Entire transitTransit at sunriseTransit at sunsetNot visibleEquirectangular projection; regions follow whether the Sun's centre is above the apparent horizon (standard refraction) at C1, greatest and C4.All times are UTC (shown in CST, UTC+08:00) \ No newline at end of file +2032-11-13 Transit of Mercury Global VisibilityC1 14:41:13 | C2 14:43:18 | Greatest 16:54:13 | C3 19:05:11 | C4 19:07:16 (CST)γ 139.0760 | Series 4 (46 yr) | least separation 572.1″ | solar semidiameter 969.4″ | planet semidiameter 4.97″ (geocentric)Entire transitTransit at sunriseTransit at sunsetNot visibleEquirectangular projection; regions follow whether the Sun's centre is above the apparent horizon (standard refraction) at C1, greatest and C4.All times are UTC (shown in CST, UTC+08:00) \ No newline at end of file diff --git a/doc/img/transit-mercury-beijing-2032-global.svg b/doc/img/transit-mercury-beijing-2032-global.svg index 28af4dd..564ad62 100644 --- a/doc/img/transit-mercury-beijing-2032-global.svg +++ b/doc/img/transit-mercury-beijing-2032-global.svg @@ -1 +1 @@ -2032-11-13 水星凌日全球可见图C1 14:41:13 | C2 14:43:18 | 凌甚 16:54:13 | C3 19:05:11 | C4 19:07:16 (CST)最小角距 572.1″ | 太阳视半径 969.4″ | 行星视半径 4.97″（地心）全程可见日出时已在凌日日落时仍在凌日不可见等经纬投影；按一触、凌甚、四触太阳中心是否在视地平上分区（计入标准大气折射）。图中时刻为 UTC（显示时区 CST，UTC+08:00） \ No newline at end of file +2032-11-13 水星凌日全球可见图C1 14:41:13 | C2 14:43:18 | 凌甚 16:54:13 | C3 19:05:11 | C4 19:07:16 (CST)γ 139.0760 | 周期族 4（46 年） | 最小角距 572.1″ | 太阳视半径 969.4″ | 行星视半径 4.97″（地心）全程可见日出时已在凌日日落时仍在凌日不可见等经纬投影；按一触、凌甚、四触太阳中心是否在视地平上分区（计入标准大气折射）。图中时刻为 UTC（显示时区 CST，UTC+08:00） \ No newline at end of file diff --git a/doc/img/transit-mercury-beijing-2032.svg b/doc/img/transit-mercury-beijing-2032.svg index 3c762b9..facb623 100644 --- a/doc/img/transit-mercury-beijing-2032.svg +++ b/doc/img/transit-mercury-beijing-2032.svg @@ -1 +1 @@ -2032-11-13 水星凌日站心图观测地 东经 116.4074° 北纬 39.9042°，海拔 0 米站心完整凌日，可见 | 最小角距 568.3″ | 太阳视半径 969.4″ | 行星视半径 4.97″北南东西一触 14:40:40二触 14:42:44凌甚 16:53:54三触 19:05:16四触 19:07:20接触时刻时刻太阳高度方位一触14:40:40+20.8°221.3°二触14:42:44+20.5°221.8°凌甚16:53:54+0.2°245.9°三触19:05:16-24.1°265.9°四触19:07:20-24.5°266.3°一触至四触历时 4h26m41s凌日中日落 16:58:19天空视向，上北左东；站心接触时刻，太阳高度为日心几何高度，灰色部分太阳中心在视地平下（计入标准大气折射）。图中时刻为 UTC（显示时区 CST，UTC+08:00） \ No newline at end of file +2032-11-13 水星凌日站心图观测地 东经 116.4074° 北纬 39.9042°，海拔 0 米站心完整凌日，可见 | 最小角距 568.3″ | 太阳视半径 969.4″ | 行星视半径 4.97″北南东西一触 14:40:40二触 14:42:44凌甚 16:53:54三触 19:05:16四触 19:07:20接触时刻时刻太阳高度方位一触14:40:40+20.8°221.3°二触14:42:44+20.5°221.8°凌甚16:53:54+0.2°245.9°三触19:05:16-24.1°265.9°四触19:07:20-24.5°266.3°一触至四触历时 4h26m41s凌日中日落 16:58:19天空视向，上北左东；站心接触时刻，太阳高度为太阳中心几何高度，灰色部分太阳中心在视地平下（计入标准大气折射）。图中时刻为 UTC（显示时区 CST，UTC+08:00） \ No newline at end of file diff --git a/doc/img/transit-venus-beijing-2012-en.svg b/doc/img/transit-venus-beijing-2012-en.svg index 532f534..739957e 100644 --- a/doc/img/transit-venus-beijing-2012-en.svg +++ b/doc/img/transit-venus-beijing-2012-en.svg @@ -1 +1 @@ -2012-06-06 Transit of Venus, Local CircumstancesSite 116.4074°E 39.9042°N, height 0 mTopocentric full transit, visible | minimum separation 544.0″ | Sun 945.3″ | planet 28.90″NSEWC1 06:10:13C2 06:27:52Greatest 09:30:40C3 12:31:56C4 12:49:17ContactsTimeSun alt.Az.C106:10:13+13.8°71.8°C206:27:52+17.0°74.4°Greatest09:30:40+51.7°104.4°C312:31:56+72.3°194.5°C412:49:17+71.2°206.7°Duration C1–C4 6h39m4sSky view, north up and east left. Topocentric contacts; Sun altitudes are geometric centre altitudes, grey parts havethe Sun's centre below the apparent horizon (standard refraction).All times are UTC (shown in CST, UTC+08:00) \ No newline at end of file +2012-06-06 Transit of Venus, Local CircumstancesSite 116.4074°E 39.9042°N, height 0 mTopocentric full transit, visible | minimum separation 544.0″ | Sun 945.3″ | planet 28.90″NSEWC1 06:10:13C2 06:27:52Greatest 09:30:40C3 12:31:56C4 12:49:17ContactsTimeSun alt.Az.C106:10:13+13.8°71.8°C206:27:52+17.0°74.4°Greatest09:30:40+51.7°104.4°C312:31:56+72.3°194.5°C412:49:17+71.2°206.7°Duration C1–C4 6h39m4sSky view, north up and east left. Topocentric contacts; Sun altitudes are geometric centre altitudes, grey parts havethe Sun's centre below the apparent horizon (standard refraction).All times are UTC (shown in CST, UTC+08:00) \ No newline at end of file diff --git a/doc/img/transit-venus-beijing-2012-global-en.svg b/doc/img/transit-venus-beijing-2012-global-en.svg index 1dd0ed6..7eec6f0 100644 --- a/doc/img/transit-venus-beijing-2012-global-en.svg +++ b/doc/img/transit-venus-beijing-2012-global-en.svg @@ -1 +1 @@ -2012-06-06 Transit of Venus Global VisibilityC1 06:09:48 | C2 06:27:37 | Greatest 09:29:36 | C3 12:31:35 | C4 12:49:24 (CST)least separation 554.4″ | solar semidiameter 945.3″ | planet semidiameter 28.90″ (geocentric)Entire transitTransit at sunriseTransit at sunsetNot visibleEquirectangular projection; regions follow whether the Sun's centre is above the apparent horizon (standard refraction) at C1, greatest and C4.All times are UTC (shown in CST, UTC+08:00) \ No newline at end of file +2012-06-06 Transit of Venus Global VisibilityC1 06:09:48 | C2 06:27:37 | Greatest 09:29:36 | C3 12:31:35 | C4 12:49:24 (CST)γ 25.4361 | Series 5 (243 yr) | least separation 554.4″ | solar semidiameter 945.3″ | planet semidiameter 28.90″ (geocentric)Entire transitTransit at sunriseTransit at sunsetNot visibleEquirectangular projection; regions follow whether the Sun's centre is above the apparent horizon (standard refraction) at C1, greatest and C4.All times are UTC (shown in CST, UTC+08:00) \ No newline at end of file diff --git a/doc/img/transit-venus-beijing-2012-global.svg b/doc/img/transit-venus-beijing-2012-global.svg index 7750fa1..72bc097 100644 --- a/doc/img/transit-venus-beijing-2012-global.svg +++ b/doc/img/transit-venus-beijing-2012-global.svg @@ -1 +1 @@ -2012-06-06 金星凌日全球可见图C1 06:09:48 | C2 06:27:37 | 凌甚 09:29:36 | C3 12:31:35 | C4 12:49:24 (CST)最小角距 554.4″ | 太阳视半径 945.3″ | 行星视半径 28.90″（地心）全程可见日出时已在凌日日落时仍在凌日不可见等经纬投影；按一触、凌甚、四触太阳中心是否在视地平上分区（计入标准大气折射）。图中时刻为 UTC（显示时区 CST，UTC+08:00） \ No newline at end of file +2012-06-06 金星凌日全球可见图C1 06:09:48 | C2 06:27:37 | 凌甚 09:29:36 | C3 12:31:35 | C4 12:49:24 (CST)γ 25.4361 | 周期族 5（243 年） | 最小角距 554.4″ | 太阳视半径 945.3″ | 行星视半径 28.90″（地心）全程可见日出时已在凌日日落时仍在凌日不可见等经纬投影；按一触、凌甚、四触太阳中心是否在视地平上分区（计入标准大气折射）。图中时刻为 UTC（显示时区 CST，UTC+08:00） \ No newline at end of file diff --git a/doc/img/transit-venus-beijing-2012.svg b/doc/img/transit-venus-beijing-2012.svg index 6f7c2ab..0035c01 100644 --- a/doc/img/transit-venus-beijing-2012.svg +++ b/doc/img/transit-venus-beijing-2012.svg @@ -1 +1 @@ -2012-06-06 金星凌日站心图观测地 东经 116.4074° 北纬 39.9042°，海拔 0 米站心完整凌日，可见 | 最小角距 544.0″ | 太阳视半径 945.3″ | 行星视半径 28.90″北南东西一触 06:10:13二触 06:27:52凌甚 09:30:40三触 12:31:56四触 12:49:17接触时刻时刻太阳高度方位一触06:10:13+13.8°71.8°二触06:27:52+17.0°74.4°凌甚09:30:40+51.7°104.4°三触12:31:56+72.3°194.5°四触12:49:17+71.2°206.7°一触至四触历时 6h39m4s天空视向，上北左东；站心接触时刻，太阳高度为日心几何高度，灰色部分太阳中心在视地平下（计入标准大气折射）。图中时刻为 UTC（显示时区 CST，UTC+08:00） \ No newline at end of file +2012-06-06 金星凌日站心图观测地 东经 116.4074° 北纬 39.9042°，海拔 0 米站心完整凌日，可见 | 最小角距 544.0″ | 太阳视半径 945.3″ | 行星视半径 28.90″北南东西一触 06:10:13二触 06:27:52凌甚 09:30:40三触 12:31:56四触 12:49:17接触时刻时刻太阳高度方位一触06:10:13+13.8°71.8°二触06:27:52+17.0°74.4°凌甚09:30:40+51.7°104.4°三触12:31:56+72.3°194.5°四触12:49:17+71.2°206.7°一触至四触历时 6h39m4s天空视向，上北左东；站心接触时刻，太阳高度为太阳中心几何高度，灰色部分太阳中心在视地平下（计入标准大气折射）。图中时刻为 UTC（显示时区 CST，UTC+08:00） \ No newline at end of file diff --git a/doc/manual/eclipse.md b/doc/manual/eclipse.md index 9437645..1deaaa5 100644 --- a/doc/manual/eclipse.md +++ b/doc/manual/eclipse.md @@ -33,6 +33,7 @@ - [凌日出图](#凌日出图) - [水星：2032-11-13（北京）](#水星2032-11-13北京) - [金星：2012-06-06（北京）](#金星2012-06-06北京) + - [特例：地心无接触的部分凌日](#特例地心无接触的部分凌日) ## 简单示例：2009长江大日食，上海某地见食情况 @@ -990,10 +991,14 @@ ut1, ok := eclipsesvg.SolarEclipseMapSVG( | 入口 | 画面 | 画布建议 | | --- | --- | --- | -| `TransitMapSVG` | 地心几何的全球可见图：全程可见、日出时已在凌日、日落时仍在凌日、不可见四区 | 960×640 起 | +| `TransitMapSVG` | 全球可见图；地心无接触时绘制站心部分凌日可见区 | 960×640 起 | | `TransitLocalSVG` | 站心日面图：上北左东，行星中心轨迹、各接触与凌甚的行星圆面、接触时刻表 | 960×680 起 | -可见区按一触、凌甚、四触三刻太阳中心是否在视地平上划分，计入标准大气折射（地平圈半径 90°34′）。等时线默认不画，用 `Isochrones` 按位开启：`TransitIsochroneC1` … `TransitIsochroneC4`，或组合 `TransitIsochroneExternal`（一触、四触）、`TransitIsochroneAll`；掠凌没有内切接触时 C2/C3 自动略去。等时线是站心接触时刻，只画太阳在视地平上的部分；`IsochroneStep` 非正时水星取 20 秒、金星取 1 分钟。凌日可见范围接近半个地球，地图始终画全球，`Projection` 可选四种投影。 +普通凌日可见区按一触、凌甚、四触三刻太阳中心是否在视地平上划分，计入标准大气折射（地平圈半径 90°34′）；地心无接触事件改用实际站心可见区。等时线默认不画，用 `Isochrones` 按位开启：`TransitIsochroneC1` … `TransitIsochroneC4`，或组合 `TransitIsochroneExternal`（一触、四触）、`TransitIsochroneAll`；地表没有站心内切接触时 C2/C3 不画，图例也不列出；地心只有外切的掠凌仍可能在部分地区出现站心二、三触，此时照常绘制。等时线是站心接触时刻，只画太阳在视地平上的部分；`IsochroneStep` 非正时水星取 20 秒、金星取 1 分钟。 + +地心无接触事件的部分凌日可见区按 2° 经度、1° 纬度探测并细化已检出的边界，这只是制图采样，不是站心接触求解精度；比采样间距更窄的可见区可能未被检出，需要精确的地点判断时直接用 `LocalTransit`。 + +`Projection` 可选四种投影，默认等经纬图展示全球。正射图只展示朝向视点的半个地球：水星、金星在地心无接触时以站心可见区的球面面积加权中心取景，例如 `1937-05-11` 水星凌日和 `2854-12-14` 金星凌日；普通凌日仍以凌甚太阳直射点取景。极区事件也可显式选择南、北极投影。 ```go venus := time.Date(2012, 6, 6, 0, 0, 0, 0, time.UTC) @@ -1022,3 +1027,9 @@ local, ok := eclipsesvg.TransitLocalSVG(eclipsesvg.TransitVenus, venus, 116.4, 3 ![2012 金星凌日北京站心图](../img/transit-venus-beijing-2012.svg) ![2012 金星凌日全球可见图](../img/transit-venus-beijing-2012-global.svg) + +### 特例：地心无接触的部分凌日 + +1937-05-11 水星凌日即为一例，可见区横跨南大西洋、非洲南端与印度洋，东至澳大利亚东侧。`2854-12-14` 金星凌日同属此类，可见区跨越南极圈，集中在南纬 50° 以南。 + +![1937 水星凌日（地心无接触）全球可见图](../img/transit-mercury-1937-partial.svg) diff --git a/doc/manual/en/eclipse.md b/doc/manual/en/eclipse.md index e017d4b..3b231b9 100644 --- a/doc/manual/en/eclipse.md +++ b/doc/manual/en/eclipse.md @@ -33,6 +33,7 @@ This manual covers solar- and lunar-eclipse queries and charts, plus the two Mer - [Transit charts](#transit-charts) - [Mercury: 2032-11-13 (Beijing)](#mercury-2032-11-13-beijing) - [Venus: 2012-06-06 (Beijing)](#venus-2012-06-06-beijing) + - [Special case: partial transits without geocentric contacts](#special-case-partial-transits-without-geocentric-contacts) ## Simple example: the 2009 Great Yangtze Eclipse at a site in Shanghai @@ -1064,10 +1065,14 @@ Transit contact times, topocentric visibility and `LocalTransit` live in [Planet | Entry | Chart | Suggested canvas | | --- | --- | --- | -| `TransitMapSVG` | Geocentric global visibility: entire transit, transit at sunrise, transit at sunset, not visible | 960×640 and up | +| `TransitMapSVG` | Global visibility, including topocentric partial-only events without geocentric contacts | 960×640 and up | | `TransitLocalSVG` | Topocentric solar disk: north up, east left, planet-center track, planet disk at each contact and at greatest, contact table | 960×680 and up | -Regions follow whether the Sun's center is above the apparent horizon at C1, greatest and C4, with standard refraction (a horizon circle of 90°34′). Isochrones are off by default; enable them bitwise with `Isochrones`: `TransitIsochroneC1` … `TransitIsochroneC4`, or the combinations `TransitIsochroneExternal` (C1 and C4) and `TransitIsochroneAll`. C2/C3 are dropped for a grazing transit without internal contacts. Isochrones are topocentric contact times drawn only where the Sun is above the apparent horizon; a non-positive `IsochroneStep` uses 20 s for Mercury and 1 min for Venus. A transit is visible from roughly half the globe, so the map always shows the whole world; `Projection` selects one of the four projections. +For ordinary transits, regions follow whether the Sun's center is above the apparent horizon at C1, greatest and C4, with standard refraction (a horizon circle of 90°34′); events without geocentric contacts use the actual topocentric visible footprint. Isochrones are off by default; enable them bitwise with `Isochrones`: `TransitIsochroneC1` … `TransitIsochroneC4`, or the combinations `TransitIsochroneExternal` (C1 and C4) and `TransitIsochroneAll`. C2/C3 are omitted, including from the legend, when no site has an internal contact; a graze with only geocentric external contacts can still have topocentric C2/C3 in some regions, and those are drawn. Isochrones are topocentric contact times drawn only where the Sun is above the apparent horizon; a non-positive `IsochroneStep` uses 20 s for Mercury and 1 min for Venus. + +For events without geocentric contacts the partial footprint is sampled on a 2° longitude by 1° latitude grid and detected boundaries are refined. That is chart sampling, not contact accuracy: a footprint narrower than the grid may go undetected, so use `LocalTransit` for exact site decisions. + +`Projection` selects one of four projections; the default equirectangular map shows the whole world. Orthographic maps show only the facing hemisphere. For both Mercury and Venus, events without geocentric contacts centre on the spherical area-weighted visible footprint, for example Mercury on `1937-05-11` and Venus on `2854-12-14`. Ordinary transits retain the subsolar centre at greatest. Polar events may also use an explicit north- or south-polar projection. ```go venus := time.Date(2012, 6, 6, 0, 0, 0, 0, time.UTC) @@ -1097,3 +1102,9 @@ The Venus transit of 2012-06-06 was visible in full from Beijing, with `6h39m4s` ![Venus transit of 2012-06-06 from Beijing](../../img/transit-venus-beijing-2012-en.svg) ![Venus transit of 2012-06-06 global visibility](../../img/transit-venus-beijing-2012-global-en.svg) + +### Special case: partial transits without geocentric contacts + +The Mercury transit of 1937-05-11 is one such event: with no geocentric contacts, its footprint spans the South Atlantic, southern Africa and the Indian Ocean east to eastern Australia. The Venus transit of `2854-12-14` belongs to the same class, with a footprint that crosses the Antarctic Circle, all south of 50°S. + +![Mercury transit of 1937-05-11 without geocentric contacts, global visibility](../../img/transit-mercury-1937-partial-en.svg) diff --git a/doc/manual/en/map-geojson.md b/doc/manual/en/map-geojson.md index 9f3054f..56b6289 100644 --- a/doc/manual/en/map-geojson.md +++ b/doc/manual/en/map-geojson.md @@ -73,7 +73,11 @@ A partial eclipse may have no central path, so `centralPath` may be `nil`. This ## Exporting transit GeoJSON -The transit entry points take a geocentric `TransitInfo` and partition the globe exactly like the transit visibility map: +The transit entry points take a global-event `TransitInfo`; contact fields remain geocentric, and regions match the transit visibility map. + +Events without geocentric contacts, such as the Venus transit of 2854-12-14, export the actual topocentric +`transit-partial` footprint and omit nonexistent geocentric C1/C4 horizons. Properties also carry +`has_geocentric_contacts`, `gamma`, `penumbral_radius`, `series` and `series_period_years`. The footprint is sampled on a 2° longitude by 1° latitude grid with boundary refinement; that is chart sampling rather than contact accuracy, so a footprint narrower than the grid may go undetected — use `LocalTransit` for exact site decisions. ```go m := mercury.NextTransit(time.Date(2026, 10, 12, 0, 0, 0, 0, cst)) @@ -194,7 +198,7 @@ Grouped by figure family, the `data-source` values are listed below; the `eclips - Solar-eclipse instants and magnitude lines: `rise-set-phase-lines`, `magnitude-contours`, `greatest-time-isochrones` - Instantaneous outlines: `penumbral-outlines`, `central-shadow-outlines` - Lunar eclipses: `p1-p4-visibility-regions`, `p1-p4-horizon-boundaries` -- Transits: `c1-c4-visibility-regions`, `c1-c4-horizon-boundaries`, `transit-contact-isochrones` +- Transits: `c1-c4-visibility-regions`, `c1-c4-horizon-boundaries`, `transit-contact-isochrones`; events without geocentric contacts fill the footprint from `topocentric-transit-visibility` ## Map Projections @@ -312,7 +316,7 @@ Each event type has plain and time-marker variants: - `MarshalLunarEclipseWithOptions`: one `LunarEclipseOptions` value carrying the time markers, `SkipRoles` and the time-envelope sampling parameters. `SkipRoles` can omit envelopes and penumbra-only bands; skipping both envelopes also skips their sampling. `EnvelopeSweepSamples` defaults to 48 and is clamped to `[2, 192]`; `EnvelopeLongitudePoints` defaults to `max(360, boundaryPoints)` and is clamped to `[12, 720]`. - `MarshalStarOccultation` / `MarshalStarOccultationWithTimeMarkers` - `MarshalPlanetOccultation` / `MarshalPlanetOccultationWithTimeMarkers` -- `MarshalMercuryTransit` / `MarshalVenusTransit`: take `mercury.TransitInfo` / `venus.TransitInfo` and partition the globe like the transit SVG map (whether the Sun's center is above the apparent horizon, standard refraction, at C1, greatest and C4). `BoundaryPoints` is the horizon sample count and also sets the region longitude columns to `max(360, BoundaryPoints)`, capped at 720. Roles are the three regions `transit-entire`, `transit-at-sunrise` and `transit-at-sunset`, the visible hemispheres and horizons `visible-at-c1` / `visible-at-c4` and `c1-horizon` / `c4-horizon`, `greatest` (the subsolar point) and the time markers `time-marker`. A positive `IsochroneStep` (at least 1 s) adds topocentric `c1-isochrone` through `c4-isochrone` (no C2/C3 for a grazing transit; the times are topocentric contacts and parallax can place them outside the geocentric C1–C4 interval), and a contact whose isochrone does not reach the visible hemisphere omits that role; skipping all three region roles skips the region computation. +- `MarshalMercuryTransit` / `MarshalVenusTransit`: take `mercury.TransitInfo` / `venus.TransitInfo` and partition the globe like the transit SVG map (whether the Sun's center is above the apparent horizon, standard refraction, at C1, greatest and C4). `BoundaryPoints` is the horizon sample count and also sets the region longitude columns to `max(360, BoundaryPoints)`, capped at 720. Roles are the three regions `transit-entire`, `transit-at-sunrise` and `transit-at-sunset`, the visible hemispheres and horizons `visible-at-c1` / `visible-at-c4` and `c1-horizon` / `c4-horizon`, `greatest` (the subsolar point) and the time markers `subsolar-track`. Events without geocentric contacts switch to a `transit-partial` region and omit the three visibility regions, the horizons and the visible hemispheres. A positive `IsochroneStep` (at least 1 s) adds topocentric `c1-isochrone` through `c4-isochrone` (no C2/C3 when no site has an internal contact, though a graze with only geocentric external contacts can still have them in some regions; the times are topocentric contacts and parallax can place them outside the geocentric C1–C4 interval), and a contact whose isochrone does not reach the visible hemisphere omits that role; skipping all three region roles skips the region computation. - `MarshalSolarEclipseWithOptions`: one `SolarEclipseOptions` value carrying both the time markers and `SkipRoles`, equivalent to the solar entry points above plus layer trimming. The usual entry in `SkipRoles` is `partial-footprint` (the instantaneous penumbral outlines); `partial-band` is the edge of the eclipse band (the magnitude-0 limit) and is normally kept. @@ -472,7 +476,7 @@ Features are grouped into Folders by `event` and `role`. Collections containing | Total and partial occultation bands | `total-band`, `partial-band`, `total-band-outline`, `occultation-band` | Yellow with translucent band fill | | Rise/set visibility edges, lunar time envelopes and solar partial-band edge | `visibility-boundary`, `p1-horizon`, `p4-horizon`, `visible-at-p1`, `visible-at-p4`, `visible-during-eclipse`, `visible-throughout-eclipse`, solar `partial-band` | Orange | | Penumbra-only moonrise/moonset bands | `penumbra-moonrise`, `penumbra-moonset` | Blue-violet for moonrise, violet for moonset, with a translucent fill | -| Transit visibility regions and contact lines | `transit-entire`, `transit-at-sunrise`, `transit-at-sunset`, `c1-horizon`, `c4-horizon`, `visible-at-c1`, `visible-at-c4`, `c1-isochrone` … `c4-isochrone` | Each of the three regions has its own fill; C1/C4 horizons and visible hemispheres are orange; the four contact isochrones each get a color | +| Transit visibility regions and contact lines | `transit-entire`, `transit-at-sunrise`, `transit-at-sunset`, `transit-partial`, `c1-horizon`, `c4-horizon`, `visible-at-c1`, `visible-at-c4`, `c1-isochrone` … `c4-isochrone` | Each of the three regions and the partial-transit footprint has its own fill; C1/C4 horizons and visible hemispheres are orange; the four contact isochrones each get a color | | Other limits, outlines, footprints and time markers | Other roles | Grey; polygons are outline-only by default | Zero-valued `Style` fields retain palette defaults: diff --git a/doc/manual/en/planets.md b/doc/manual/en/planets.md index f57f5c1..bed91bc 100644 --- a/doc/manual/en/planets.md +++ b/doc/manual/en/planets.md @@ -38,7 +38,7 @@ The inner planets (Mercury and Venus) add superior/inferior conjunction, greates - [Position and coordinates](#position-and-coordinates) - [Rise, set and culmination](#rise-set-and-culmination) - [Conjunctions, oppositions, stations and quadratures](#conjunctions-oppositions-stations-and-quadratures) - - [Greatest elongation and geocentric transits](#greatest-elongation-and-geocentric-transits) + - [Greatest elongation and transits](#greatest-elongation-and-transits) - [Nodes, phase, magnitude, diameter and parallactic angle](#nodes-phase-magnitude-diameter-and-parallactic-angle) - [Division of labour with other manuals](#division-of-labour-with-other-manuals) - [Physical ephemerides](#physical-ephemerides) @@ -188,22 +188,22 @@ fmt.Println(saturn.NextProgradeToRetrograde(date)) // Saturn's next prograde-to- Event searches return the nearest event at or after the input instant and keep its time zone; the paired `Last...` functions and the direction-agnostic `NextRetrograde` are in [Conjunctions, oppositions, stations and quadratures](#conjunctions-oppositions-stations-and-quadratures); -Mercury's and Venus's `NextGreatestElongationEast` / `...West` are in [Greatest elongation and geocentric transits](#greatest-elongation-and-geocentric-transits). +Mercury's and Venus's `NextGreatestElongationEast` / `...West` are in [Greatest elongation and transits](#greatest-elongation-and-transits). ### Mercury and Venus transits ```go -transit := mercury.NextTransit(date) // next geocentric Mercury transit after 2020 +transit := mercury.NextTransit(date) // next Mercury transit after 2020 fmt.Println(transit.Valid, transit.Start, transit.Greatest) // whether one was found, first contact, greatest transit fmt.Println(transit.Duration, transit.MinimumSeparationArcsec) // duration and minimum separation at greatest transit ``` ```text -true 2032-11-13 14:41:13.161198198 +0800 CST 2032-11-13 16:54:12.821315824 +0800 CST -4h26m2.695272267s 572.0643215495325 +true 2032-11-13 14:41:13.691189289 +0800 CST 2032-11-13 16:54:13.351347148 +0800 CST +4h26m2.484410405s 572.0831814505963 ``` -`TransitInfo.Valid` false means no transit inside the search window and every other field is the zero value; `TransitInfo` is the geocentric geometry. For the contact times seen from one site, `LocalTransit` returns the topocentric result: +`TransitInfo.Valid` false means no transit inside the search window and every other field is the zero value; otherwise `TransitInfo` describes a global event whose contact fields are geocentric reference times. Without geocentric contacts `HasExternal` is false, `Start` / `End` are zero, and only the `WindowStart` / `WindowEnd` search window is set. For the contact times seen from one site, `LocalTransit` returns the topocentric result: ```go transit := mercury.NextTransit(time.Date(2019, 1, 1, 0, 0, 0, 0, time.UTC)) @@ -214,12 +214,14 @@ fmt.Println(local.Start.Time, local.Start.SunAltitude) // topocentric C1 ```text true true true -2019-11-11 12:36:09.037762284 +0000 UTC 7.933256722170642 +2019-11-11 12:36:08.979182839 +0000 UTC 7.933092450591312 ``` Topocentric contacts are re-solved after applying parallax to the Sun and the planet separately; Venus contacts can shift by about ±7 minutes from the geocentric ones across the globe, Mercury by about ±2 minutes. `SunAltitude` is the geometric altitude of the Sun's center; `Visible`, `Sunrise` and `Sunset` test whether the Sun's center is above the apparent horizon with standard refraction (1010 hPa, 10 °C), a geometric threshold of about −34′. A site may see only a grazing transit even when the geocentric one has internal contacts (the 1999-11-15 Mercury transit at high southern latitudes); then `HasInternal` is false and `InternalStart` / `InternalEnd` are zero. Global maps and topocentric disk charts are in [Eclipses and transits](eclipse.md#transit-charts). -The four contacts, partial transits and the internal-contact fields are in [Greatest elongation and geocentric transits](#greatest-elongation-and-geocentric-transits). +The four contacts, partial transits and the internal-contact fields are in [Greatest elongation and transits](#greatest-elongation-and-transits). + +`Gamma` is the signed distance from Earth's centre to the shadow axis at greatest (positive to the north) and `PenumbralRadius` is the penumbral radius on the fundamental plane, both in Earth equatorial radii. The transit penumbra is far larger than Earth, so the solar-eclipse test `|Gamma| < 1` does not apply. `Series` is the transit family: Venus families average 243 years and are numbered for astronomical years −2000 to +4000, Mercury families average 46 years and are numbered for +1601 to +2300; outside those spans, or for an unnumbered family, `Number` is zero while `PeriodYears` still gives the mean period. ### Phase, diameter, magnitude and nodes @@ -580,15 +582,15 @@ fmt.Println(saturn.LastProgradeToRetrograde(date), saturn.NextProgradeToRetrogra fmt.Println(saturn.LastRetrogradeToPrograde(date), saturn.NextRetrogradeToPrograde(date)) ``` -### Greatest elongation and geocentric transits +### Greatest elongation and transits Greatest elongation only makes sense for Mercury and Venus: `LastGreatestElongation` / `NextGreatestElongation` ignore the side, while `LastGreatestElongationEast` / `NextGreatestElongationEast` and the `...West` forms distinguish it. -The public entry points for geocentric transits are likewise only in `mercury` / `venus`: `LastTransit` / `NextTransit` / `ClosestTransit` return a `TransitInfo`. +The public transit entry points are likewise only in `mercury` / `venus`: `LastTransit` / `NextTransit` / `ClosestTransit` return a `TransitInfo`. -`Valid == false` means no transit inside the search window and every other field is a zero value; a partial transit has no internal contacts, so `HasInternal` is false and `InternalStart` / `InternalEnd` are zero. +`Valid == false` means no transit inside the search window and every other field is a zero value. Without geocentric contacts `HasExternal` is false and `Start` / `End` are zero; the event is then accepted or rejected by topocentric surface visibility. Without a topocentric internal contact `HasInternal` is false and `InternalStart` / `InternalEnd` are zero. -A transit is a purely geocentric geometry test and does not check whether the Sun is above the horizon at a given site. +Contact fields come from the geocentric geometry, but whether a site can see the transit, and which contacts it sees, is decided by `LocalTransit`. ```go // Greatest elongation: use ...East / ...West to pick a side, or the plain form for either. @@ -902,7 +904,7 @@ Geocentric, topocentric, and distance are three different conventions: `Apparent ### Zero values, out-of-range input and sentinel errors -When an event search finds nothing inside its window it returns zero values or a struct with `Valid == false` (`TransitInfo.Valid`, `GalileanPhenomenonEvent.Valid`, `GalileanPhenomenonContactEvent.Valid`) instead of an error; when `TransitInfo.HasInternal` is false, `InternalStart` / `InternalEnd` are zero values because a partial transit has no internal contacts. +When an event search finds nothing inside its window it returns zero values or a struct with `Valid == false` (`TransitInfo.Valid`, `GalileanPhenomenonEvent.Valid`, `GalileanPhenomenonContactEvent.Valid`) instead of an error; when `TransitInfo.HasExternal` is false, `Start` / `End` are zero values, and when `TransitInfo.HasInternal` is false, `InternalStart` / `InternalEnd` are zero values. `planet.WherePlanet` / `planet.WherePlanetN` return `NaN` instead of panicking when `xt` / `zn` is out of range. `n < 0` selects the full series; all nonnegative `n` select one semianalytic fast table, while Legacy truncates according to `n`. diff --git a/doc/manual/map-geojson.md b/doc/manual/map-geojson.md index 732cd2e..b134df6 100644 --- a/doc/manual/map-geojson.md +++ b/doc/manual/map-geojson.md @@ -73,7 +73,11 @@ func main() { ## 导出凌日 GeoJSON -凌日入口接收地心 `TransitInfo`，分区口径与凌日全球图一致： +凌日入口接收全球事件的 `TransitInfo`，接触字段仍为地心参考时刻，分区口径与凌日全球图一致： + +地心无接触的边缘事件（如 2854-12-14 金星凌日）导出站心实际可见的 `transit-partial` 区域， +不生成不存在的地心 C1/C4 地平线。属性还包含 `has_geocentric_contacts`、`gamma`、`penumbral_radius`、 +`series`、`series_period_years`。足迹按 2° 经度、1° 纬度采样并细化边界，这是制图采样而不是接触求解精度，比采样间距更窄的可见区可能未检出；精确判断请用 `LocalTransit`。 ```go m := mercury.NextTransit(time.Date(2026, 10, 12, 0, 0, 0, 0, cst)) @@ -190,7 +194,7 @@ true true 4217 - 日食时刻与食分线：`rise-set-phase-lines`、`magnitude-contours`、`greatest-time-isochrones` - 瞬时轮廓：`penumbral-outlines`、`central-shadow-outlines` - 月食：`p1-p4-visibility-regions`、`p1-p4-horizon-boundaries` -- 凌日：`c1-c4-visibility-regions`、`c1-c4-horizon-boundaries`、`transit-contact-isochrones` +- 凌日：`c1-c4-visibility-regions`、`c1-c4-horizon-boundaries`、`transit-contact-isochrones`；地心无接触时可见区填充改用 `topocentric-transit-visibility` ## 地图投影 @@ -294,7 +298,7 @@ func main() { - `MarshalLunarEclipseWithOptions`：用 `LunarEclipseOptions` 一次给出时间标记、`SkipRoles` 与时间包络的采样参数。`SkipRoles` 可跳过时间包络和仅见半影带；跳过两个时间包络时不会执行包络采样。`EnvelopeSweepSamples` 默认 48，限制在 `[2, 192]`；`EnvelopeLongitudePoints` 默认 `max(360, boundaryPoints)`，限制在 `[12, 720]`。 - `MarshalStarOccultation` / `MarshalStarOccultationWithTimeMarkers` - `MarshalPlanetOccultation` / `MarshalPlanetOccultationWithTimeMarkers` -- `MarshalMercuryTransit` / `MarshalVenusTransit`：接收 `mercury.TransitInfo` / `venus.TransitInfo`，分区口径同凌日 SVG 全球图（一触、凌甚、四触太阳中心是否在视地平上，计入标准大气折射）。`BoundaryPoints` 是地平线采样点数，同时决定分区的经度列数 `max(360, BoundaryPoints)`（上限 720）。role 有 `transit-entire`、`transit-at-sunrise`、`transit-at-sunset` 三个分区，`visible-at-c1` / `visible-at-c4` 与 `c1-horizon` / `c4-horizon`，`greatest`（太阳直射点）与时间标记 `time-marker`；`IsochroneStep` 为正（至少 1 秒）时追加 `c1-isochrone` 到 `c4-isochrone` 站心等时线（掠凌无二、三触；时刻是站心接触时刻，受视差影响可落在地心一触到四触之外），某个接触的等时线不落在可见半球上时该 role 不出现；`SkipRoles` 跳过全部三个分区时不做分区计算。 +- `MarshalMercuryTransit` / `MarshalVenusTransit`：接收 `mercury.TransitInfo` / `venus.TransitInfo`，分区口径同凌日 SVG 全球图（一触、凌甚、四触太阳中心是否在视地平上，计入标准大气折射）。`BoundaryPoints` 是地平线采样点数，同时决定分区的经度列数 `max(360, BoundaryPoints)`（上限 720）。role 有 `transit-entire`、`transit-at-sunrise`、`transit-at-sunset` 三个分区，`visible-at-c1` / `visible-at-c4` 与 `c1-horizon` / `c4-horizon`，`greatest`（太阳直射点）与时间标记 `subsolar-track`；地心无接触事件改用 `transit-partial` 分区，不生成三个可见分区、视地平与可见半球 role；`IsochroneStep` 为正（至少 1 秒）时追加 `c1-isochrone` 到 `c4-isochrone` 站心等时线（地表没有站心内切接触时无二、三触，地心只有外切的掠凌仍可能在部分地区出现；时刻是站心接触时刻，受视差影响可落在地心一触到四触之外），某个接触的等时线不落在可见半球上时该 role 不出现；`SkipRoles` 跳过全部三个分区时不做分区计算。 - `MarshalSolarEclipseWithOptions`：用 `SolarEclipseOptions` 一次给出时间标记与 `SkipRoles`，等价于上面的日食入口再加图层裁剪。`SkipRoles` 里最常用的是 `partial-footprint`（瞬时半影轮廓）；`partial-band` 是日食带边缘（食分 0 界限），通常要留着。 ### 图层筛选 @@ -434,7 +438,7 @@ func main() { | 掩星全掩与偏掩带 | `total-band`、`partial-band`、`total-band-outline`、`occultation-band` | 黄色，带半透明填充 | | 升落可见边界、月食时间包络、日食偏食带边缘 | `visibility-boundary`、`p1-horizon`、`p4-horizon`、`visible-at-p1`、`visible-at-p4`、`visible-during-eclipse`、`visible-throughout-eclipse`、日食的 `partial-band` | 橙色 | | 仅见半影的月出/月落带 | `penumbra-moonrise`、`penumbra-moonset` | 月出蓝紫、月落紫红，带半透明填充 | -| 凌日可见分区与接触线 | `transit-entire`、`transit-at-sunrise`、`transit-at-sunset`、`c1-horizon`、`c4-horizon`、`visible-at-c1`、`visible-at-c4`、`c1-isochrone` … `c4-isochrone` | 三个分区各有填充色；C1/C4 视地平与可见半球橙色；四触等时线各一色 | +| 凌日可见分区与接触线 | `transit-entire`、`transit-at-sunrise`、`transit-at-sunset`、`transit-partial`、`c1-horizon`、`c4-horizon`、`visible-at-c1`、`visible-at-c4`、`c1-isochrone` … `c4-isochrone` | 三个分区与部分凌日可见区各有填充色；C1/C4 视地平与可见半球橙色；四触等时线各一色 | | 其他限界、轮廓、足迹与时间标记 | 其他 `role` | 灰色，面默认只画轮廓 | `Style` 的零值字段沿用默认样式： diff --git a/doc/manual/planets.md b/doc/manual/planets.md index fcb0d92..3b1e360 100644 --- a/doc/manual/planets.md +++ b/doc/manual/planets.md @@ -32,7 +32,7 @@ - [位置与坐标](#位置与坐标) - [升落与中天](#升落与中天) - [合冲留与方照](#合冲留与方照) - - [大距与地心凌日](#大距与地心凌日) + - [大距与凌日](#大距与凌日) - [节点、相位、视星等、视直径与视差角](#节点相位视星等视直径与视差角) - [与其他手册的分工](#与其他手册的分工) - [物理星历](#物理星历) @@ -174,22 +174,22 @@ fmt.Println(saturn.NextProgradeToRetrograde(date)) // 土星下次顺转逆的 2020-05-11 17:26:53.961271941 +0800 CST ``` -事件搜索取"当前或之后最近一次"并保持输入时区，成对的 `Last...` 与不区分方向的 `NextRetrograde` 见[合冲留与方照](#合冲留与方照)；水星、金星的 `NextGreatestElongationEast` / `...West` 见[大距与地心凌日](#大距与地心凌日)。 +事件搜索取"当前或之后最近一次"并保持输入时区，成对的 `Last...` 与不区分方向的 `NextRetrograde` 见[合冲留与方照](#合冲留与方照)；水星、金星的 `NextGreatestElongationEast` / `...West` 见[大距与凌日](#大距与凌日)。 ### 水星与金星凌日 ```go -transit := mercury.NextTransit(date) // 2020 年之后下一场地心水星凌日 +transit := mercury.NextTransit(date) // 2020 年之后下一场水星凌日 fmt.Println(transit.Valid, transit.Start, transit.Greatest) // 是否有凌日、一触与凌甚 fmt.Println(transit.Duration, transit.MinimumSeparationArcsec) // 历时与凌甚最小角距 ``` ```text -true 2032-11-13 14:41:13.161198198 +0800 CST 2032-11-13 16:54:12.821315824 +0800 CST -4h26m2.695272267s 572.0643215495325 +true 2032-11-13 14:41:13.691189289 +0800 CST 2032-11-13 16:54:13.351347148 +0800 CST +4h26m2.484410405s 572.0831814505963 ``` -`TransitInfo.Valid` 为假表示搜索窗口内没有凌日、其余字段是零值；`TransitInfo` 是地心几何结果。某个地点看到的接触时刻用 `LocalTransit` 求站心结果： +`TransitInfo.Valid` 为假表示搜索窗口内没有凌日、其余字段是零值；其余情况下 `TransitInfo` 表示全球事件，接触字段是地心参考时刻。地心无接触时 `HasExternal` 为假，`Start` / `End` 也是零值，只有 `WindowStart` / `WindowEnd` 站心求解窗口有值。某个地点看到的接触时刻用 `LocalTransit` 求站心结果： ```go transit := mercury.NextTransit(time.Date(2019, 1, 1, 0, 0, 0, 0, time.UTC)) @@ -200,15 +200,17 @@ fmt.Println(local.Start.Time, local.Start.SunAltitude) // 站心一触 ```text true true true -2019-11-11 12:36:09.037762284 +0000 UTC 7.933256722170642 +2019-11-11 12:36:08.979182839 +0000 UTC 7.933092450591312 ``` 站心接触对太阳与行星分别做视差修正后重解，所以各地接触时刻与地心值不同：金星可达约 ±7 分钟，水星约 ±2 分钟。 -- `SunAltitude` 是日心几何高度。 +- `SunAltitude` 是太阳中心的几何高度，不含折射。 - `Visible`、`Sunrise`、`Sunset` 按太阳中心是否在视地平上判定，计入标准大气（1010 hPa、10 °C）折射，几何阈值约 −34′。 - 地心有内切接触、站心可能只是掠凌（例如 1999-11-15 水星凌日在南半球高纬），此时 `HasInternal` 为假，`InternalStart` / `InternalEnd` 是零值。 -- 全球可见图与站心日面图见[日月食与凌日](eclipse.md#凌日出图)；四触、偏凌与内切口径见[大距与地心凌日](#大距与地心凌日)。 +- 全球可见图与站心日面图见[日月食与凌日](eclipse.md#凌日出图)；四触、偏凌与内切口径见[大距与凌日](#大距与凌日)。 + +`Gamma` 是凌甚时影轴到地心的有符号距离（北侧为正），`PenumbralRadius` 是地心基本面上的半影半径，两者都以地球赤道半径为单位。凌日半影远大于地球，不要套用日食的 `|Gamma| < 1` 判断。`Series` 是周期族：金星平均 243 年、编号适用于天文年 −2000 至 +4000，水星平均 46 年、适用于 +1601 至 +2300；超出该范围或不属于已编号族时 `Number` 为 0，`PeriodYears` 仍给出平均周期。 ### 相位、视直径、视星等与节点 @@ -348,7 +350,7 @@ fmt.Println(venus.AscendingNode(date), venus.DescendingNode(date)) 76.86008484515058 256.8600848451506 // 金星升交点、降交点黄经，单位度 ``` -水星和金星还提供 `NextTransit` / `LastTransit` / `ClosestTransit` 地心凌日查询。这里的“地心”指从地球中心看到的行星圆面经过太阳圆面，不判断某个地点当时太阳是否在地平线上；观测计划用 `LocalTransit` 求站心接触与当时太阳高度（见[水星与金星凌日](#水星与金星凌日)）。 +水星和金星还提供 `NextTransit` / `LastTransit` / `ClosestTransit` 凌日查询。搜索返回全球事件，接触字段是地心参考时刻；地心无接触、只有部分地区可见的掠凌也会返回，这类事件是否成立按地表站心可见性判定。某个地点当时太阳是否在地平线上、当地看到哪些接触，用 `LocalTransit` 求站心结果（见[水星与金星凌日](#水星与金星凌日)）。 ```go package main @@ -553,11 +555,11 @@ fmt.Println(saturn.LastProgradeToRetrograde(date), saturn.NextProgradeToRetrogra fmt.Println(saturn.LastRetrogradeToPrograde(date), saturn.NextRetrogradeToPrograde(date)) ``` -### 大距与地心凌日 +### 大距与凌日 大距只对水星、金星有意义：`LastGreatestElongation` / `NextGreatestElongation` 不区分东西，`LastGreatestElongationEast` / `NextGreatestElongationEast` 与 `...West` 区分。 -地心凌日的公开入口同样只在 `mercury` / `venus`：`LastTransit` / `NextTransit` / `ClosestTransit` 返回 `TransitInfo`。`Valid` 为假表示搜索窗口内没有凌日，其余字段都是零值；偏凌日没有内切，此时 `HasInternal` 为假、`InternalStart` / `InternalEnd` 是零值。凌日只判断地心几何，不判断某地当时太阳是否在地平线上。 +凌日的公开入口同样只在 `mercury` / `venus`：`LastTransit` / `NextTransit` / `ClosestTransit` 返回 `TransitInfo`。`Valid` 为假表示搜索窗口内没有凌日，其余字段都是零值。地心无接触时 `HasExternal` 为假、`Start` / `End` 为零值，事件是否成立按地表站心可见性判定；站心内切不存在时 `HasInternal` 为假、`InternalStart` / `InternalEnd` 是零值。某个地点是否可见、看到哪些接触用 `LocalTransit` 判定。 ```go // 大距：区分东西用 ...East / ...West，不区分时用 NextGreatestElongation。 @@ -566,7 +568,7 @@ fmt.Println(venus.NextGreatestElongation(date)) ``` ```go -// 地心凌日：Valid 为假时其余字段为零值。 +// 凌日：Valid 为假时其余字段为零值。 transit := mercury.NextTransit(date) if transit.Valid { fmt.Println(transit.Start, transit.InternalStart, transit.Greatest, transit.InternalEnd, transit.End) @@ -859,7 +861,7 @@ fmt.Println(planet.WherePlanet(-1, 0, 2460000.5), planet.WherePlanet(4, 0, 24600 ### 零值、越界与哨兵错误 -事件搜索在窗口内没有事件时返回零值或 `Valid=false` 的结构（`TransitInfo.Valid`、`GalileanPhenomenonEvent.Valid`、`GalileanPhenomenonContactEvent.Valid`），不返回错误；`TransitInfo` 在 `HasInternal` 为假时 `InternalStart` / `InternalEnd` 为零值，偏凌日没有内切。 +事件搜索在窗口内没有事件时返回零值或 `Valid=false` 的结构（`TransitInfo.Valid`、`GalileanPhenomenonEvent.Valid`、`GalileanPhenomenonContactEvent.Valid`），不返回错误；`TransitInfo` 在 `HasExternal` 为假时 `Start` / `End` 为零值，`HasInternal` 为假时 `InternalStart` / `InternalEnd` 为零值。 `planet.WherePlanet` / `planet.WherePlanetN` 的 `xt` / `zn` 越界返回 `NaN` 而不 panic；`...N` 的 `n < 0` 用完整级数，半解析的非负 `n` 用同一张快速表，Legacy 按 `n` 截断。`RiseTime` / `SetTime` / `DownTime`（含 `N` 版）当天没有几何升/落时返回下面的哨兵错误，时刻为零值。 diff --git a/eclipse/svg/data_source_test.go b/eclipse/svg/data_source_test.go index c8f550f..c92d797 100644 --- a/eclipse/svg/data_source_test.go +++ b/eclipse/svg/data_source_test.go @@ -18,6 +18,7 @@ func TestEclipseMapDataSourcesUseDocumentedVocabulary(t *testing.T) { "central-path-limits", "penumbral-outlines", "central-shadow-outlines", "p1-p4-visibility-regions", "p1-p4-horizon-boundaries", "c1-c4-visibility-regions", "c1-c4-horizon-boundaries", "transit-contact-isochrones", + "topocentric-transit-visibility", } if len(documented) != len(eclipseMapDataSources) { t.Fatalf("vocabulary has %d entries, want %d", len(eclipseMapDataSources), len(documented)) diff --git a/eclipse/svg/doc.go b/eclipse/svg/doc.go index 81ee315..4143474 100644 --- a/eclipse/svg/doc.go +++ b/eclipse/svg/doc.go @@ -19,6 +19,7 @@ // c1-c4-visibility-regions 凌日可见性分区 // c1-c4-horizon-boundaries 凌日 C1、C4 视地平边界 // transit-contact-isochrones 凌日站心接触等时线 +// topocentric-transit-visibility 地心无接触时按站心求得的部分凌日可见区 // // Package svg renders solar-eclipse, lunar-eclipse and Mercury/Venus transit SVG charts. Degradable layers tag the geometry // source they actually used in data-source, using exactly the values above. @@ -43,6 +44,7 @@ const ( eclipseMapSourceTransitVisibilityRegions = "c1-c4-visibility-regions" eclipseMapSourceTransitHorizonBoundaries = "c1-c4-horizon-boundaries" eclipseMapSourceTransitContactIsochrones = "transit-contact-isochrones" + eclipseMapSourceTransitPartialVisibility = "topocentric-transit-visibility" ) // eclipseMapDataSources 是全部合法取值，测试用它封闭词表。 @@ -64,4 +66,5 @@ var eclipseMapDataSources = []string{ eclipseMapSourceTransitVisibilityRegions, eclipseMapSourceTransitHorizonBoundaries, eclipseMapSourceTransitContactIsochrones, + eclipseMapSourceTransitPartialVisibility, } diff --git a/eclipse/svg/transit_local.go b/eclipse/svg/transit_local.go index 4b9b766..930248e 100644 --- a/eclipse/svg/transit_local.go +++ b/eclipse/svg/transit_local.go @@ -175,13 +175,31 @@ func renderTransitLocalSVG( fmt.Fprintf(&b, `%s`, width/2, html.EscapeString(titleText)) for index, line := range []string{site, summary} { + line = svgchart.EllipsizeText(line, width-80, 13) fmt.Fprintf(&b, `%s`, width/2, 72+float64(index)*19, html.EscapeString(line)) } - // 日面区占左侧，接触表占右侧。 + // 日面区占左侧，接触表占右侧；日面半径要容下盘外的东南西北轴标签，所以上下界由副标题与页脚算出。 panelWidth := 300.0 - top, bottom := 112.0, height-86 + footerLines := transitLocalFooterLines(local, options, width) + footerBaselines := transitLocalFooterBaselines(len(footerLines), height) + footerTop := height - 30.0 + if len(footerBaselines) > 0 { + footerTop = footerBaselines[0] - 10.2 + } + const ( + headerBottom = 72 + 19 + 10 // 经纬度行、摘要行与留白 + axisLabelAbove = 33.0 // 盘顶到“北”字上沿 + axisLabelBelow = 36.0 // 盘底到“南”字下沿 + minimumDiskHeight = 150.0 + ) + top, bottom := headerBottom+axisLabelAbove, footerTop-axisLabelBelow + if bottom-top < minimumDiskHeight { + shrink := (minimumDiskHeight - (bottom - top)) / 2 + top -= shrink + bottom += shrink + } diskLeft, diskRight := 46.0, width-panelWidth-60 cx, cy := (diskLeft+diskRight)/2, (top+bottom)/2 extent := local.SunSemidiameterArcsec + 2.2*local.PlanetSemidiameterArcsec + 60 @@ -196,6 +214,9 @@ func renderTransitLocalSVG( events := transitLocalEvents(local, options.Language) planetRadius := math.Max(local.PlanetSemidiameterArcsec*scale, 2) table := &svgchart.LabelTable{} + table.Reserve(30, footerTop-4, width-60, height-footerTop+4) + table.Reserve(cx-14, cy+sunRadius+16, 28, 18) + table.Reserve(cx-14, cy-sunRadius-32, 28, 20) for _, event := range events { x, y := transitLocalEventPosition(frames, event.contact.JD, mapX, mapY) fill, stroke := "#1d1d1d", "#000000" @@ -234,7 +255,7 @@ func renderTransitLocalSVG( // 标签平移进白色图框内、接触表左侧，不能伸出画面。 box := svgchart.LabelTextBox(candidate.X, candidate.Y, 11, text, anchor) candidate.X += math.Max(0, 30-box.X) - math.Max(0, box.X+box.Width-(width-panelWidth-44)) - candidate.Y += math.Max(0, top-8-box.Y) - math.Max(0, box.Y+box.Height-(height-30)) + candidate.Y += math.Max(0, top-8-box.Y) - math.Max(0, box.Y+box.Height-(footerTop-4)) candidates = append(candidates, candidate) } placement, placed = table.Place(text, 11, candidates) @@ -254,8 +275,14 @@ func renderTransitLocalSVG( placement.X, placement.Y, placement.Anchor, html.EscapeString(text)) } - writeTransitLocalContactTable(&b, local, events, at, width-panelWidth-36, top+10, panelWidth, options.Language) - writeTransitLocalFooter(&b, local, options, width, height) + tableWidth := panelWidth + tableX := width - panelWidth - 36 + if minimum := cx + sunRadius + 40; tableX < minimum { + tableX = minimum + tableWidth = width - tableX - 36 + } + writeTransitLocalContactTable(&b, local, events, at, tableX, top+10, tableWidth, options.Language) + writeTransitLocalFooter(&b, footerLines, footerBaselines) b.WriteString(``) return b.String() } @@ -301,7 +328,7 @@ func transitLocalSummaryText(local basic.PlanetTransitLocalResult, language stri // writeTransitLocalAxes 画过日心的南北、东西虚线与方向字。 func writeTransitLocalAxes(b *strings.Builder, cx, cy, radius float64, language string) { - reach := radius + 18 + reach := radius + 12 fmt.Fprintf(b, ``, cx, cy-reach, cx, cy+reach) fmt.Fprintf(b, ``, cx-reach, cy, cx+reach, cy) labels := [4]string{"北", "南", "东", "西"} @@ -419,18 +446,25 @@ func writeTransitLocalContactTable( } } -func writeTransitLocalFooter(b *strings.Builder, local basic.PlanetTransitLocalResult, options TransitLocalSVGOptions, width, height float64) { +func transitLocalFooterLines(local basic.PlanetTransitLocalResult, options TransitLocalSVGOptions, width float64) []string { note := options.FooterNote if note == "" { if options.Language == "en" { note = "Sky view, north up and east left. Topocentric contacts; Sun altitudes are geometric centre altitudes, grey parts have the Sun's centre below the apparent horizon (standard refraction)." } else { - note = "天空视向，上北左东；站心接触时刻，太阳高度为日心几何高度，灰色部分太阳中心在视地平下（计入标准大气折射）。" + note = "天空视向，上北左东；站心接触时刻，太阳高度为太阳中心几何高度，灰色部分太阳中心在视地平下（计入标准大气折射）。" } } scaleNote := timenote.Scale(options.TimeScale, basic.JD2DateByZone(local.Greatest.JD, time.UTC, false), options.Location, options.Language) - lines := svgchart.FooterLinesWithScale(svgchart.WrapText(note, width-80, 12), scaleNote, width-80, 12, 3) - baselines, _ := svgchart.FooterBlock(len(lines), height, 12, svgchart.FooterLineHeight(12), svgchart.FooterBottomPadding(12)+18) + return svgchart.FooterLinesWithScale(svgchart.WrapText(note, width-80, 12), scaleNote, width-80, 12, 3) +} + +func transitLocalFooterBaselines(lineCount int, height float64) []float64 { + baselines, _ := svgchart.FooterBlock(lineCount, height, 12, svgchart.FooterLineHeight(12), svgchart.FooterBottomPadding(12)+18) + return baselines +} + +func writeTransitLocalFooter(b *strings.Builder, lines []string, baselines []float64) { for index, line := range lines { fmt.Fprintf(b, `%s`, baselines[index], html.EscapeString(line)) diff --git a/eclipse/svg/transit_map.go b/eclipse/svg/transit_map.go index 8306a3e..075d043 100644 --- a/eclipse/svg/transit_map.go +++ b/eclipse/svg/transit_map.go @@ -36,9 +36,9 @@ type TransitIsochrones uint8 const ( // TransitIsochroneC1 画一触等时线 / draw C1 isochrones. TransitIsochroneC1 TransitIsochrones = 1 << iota - // TransitIsochroneC2 画二触等时线；掠凌没有内切接触时忽略 / draw C2 isochrones, ignored for a grazing transit. + // TransitIsochroneC2 画二触等时线；地表没有站心内切时忽略 / draw C2 isochrones, ignored where no site has internal contact. TransitIsochroneC2 - // TransitIsochroneC3 画三触等时线；掠凌没有内切接触时忽略 / draw C3 isochrones, ignored for a grazing transit. + // TransitIsochroneC3 画三触等时线；地表没有站心内切时忽略 / draw C3 isochrones, ignored where no site has internal contact. TransitIsochroneC3 // TransitIsochroneC4 画四触等时线 / draw C4 isochrones. TransitIsochroneC4 @@ -64,8 +64,7 @@ type TransitMapSVGOptions struct { // TimeScale 选择图中时刻的时标：零值 UTC；TimeScaleUT1 要求 Location 为 nil 或 UTC，否则返回 false。 // TimeScale selects the label scale: the zero value is UTC; TimeScaleUT1 requires Location to be nil or UTC, otherwise the renderer returns false. TimeScale astro.TimeScale - // Projection 选择地图投影；零值使用等经纬投影，正射投影视点取凌甚时的太阳直射点。 - // Projection selects the map projection. The zero value uses equirectangular; orthographic centres on the subsolar point at greatest. + // Projection 选择地图投影，默认等经纬；正射视点在地心无接触时取站心可见区中心，否则取凌甚太阳直射点 / projection, default equirectangular; orthographic centres on the visible footprint without geocentric contacts, otherwise the subsolar point at greatest. Projection EclipseMapProjection // Isochrones 选择要画的站心接触等时线，可按位组合；零值不画（默认关闭，需显式开启）。 // Isochrones selects the topocentric contact isochrones to draw and may combine flags. The zero value draws none (off by default). @@ -79,12 +78,14 @@ type TransitMapSVGOptions struct { FooterNote string } -// TransitMapSVG 绘制距 date 最近一次水星或金星凌日的全球可见图：按一触、凌甚、四触太阳中心在视地平上 -// （计入标准大气折射）的区域分区；Options.Isochrones 显式开启时叠加所选接触的站心等时线。找不到凌日或选项无效时返回 false。 -// TransitMapSVG renders the global visibility of the Mercury or Venus transit closest to date: regions follow where the -// Sun's centre is above the apparent horizon (standard refraction) at C1, greatest and C4. Topocentric isochrones of the -// selected contacts are overlaid only when Options.Isochrones enables them. It returns false when no transit is found -// or the options are invalid. +// TransitMapSVG 绘制距 date 最近一次水星或金星凌日的全球可见图：普通凌日按一触、凌甚、四触太阳中心在视地平上 +// （计入标准大气折射）的区域分区，地心无接触时绘制站心部分凌日可见区；Options.Isochrones 显式开启时叠加所选接触的 +// 站心等时线，图例只列实际画出的接触。找不到凌日或选项无效时返回 false。 +// TransitMapSVG renders the global visibility of the Mercury or Venus transit closest to date: ordinary transits are +// partitioned by whether the Sun's centre is above the apparent horizon (standard refraction) at C1, greatest and C4; +// events without geocentric contacts show the topocentric partial-transit footprint. Topocentric isochrones of the +// selected contacts are overlaid only when Options.Isochrones enables them, and the legend lists only contacts that +// produced lines. It returns false when no transit is found or the options are invalid. func TransitMapSVG(planet TransitPlanet, date time.Time, options TransitMapSVGOptions) (string, bool) { transit, ok := closestPlanetTransit(planet, date) if !ok || !validEclipseMapProjection(options.Projection) { @@ -97,7 +98,7 @@ func TransitMapSVG(planet TransitPlanet, date time.Time, options TransitMapSVGOp options.Location = time.UTC } options = normalizeTransitMapSVGOptions(date, options) - if !transit.HasInternal { + if !transit.HasExternal { options.Isochrones &^= TransitIsochroneC2 | TransitIsochroneC3 } options.Isochrones &= TransitIsochroneAll @@ -114,7 +115,7 @@ func TransitMapSVG(planet TransitPlanet, date time.Time, options TransitMapSVGOp return renderTransitMapSVG(planet, transit, options, projection), true } -// closestPlanetTransit 返回距 date 最近的一次地心凌日；行星无效或找不到时第二个返回值为 false。 +// closestPlanetTransit 返回距 date 最近的一次全球凌日。 func closestPlanetTransit(planet TransitPlanet, date time.Time) (basic.PlanetTransitResult, bool) { jd := basic.Date2JD(date.UTC()) var transit basic.PlanetTransitResult @@ -155,7 +156,12 @@ type transitMapContacts struct { } func transitMapContactsOf(transit basic.PlanetTransitResult) transitMapContacts { - at := func(jd float64) time.Time { return basic.JD2DateByZone(jd, time.UTC, false) } + at := func(jd float64) time.Time { + if jd == 0 { + return time.Time{} + } + return basic.JD2DateByZone(jd, time.UTC, false) + } contacts := transitMapContacts{ start: at(transit.ExternalIngress), greatest: at(transit.Greatest), @@ -178,7 +184,7 @@ func renderTransitMapSVG( // 底部比月食图多留两行：图例下要放可折行的说明与时标声明。 frame := eclipseMapFrame(options.Width, options.Height, projection, 142, 122) if projection == svgmap.ProjectionOrthographic { - center := solarEclipseSubsolarPoint(contacts.greatest) + center := transitMapOrthographicCenter(transit) frame.CenterLongitude, frame.CenterLatitude = center.Longitude, center.Latitude } title := options.Title @@ -192,8 +198,12 @@ func renderTransitMapSVG( } var b strings.Builder - fmt.Fprintf(&b, ``, + fmt.Fprintf(&b, ``) b.WriteString(``) writeTransitMapDefinitions(&b, frame, contacts) b.WriteString(``) @@ -211,13 +221,15 @@ func renderTransitMapSVG( frame.WriteOcean(&b) frame.WriteGraticule(&b, "transit-map-clip") frame.WriteLand(&b, "transit-map-clip") - writeTransitMapRegions(&b, frame, contacts) - if options.Isochrones != 0 { - writeTransitMapIsochrones(&b, frame, transit, options) + if transit.HasExternal { + writeTransitMapRegions(&b, frame, contacts) + } else { + writeTransitMapPartialRegions(&b, frame, transit) } + options.Isochrones = writeTransitMapIsochrones(&b, frame, transit, options) frame.WriteFrame(&b) - legendRows := writeTransitMapLegend(&b, frame, float64(options.Width), options) - writeTransitMapFooter(&b, frame, options, projection, contacts.greatest, legendRows) + legendRows := writeTransitMapLegend(&b, frame, float64(options.Width), options, transit.HasExternal) + writeTransitMapFooter(&b, frame, options, projection, contacts.greatest, legendRows, transit.HasExternal) b.WriteString(``) return b.String() } @@ -264,6 +276,9 @@ const ( // writeTransitMapDefinitions 输出一触、凌甚、四触三个太阳在地平上半球的形状与排除掩膜，分区口径同月食可见图。 func writeTransitMapDefinitions(b *strings.Builder, frame svgmap.Frame, contacts transitMapContacts) { b.WriteString(frame.ClipDefinition("transit-map-clip")) + if contacts.start.IsZero() || contacts.end.IsZero() { + return + } startPath, _ := transitMapSunUpPath(contacts.start, frame) endPath, _ := transitMapSunUpPath(contacts.end, frame) greatestPath, _ := transitMapSunUpPath(contacts.greatest, frame) @@ -285,10 +300,19 @@ const ( transitMapColorSunrise = "#4e9da0" transitMapColorSunset = "#e2aa4b" transitMapColorNone = "#747b7d" + transitMapColorPartial = "#f7dd93" transitMapColorC1 = "#a3392f" transitMapColorC2 = "#c0702a" transitMapColorC3 = "#7a4ca0" transitMapColorC4 = "#2f5ea3" + + transitMapFooterSideInset = 120.0 + + transitMapRegionFillOpacity = 0.34 + // 部分凌日只有一块可见区，要比普通凌日分区更浓才能从底图上辨认。 + transitMapPartialFillOpacity = 0.55 + // 图例色块直接画在白底上，比叠在海陆底色上的分区取更浓的值，同月食图。 + transitMapLegendFillOpacity = 0.78 ) // writeTransitMapRegions 着色四类可见区：全程可见、日出时凌日已开始、日落时凌日未结束、不可见。 @@ -298,25 +322,27 @@ func writeTransitMapRegions(b *strings.Builder, frame svgmap.Frame, contacts tra _, endBoundary := transitMapSunUpPath(contacts.end, frame) fmt.Fprintf(b, ``, eclipseMapSourceTransitVisibilityRegions) if frame.IsPolar() { - fmt.Fprintf(b, ``, - frame.X+frame.Width/2, frame.Y+frame.Height/2, frame.Width/2, transitMapColorNone) + fmt.Fprintf(b, ``, + frame.X+frame.Width/2, frame.Y+frame.Height/2, frame.Width/2, transitMapColorNone, transitMapRegionFillOpacity) } else { - fmt.Fprintf(b, ``, - frame.X, frame.Y, frame.Width, frame.Height, transitMapColorNone) + fmt.Fprintf(b, ``, + frame.X, frame.Y, frame.Width, frame.Height, transitMapColorNone, transitMapRegionFillOpacity) } - fmt.Fprintf(b, ``, transitMapColorSunset) - fmt.Fprintf(b, ``, transitMapColorSunrise) - fmt.Fprintf(b, ``, transitMapColorEntire) - fmt.Fprintf(b, ``, transitMapColorSunset) - fmt.Fprintf(b, ``, transitMapColorSunset) + fmt.Fprintf(b, ``, transitMapColorSunset, transitMapRegionFillOpacity) + fmt.Fprintf(b, ``, transitMapColorSunrise, transitMapRegionFillOpacity) + fmt.Fprintf(b, ``, transitMapColorEntire, transitMapRegionFillOpacity) + fmt.Fprintf(b, ``, transitMapColorSunset, transitMapRegionFillOpacity) + fmt.Fprintf(b, ``, transitMapColorSunset, transitMapRegionFillOpacity) b.WriteString(``) writeEclipseMapGeoLine(b, frame, startBoundary, "c1-horizon", "#a56c16", 1.2, "4 3", "transit-map-clip", eclipseMapSourceTransitHorizonBoundaries) writeEclipseMapGeoLine(b, frame, endBoundary, "c4-horizon", "#197a82", 1.2, "4 3", "transit-map-clip", eclipseMapSourceTransitHorizonBoundaries) } // writeTransitMapIsochrones 画选中接触的站心等时线：每条线是同一时刻看到该接触的地点，只保留该时刻太阳在视地平上的部分。 +// 返回实际画出线的接触，供图例使用。 // 每隔若干条标一次时刻，标签按 LabelTable 避让。 -func writeTransitMapIsochrones(b *strings.Builder, frame svgmap.Frame, transit basic.PlanetTransitResult, options TransitMapSVGOptions) { +func writeTransitMapIsochrones(b *strings.Builder, frame svgmap.Frame, transit basic.PlanetTransitResult, options TransitMapSVGOptions) TransitIsochrones { + var drawn TransitIsochrones table := &svgchart.LabelTable{} format := "15:04" if options.IsochroneStep%time.Minute != 0 { @@ -337,6 +363,9 @@ func writeTransitMapIsochrones(b *strings.Builder, frame svgmap.Frame, transit b continue } lines := basic.PlanetTransitContactIsochrones(transit, contact.kind, options.IsochroneStep) + if len(lines) > 0 { + drawn |= contact.flag + } labelEvery := transitMapIsochroneLabelEvery(len(lines)) for index, line := range lines { segments := make([][]svgmap.GeoPoint, 0, len(line.Segments)) @@ -360,6 +389,7 @@ func writeTransitMapIsochrones(b *strings.Builder, frame svgmap.Frame, transit b contact.class, placement.X, placement.Y, contact.color, placement.Anchor, html.EscapeString(text)) } } + return drawn } // transitMapIsochroneLabelEvery 让每个接触最多标八条线。 @@ -413,7 +443,10 @@ func writeTransitMapSummary(b *strings.Builder, transit basic.PlanetTransitResul if options.Language == "en" { greatest, separation = "Greatest", "least separation" } - parts := []string{"C1 " + clock(contacts.start)} + parts := []string{} + if !contacts.start.IsZero() { + parts = append(parts, "C1 "+clock(contacts.start)) + } if contacts.hasInternal { parts = append(parts, "C2 "+clock(contacts.internalStart)) } @@ -421,7 +454,13 @@ func writeTransitMapSummary(b *strings.Builder, transit basic.PlanetTransitResul if contacts.hasInternal { parts = append(parts, "C3 "+clock(contacts.internalEnd)) } - parts = append(parts, "C4 "+clock(contacts.end)) + if !contacts.end.IsZero() { + parts = append(parts, "C4 "+clock(contacts.end)) + } else if options.Language == "en" { + parts = append(parts, "Partial transit; no geocentric contacts") + } else { + parts = append(parts, "仅部分凌日，地心无接触") + } line := fmt.Sprintf("%s (%s)", strings.Join(parts, " | "), zone) fmt.Fprintf(b, `%s`, float64(options.Width)/2, html.EscapeString(svgchart.EllipsizeText(line, float64(options.Width)-60, 13))) @@ -432,23 +471,33 @@ func writeTransitMapSummary(b *strings.Builder, transit basic.PlanetTransitResul transit.MinimumSeparationArcsec, transit.SunSemidiameterArcsec, transit.PlanetSemidiameterArcsec) } fmt.Fprintf(b, `%s`, - float64(options.Width)/2, html.EscapeString(note)) + float64(options.Width)/2, html.EscapeString(svgchart.EllipsizeText(transitMapGeometryNote(transit, options.Language)+" | "+note, float64(options.Width)-60, 13))) } -// writeTransitMapLegend 画区域色块与等时线图例，按画布宽度折行，返回行数。 -func writeTransitMapLegend(b *strings.Builder, frame svgmap.Frame, canvasWidth float64, options TransitMapSVGOptions) int { - type entry struct { - text, color string - line bool - } +type transitMapLegendEntry struct { + text, color string + line bool +} + +// transitMapLegendEntries 列出图例项：地心无接触时是部分凌日可见区，否则是四个分区；两者都追加实际画出的等时线。 +func transitMapLegendEntries(options TransitMapSVGOptions, hasExternal bool) []transitMapLegendEntry { labels := []string{"全程可见", "日出时已在凌日", "日落时仍在凌日", "不可见", "一触等时线", "二触等时线", "三触等时线", "四触等时线"} if options.Language == "en" { labels = []string{"Entire transit", "Transit at sunrise", "Transit at sunset", "Not visible", "C1 isochrones", "C2 isochrones", "C3 isochrones", "C4 isochrones"} } - entries := []entry{ - {labels[0], transitMapColorEntire, false}, {labels[1], transitMapColorSunrise, false}, - {labels[2], transitMapColorSunset, false}, {labels[3], transitMapColorNone, false}, + var entries []transitMapLegendEntry + if !hasExternal { + partial := "部分凌日可见区" + if options.Language == "en" { + partial = "Partial transit visible" + } + entries = []transitMapLegendEntry{{partial, transitMapColorPartial, false}} + } else { + entries = []transitMapLegendEntry{ + {labels[0], transitMapColorEntire, false}, {labels[1], transitMapColorSunrise, false}, + {labels[2], transitMapColorSunset, false}, {labels[3], transitMapColorNone, false}, + } } for index, item := range []struct { flag TransitIsochrones @@ -458,14 +507,20 @@ func writeTransitMapLegend(b *strings.Builder, frame svgmap.Frame, canvasWidth f {TransitIsochroneC3, transitMapColorC3}, {TransitIsochroneC4, transitMapColorC4}, } { if options.Isochrones&item.flag != 0 { - entries = append(entries, entry{labels[4+index], item.color, true}) + entries = append(entries, transitMapLegendEntry{labels[4+index], item.color, true}) } } - width := func(item entry) float64 { + return entries +} + +// writeTransitMapLegend 画区域色块与等时线图例，按画布宽度折行，返回行数。 +func writeTransitMapLegend(b *strings.Builder, frame svgmap.Frame, canvasWidth float64, options TransitMapSVGOptions, hasExternal bool) int { + entries := transitMapLegendEntries(options, hasExternal) + var rows [][]transitMapLegendEntry + width := func(item transitMapLegendEntry) float64 { return lunarEclipseDetailedMapLegendIconWidth + svgchart.EstimatedTextWidth(item.text, lunarEclipseDetailedMapLegendFontSize) } available := canvasWidth - 2*lunarEclipseDetailedPageInset - var rows [][]entry used := 0.0 for _, item := range entries { if len(rows) == 0 || used+lunarEclipseDetailedMapLegendGap+width(item) > available { @@ -487,7 +542,7 @@ func writeTransitMapLegend(b *strings.Builder, frame svgmap.Frame, canvasWidth f if item.line { fmt.Fprintf(b, ``, x, y-4, x+18, y-4, item.color) } else { - fmt.Fprintf(b, ``, x, y-8, item.color) + fmt.Fprintf(b, ``, x, y-8, item.color, transitMapLegendFillOpacity) } fmt.Fprintf(b, `%s`, x+lunarEclipseDetailedMapLegendIconWidth, y, lunarEclipseDetailedMapLegendFontSize, html.EscapeString(item.text)) @@ -504,17 +559,25 @@ func writeTransitMapFooter( projection svgmap.Projection, instant time.Time, legendRows int, + hasExternal bool, ) { text := options.FooterNote if text == "" { + projectionLabel := transitMapProjectionLabel(projection, hasExternal, options.Language) if options.Language == "en" { - text = eclipseMapProjectionLabel(projection, "en") + "; regions follow whether the Sun's centre is above the apparent horizon (standard refraction) at C1, greatest and C4" + text = projectionLabel + "; regions follow whether the Sun's centre is above the apparent horizon (standard refraction) at C1, greatest and C4" + if !hasExternal { + text = projectionLabel + "; shaded footprint uses topocentric disk overlap and apparent-horizon visibility (standard refraction)" + } if options.Isochrones != 0 { text += "; isochrones are topocentric contact times" } text += "." } else { - text = eclipseMapProjectionLabel(projection, "zh") + "；按一触、凌甚、四触太阳中心是否在视地平上分区（计入标准大气折射）" + text = projectionLabel + "；按一触、凌甚、四触太阳中心是否在视地平上分区（计入标准大气折射）" + if !hasExternal { + text = projectionLabel + "；着色区域按站心圆面相交与视地平可见性计算（计入标准大气折射）" + } if options.Isochrones != 0 { text += "；等时线为站心接触时刻" } @@ -525,12 +588,18 @@ func writeTransitMapFooter( // 页脚紧跟图例末行，不与图例重叠；行数按画布底边截断，时标声明独占最后一行。 baseline := frame.Y + frame.Height + 31 + float64(legendRows-1)*lunarEclipseDetailedMapLegendLineStep + 18 + x, anchor := frame.X, "start" maxWidth := float64(options.Width) - frame.X - 24 + // 球面投影的盘面居中而左缘靠里，页脚改为随图例居中，否则文字整体右偏。 + if frame.IsDisk() { + x, anchor = float64(options.Width)/2, "middle" + maxWidth = float64(options.Width) - 2*transitMapFooterSideInset + } lines := svgchart.WrapText(text, maxWidth, 11) lines = svgchart.FooterLinesWithScale(lines, timenote.Scale(options.TimeScale, instant, options.Location, options.Language), maxWidth, 11, svgchart.BaselineLineLimit(11, svgchart.FooterLineHeight(11), baseline, float64(options.Height)-footerBottomPadding)) for index, line := range lines { - fmt.Fprintf(b, `%s`, - frame.X, baseline+float64(index)*svgchart.FooterLineHeight(11), html.EscapeString(line)) + fmt.Fprintf(b, `%s`, + x, baseline+float64(index)*svgchart.FooterLineHeight(11), anchor, html.EscapeString(line)) } } diff --git a/eclipse/svg/transit_orthographic_test.go b/eclipse/svg/transit_orthographic_test.go new file mode 100644 index 0000000..ebb841f --- /dev/null +++ b/eclipse/svg/transit_orthographic_test.go @@ -0,0 +1,158 @@ +package svg + +import ( + "fmt" + "math" + "regexp" + "strings" + "testing" + "time" + + "b612.me/astro/basic" +) + +func TestTransitMapVisibilityCenter(tester *testing.T) { + for _, test := range []struct { + name string + columns []basic.PlanetTransitVisibilityColumn + longitude float64 + latitude float64 + valid bool + }{ + {name: "empty"}, + {name: "no-visible-interval", columns: []basic.PlanetTransitVisibilityColumn{{Longitude: -180}, {Longitude: 180}}}, + {name: "antimeridian", columns: []basic.PlanetTransitVisibilityColumn{ + {Longitude: -180, Intervals: [][2]float64{{-10, 10}}}, + {Longitude: -178, Intervals: [][2]float64{{-10, 10}}}, + {Longitude: -176}, {Longitude: 176}, + {Longitude: 178, Intervals: [][2]float64{{-10, 10}}}, + {Longitude: 180, Intervals: [][2]float64{{-10, 10}}}, + }, longitude: 180, valid: true}, + {name: "north", columns: []basic.PlanetTransitVisibilityColumn{ + {Longitude: -180, Intervals: [][2]float64{{70, 90}}}, + {Longitude: -90, Intervals: [][2]float64{{70, 90}}}, + {Longitude: 0, Intervals: [][2]float64{{70, 90}}}, + {Longitude: 90, Intervals: [][2]float64{{70, 90}}}, + {Longitude: 180, Intervals: [][2]float64{{70, 90}}}, + }, latitude: 90, valid: true}, + {name: "south", columns: []basic.PlanetTransitVisibilityColumn{ + {Longitude: -180, Intervals: [][2]float64{{-90, -70}}}, + {Longitude: -90, Intervals: [][2]float64{{-90, -70}}}, + {Longitude: 0, Intervals: [][2]float64{{-90, -70}}}, + {Longitude: 90, Intervals: [][2]float64{{-90, -70}}}, + {Longitude: 180, Intervals: [][2]float64{{-90, -70}}}, + }, latitude: -90, valid: true}, + } { + tester.Run(test.name, func(tester *testing.T) { + center, valid := transitMapVisibilityCenter(test.columns) + if valid != test.valid { + tester.Fatalf("valid=%v, want %v", valid, test.valid) + } + if valid && (math.Abs(center.Latitude-test.latitude) > 1e-8 || + math.Abs(test.latitude) < 90 && math.Abs(math.Remainder(center.Longitude-test.longitude, 360)) > 1e-8) { + tester.Fatalf("center=%+v, want longitude=%v latitude=%v", center, test.longitude, test.latitude) + } + }) + } +} + +func TestTransitMapOrthographicPartialCenter(tester *testing.T) { + for _, test := range []struct { + planet TransitPlanet + year int + month time.Month + day int + longitude float64 + latitude float64 + }{ + {TransitVenus, 2854, time.December, 14, -135, -60}, + {TransitMercury, 1937, time.May, 11, 0, -60}, + } { + tester.Run(fmt.Sprintf("%d-%d", test.planet, test.year), func(tester *testing.T) { + date := time.Date(test.year, test.month, test.day, 0, 0, 0, 0, time.UTC) + transit, valid := closestPlanetTransit(test.planet, date) + if !valid || transit.HasExternal { + tester.Fatal("expected a transit without geocentric contacts") + } + center := transitMapOrthographicCenter(transit) + local := basic.PlanetTransitLocal(transit, basic.Observer{Longitude: center.Longitude, Latitude: center.Latitude}) + if !local.Valid || !local.Visible { + tester.Fatalf("view centre must lie in the visible footprint: %+v", center) + } + for _, projection := range []EclipseMapProjection{EclipseMapProjectionEquirectangular, EclipseMapProjectionSouthPolar, EclipseMapProjectionOrthographic} { + chart, valid := TransitMapSVG(test.planet, date, TransitMapSVGOptions{Projection: projection, Language: "en"}) + if !valid || strings.Contains(chart, "NaN") || strings.Contains(chart, "Inf") { + tester.Fatalf("%s: invalid chart", projection) + } + paths := regexp.MustCompile(`]+>`).FindAllStringSubmatch(chart, -1) + if len(paths) != 1 || !strings.Contains(paths[0][0], `fill-rule="nonzero"`) { + tester.Fatalf("%s: footprint must use a single nonzero fill, got %d paths", projection, len(paths)) + } + frame := eclipseMapFrame(960, 640, internalEclipseMapProjection(projection), 142, 122) + if projection == EclipseMapProjectionOrthographic { + frame.CenterLongitude, frame.CenterLatitude = center.Longitude, center.Latitude + if !strings.Contains(chart, fmt.Sprintf(`data-center-longitude="%.6f" data-center-latitude="%.6f"`, center.Longitude, center.Latitude)) || + !strings.Contains(chart, "centred on the topocentric visible footprint") { + tester.Fatal("orthographic view centre or explanation missing") + } + } + x, y, projectable := frame.Project(test.longitude, test.latitude) + if !projectable || !lunarPenumbralTestInside(lunarPenumbralTestRings(tester, paths[0][1]), x, y) { + tester.Fatalf("%s: known visible site is outside the painted footprint", projection) + } + } + }) + } +} + +func TestTransitMapOrthographicCenterFallback(tester *testing.T) { + for _, test := range []struct { + planet TransitPlanet + date time.Time + }{ + {TransitMercury, time.Date(2019, 11, 11, 0, 0, 0, 0, time.UTC)}, + {TransitVenus, time.Date(2012, 6, 6, 0, 0, 0, 0, time.UTC)}, + } { + transit, valid := closestPlanetTransit(test.planet, test.date) + if !valid || !transit.HasExternal { + tester.Fatal("expected an ordinary geocentric transit") + } + want := solarEclipseSubsolarPoint(transitMapContactsOf(transit).greatest) + if got := transitMapOrthographicCenter(transit); got != want { + tester.Fatalf("ordinary transit view changed: got %+v want %+v", got, want) + } + transit.Valid, transit.HasExternal = false, false + if got := transitMapOrthographicCenter(transit); got != want { + tester.Fatalf("empty visibility must fall back to the subsolar point: %+v", got) + } + } +} + +func TestTransitMapFooterAlignment(tester *testing.T) { + footer := regexp.MustCompile(``) + date := time.Date(2854, 12, 14, 0, 0, 0, 0, time.UTC) + for _, test := range []struct { + projection EclipseMapProjection + x, anchor string + }{ + {EclipseMapProjectionOrthographic, "480.000", "middle"}, + {"", "104.000", "start"}, + } { + chart, ok := TransitMapSVG(TransitVenus, date, TransitMapSVGOptions{Projection: test.projection, Language: "en"}) + if !ok { + tester.Fatalf("%q: render failed", test.projection) + } + lines := footer.FindAllStringSubmatch(chart, -1) + if len(lines) < 2 { + tester.Fatalf("%q: footer lines=%d", test.projection, len(lines)) + } + for _, line := range lines { + if line[1] != test.x || line[2] != test.anchor { + tester.Fatalf("%q: footer at x=%s anchor=%s, want %s %s", test.projection, line[1], line[2], test.x, test.anchor) + } + } + if !strings.Contains(chart, `class="transit-partial-region"`) || strings.Contains(chart, "#747b7d") || !strings.Contains(chart, `fill="#f7dd93"`) { + tester.Fatalf("%q: partial region fill changed", test.projection) + } + } +} diff --git a/eclipse/svg/transit_partial.go b/eclipse/svg/transit_partial.go new file mode 100644 index 0000000..14aec10 --- /dev/null +++ b/eclipse/svg/transit_partial.go @@ -0,0 +1,94 @@ +package svg + +import ( + "fmt" + "math" + "strings" + + "b612.me/astro/basic" + "b612.me/astro/internal/svgmap" +) + +func transitMapOrthographicCenter(transit basic.PlanetTransitResult) svgmap.GeoPoint { + if !transit.HasExternal { + if center, ok := transitMapVisibilityCenter(basic.PlanetTransitVisibilityColumns(transit)); ok { + return center + } + } + return solarEclipseSubsolarPoint(transitMapContactsOf(transit).greatest) +} + +func transitMapVisibilityCenter(columns []basic.PlanetTransitVisibilityColumn) (svgmap.GeoPoint, bool) { + var momentX, momentY, momentZ float64 + for index, column := range columns { + weight := 0.0 + if index > 0 { + weight += (column.Longitude - columns[index-1].Longitude) / 2 + } + if index+1 < len(columns) { + weight += (columns[index+1].Longitude - column.Longitude) / 2 + } + longitude := column.Longitude * math.Pi / 180 + for _, interval := range column.Intervals { + lower, upper := interval[0]*math.Pi/180, interval[1]*math.Pi/180 + horizontalMoment := (upper-lower)/2 + (math.Sin(2*upper)-math.Sin(2*lower))/4 + verticalMoment := (math.Sin(upper)*math.Sin(upper) - math.Sin(lower)*math.Sin(lower)) / 2 + momentX += weight * horizontalMoment * math.Cos(longitude) + momentY += weight * horizontalMoment * math.Sin(longitude) + momentZ += weight * verticalMoment + } + } + if math.Sqrt(momentX*momentX+momentY*momentY+momentZ*momentZ) < 1e-12 { + return svgmap.GeoPoint{}, false + } + return svgmap.GeoPoint{ + Longitude: math.Atan2(momentY, momentX) * 180 / math.Pi, + Latitude: math.Atan2(momentZ, math.Hypot(momentX, momentY)) * 180 / math.Pi, + }, true +} + +func transitMapProjectionLabel(projection svgmap.Projection, hasExternal bool, language string) string { + if projection != svgmap.ProjectionOrthographic { + return eclipseMapProjectionLabel(projection, language) + } + if language == "en" { + center := "the subsolar point at greatest" + if !hasExternal { + center = "the topocentric visible footprint" + } + return "Orthographic globe projection, centred on " + center + "; one hemisphere only" + } + center := "凌甚太阳直射点" + if !hasExternal { + center = "站心可见区中心" + } + return "正射球面投影，视点取" + center + "，只画朝向视点的半个地球" +} + +func writeTransitMapPartialRegions(builder *strings.Builder, frame svgmap.Frame, transit basic.PlanetTransitResult) { + fmt.Fprintf(builder, ``, eclipseMapSourceTransitPartialVisibility) + var path strings.Builder + for _, polygon := range basic.PlanetTransitVisibilityPolygons(transit) { + points := make([]svgmap.GeoPoint, len(polygon)) + for index, point := range polygon { + points[index] = svgmap.GeoPoint{Longitude: point.Longitude, Latitude: point.Latitude} + } + path.WriteString(lunarEclipseGeoPolygonPath(points, frame)) + } + if path.Len() > 0 { + fmt.Fprintf(builder, ``, path.String(), transitMapColorPartial, transitMapPartialFillOpacity) + } + builder.WriteString(``) +} + +func transitMapGeometryNote(transit basic.PlanetTransitResult, language string) string { + text := fmt.Sprintf("γ %.4f", transit.Gamma) + if transit.Series.Number != 0 { + if language == "en" { + text += fmt.Sprintf(" | Series %d (%d yr)", transit.Series.Number, transit.Series.PeriodYears) + } else { + text += fmt.Sprintf(" | 周期族 %d（%d 年）", transit.Series.Number, transit.Series.PeriodYears) + } + } + return text +} diff --git a/eclipse/svg/transit_partial_test.go b/eclipse/svg/transit_partial_test.go new file mode 100644 index 0000000..b2a4e6d --- /dev/null +++ b/eclipse/svg/transit_partial_test.go @@ -0,0 +1,106 @@ +package svg + +import ( + "regexp" + "strconv" + "strings" + "testing" + "time" +) + +func TestTransitSVGGlobalPartial(tester *testing.T) { + date := time.Date(2854, 12, 14, 0, 0, 0, 0, time.UTC) + chart, ok := TransitMapSVG(TransitVenus, date, TransitMapSVGOptions{Isochrones: TransitIsochroneAll, IsochroneStep: 30 * time.Second}) + if !ok { + tester.Fatal("global partial transit render failed") + } + for _, expected := range []string{"2854-12-14", "transit-partial-region", "地心无接触", "周期族 4", "γ ", "c1-isochrone", "c4-isochrone"} { + if !strings.Contains(chart, expected) { + tester.Fatalf("missing %s", expected) + } + } + for _, absent := range []string{"NaN", "Inf", "transit-sun-c1-shape", "c2-isochrone", "c3-isochrone"} { + if strings.Contains(chart, absent) { + tester.Fatalf("invalid SVG contains %s", absent) + } + } + local, ok := TransitLocalSVG(TransitVenus, date, -135, -60, 0, TransitLocalSVGOptions{}) + if !ok || !strings.Contains(local, "transit-planet-c1") || !strings.Contains(local, "transit-planet-c4") || strings.Contains(local, "transit-planet-c2") { + tester.Fatal("local partial-transit diagram has invalid contacts") + } + labels := regexp.MustCompile(`class="transit-event-label" x="[0-9.-]+" y="([0-9.-]+)"`).FindAllStringSubmatch(local, -1) + if len(labels) != 3 { + tester.Fatalf("grazing contact labels=%d, want 3", len(labels)) + } + footer := regexp.MustCompile(` firstFooterBaseline-8 { + tester.Fatal("contact label overlaps the footer") + } + } + if _, ok := TransitLocalSVG(TransitVenus, date, 0, 45, 0, TransitLocalSVGOptions{}); ok { + tester.Fatal("a northern site must miss the transit") + } +} + +func TestTransitMapLegendListsDrawnIsochrones(tester *testing.T) { + for _, test := range []struct { + planet TransitPlanet + date time.Time + inner bool + }{ + {TransitVenus, time.Date(554, 5, 22, 0, 0, 0, 0, time.UTC), true}, + {TransitMercury, time.Date(-513, 10, 22, 0, 0, 0, 0, time.UTC), false}, + } { + chart, ok := TransitMapSVG(test.planet, test.date, TransitMapSVGOptions{Isochrones: TransitIsochroneAll, Language: "en"}) + if !ok { + tester.Fatalf("%v: render failed", test.date) + } + for _, marker := range []string{"c2-isochrone", "C2 isochrones", "C3 isochrones"} { + if strings.Contains(chart, marker) != test.inner { + tester.Fatalf("%v: %s present=%v, want %v", test.date, marker, !test.inner, test.inner) + } + } + } +} + +// 地心无接触的部分凌日也要在图例里列出实际画出的一触、四触等时线；二触、三触不存在。 +func TestTransitPartialLegendListsDrawnIsochrones(tester *testing.T) { + for _, language := range []string{"zh", "en"} { + chart, ok := TransitMapSVG(TransitVenus, time.Date(2854, 12, 14, 0, 0, 0, 0, time.UTC), TransitMapSVGOptions{ + Isochrones: TransitIsochroneAll, Language: language, + }) + if !ok { + tester.Fatalf("%s: render failed", language) + } + partial, c1, c4 := "部分凌日可见区", "一触等时线", "四触等时线" + if language == "en" { + partial, c1, c4 = "Partial transit visible", "C1 isochrones", "C4 isochrones" + } + for _, expected := range []string{"c1-isochrone", "c4-isochrone", partial, c1, c4} { + if !strings.Contains(chart, expected) { + tester.Fatalf("%s: partial chart is missing %q", language, expected) + } + } + for _, absent := range []string{"c2-isochrone", "c3-isochrone", "二触等时线", "三触等时线", "C2 isochrones", "C3 isochrones"} { + if strings.Contains(chart, absent) { + tester.Fatalf("%s: partial chart must not contain %q", language, absent) + } + } + } + // 默认不画等时线时图例只留可见区色块，不凭空加等时线条目。 + chart, ok := TransitMapSVG(TransitVenus, time.Date(2854, 12, 14, 0, 0, 0, 0, time.UTC), TransitMapSVGOptions{Language: "en"}) + if !ok { + tester.Fatal("default render failed") + } + for _, absent := range []string{"c1-isochrone", "C1 isochrones", "C4 isochrones"} { + if strings.Contains(chart, absent) { + tester.Fatalf("default chart must not contain %q", absent) + } + } +} diff --git a/eclipse/svg/transit_text_overlap_test.go b/eclipse/svg/transit_text_overlap_test.go new file mode 100644 index 0000000..0317b6f --- /dev/null +++ b/eclipse/svg/transit_text_overlap_test.go @@ -0,0 +1,126 @@ +package svg + +import ( + "fmt" + "regexp" + "strings" + "testing" + "time" +) + +var transitLegendSwatchPattern = regexp.MustCompile(`]*width="18" height="9"[^>]*fill-opacity="([0-9.]+)"`) + +type transitOverlapCase struct { + name string + render func(width, height int) (string, bool) +} + +func transitOverlapCases() []transitOverlapCase { + partial := time.Date(1937, time.May, 11, 0, 0, 0, 0, time.UTC) + ordinary := time.Date(2012, time.June, 6, 0, 0, 0, 0, time.UTC) + return []transitOverlapCase{ + {"local-partial-zh", func(width, height int) (string, bool) { + return TransitLocalSVG(TransitMercury, partial, 0, -60, 0, TransitLocalSVGOptions{Width: width, Height: height}) + }}, + {"local-partial-en", func(width, height int) (string, bool) { + return TransitLocalSVG(TransitMercury, partial, 0, -60, 0, TransitLocalSVGOptions{Width: width, Height: height, Language: "en"}) + }}, + {"local-full-zh", func(width, height int) (string, bool) { + return TransitLocalSVG(TransitVenus, ordinary, 116.4, 39.9, 0, TransitLocalSVGOptions{Width: width, Height: height}) + }}, + {"local-full-en", func(width, height int) (string, bool) { + return TransitLocalSVG(TransitVenus, ordinary, 116.4, 39.9, 0, TransitLocalSVGOptions{Width: width, Height: height, Language: "en"}) + }}, + {"map-partial", func(width, height int) (string, bool) { + return TransitMapSVG(TransitMercury, partial, TransitMapSVGOptions{Width: width, Height: height, Isochrones: TransitIsochroneAll}) + }}, + {"map-partial-orthographic", func(width, height int) (string, bool) { + return TransitMapSVG(TransitMercury, partial, TransitMapSVGOptions{ + Width: width, Height: height, Isochrones: TransitIsochroneAll, Projection: EclipseMapProjectionOrthographic, + }) + }}, + {"map-ordinary", func(width, height int) (string, bool) { + return TransitMapSVG(TransitVenus, ordinary, TransitMapSVGOptions{Width: width, Height: height, Isochrones: TransitIsochroneAll}) + }}, + } +} + +// 凌日图与日月食图同一口径：任意两个尺寸下返回的图都必须 0 文字压叠。 +func TestTransitChartTextOverlapMatrix(t *testing.T) { + for _, size := range [][2]int{{640, 420}, {640, 480}, {800, 600}, {960, 640}, {960, 680}, {1200, 800}} { + for _, test := range transitOverlapCases() { + name := fmt.Sprintf("%s/%dx%d", test.name, size[0], size[1]) + document, ok := test.render(size[0], size[1]) + if !ok { + continue + } + assertNoTextOverlap(t, name, document) + } + } + // 文档下限尺寸必须继续出图，拒绝只允许出现在真放不下的画布上。 + if _, ok := TransitLocalSVG(TransitMercury, time.Date(1937, time.May, 11, 0, 0, 0, 0, time.UTC), 0, -60, 0, + TransitLocalSVGOptions{Width: 640, Height: 480}); !ok { + t.Error("documented 640x480 local transit canvas was rejected") + } +} + +// 地心无接触的可见区填充必须落在 data-source 词表内，消费方才能判断该层来源。 +func TestTransitPartialRegionsUseDocumentedDataSource(t *testing.T) { + chart, ok := TransitMapSVG(TransitMercury, time.Date(1937, time.May, 11, 0, 0, 0, 0, time.UTC), TransitMapSVGOptions{}) + if !ok { + t.Fatal("partial transit render failed") + } + if !strings.Contains(chart, `data-source="`+eclipseMapSourceTransitPartialVisibility+`"`) { + t.Fatalf("partial regions carry no documented data-source marker") + } + for _, match := range eclipseMapDataSourcePattern.FindAllStringSubmatch(chart, -1) { + known := false + for _, source := range eclipseMapDataSources { + known = known || source == match[1] + } + if !known { + t.Fatalf("undocumented data-source %q", match[1]) + } + } +} + +// 图例色块画在白底上取月食图同值；部分凌日单块可见区比普通分区更浓。 +func TestTransitMapFillOpacities(t *testing.T) { + for _, test := range []struct { + planet TransitPlanet + date time.Time + pattern *regexp.Regexp + region float64 + }{ + {TransitVenus, time.Date(2012, time.June, 6, 0, 0, 0, 0, time.UTC), + regexp.MustCompile(`class="transit-[a-z-]+-region"[^>]*fill-opacity="([0-9.]+)"`), transitMapRegionFillOpacity}, + {TransitMercury, time.Date(1937, time.May, 11, 0, 0, 0, 0, time.UTC), + regexp.MustCompile(`class="transit-partial-region"[^>]*fill-opacity="([0-9.]+)"`), transitMapPartialFillOpacity}, + } { + chart, ok := TransitMapSVG(test.planet, test.date, TransitMapSVGOptions{}) + if !ok { + t.Fatalf("%v: render failed", test.date) + } + regions := test.pattern.FindAllStringSubmatch(chart, -1) + if len(regions) == 0 { + t.Fatalf("%v: no region fill found", test.date) + } + for _, match := range regions { + if match[1] != fmt.Sprintf("%.2f", test.region) { + t.Fatalf("%v: region fill-opacity %s, want %.2f", test.date, match[1], test.region) + } + } + legends := transitLegendSwatchPattern.FindAllStringSubmatch(chart, -1) + if len(legends) == 0 { + t.Fatalf("%v: no legend swatch found", test.date) + } + for _, match := range legends { + if match[1] != fmt.Sprintf("%.2f", transitMapLegendFillOpacity) { + t.Fatalf("%v: legend fill-opacity %s, want %.2f", test.date, match[1], transitMapLegendFillOpacity) + } + } + } + if transitMapPartialFillOpacity <= transitMapRegionFillOpacity { + t.Fatal("partial-transit footprint must be denser than ordinary regions") + } +} diff --git a/geojson/transit.go b/geojson/transit.go index b90bff5..fa4fc04 100644 --- a/geojson/transit.go +++ b/geojson/transit.go @@ -28,16 +28,17 @@ type PlanetTransitOptions struct { // 小于等于零时使用 360，其他值限制在 [12, 1440]。 // BoundaryPoints is the horizon sample count and also sets the region longitude columns to // max(360, BoundaryPoints), capped at 720. Values <= 0 use 360; others are clamped to [12, 1440]. + // 地心无接触事件使用固定站心网格，不由 BoundaryPoints 控制 / events without geocentric contacts use a fixed topocentric grid. BoundaryPoints int - // IsochroneStep 为正时输出一触到四触的站心接触等时线（掠凌没有二触、三触），至少一秒；零值不输出。 - // A positive IsochroneStep adds topocentric C1-C4 isochrones (no C2/C3 for a grazing transit) and must be + // IsochroneStep 为正时输出一触到四触的站心接触等时线（地表没有站心内切的事件没有二触、三触），至少一秒；零值不输出。 + // A positive IsochroneStep adds topocentric C1-C4 isochrones (no C2/C3 where no site has internal contact) and must be // at least one second. The zero value adds none. IsochroneStep time.Duration - // TimeMarkers 非空时沿一触到四触的太阳直射点轨迹追加时间标记 Point 要素。 - // TimeMarkers adds time-marker Point Features along the subsolar track from C1 through C4 when non-nil. + // TimeMarkers 非空时沿太阳直射点轨迹追加时间标记；无地心接触时覆盖求解窗口。 + // TimeMarkers adds subsolar-track markers over C1-C4, or over the search window when geocentric contacts are absent. TimeMarkers *TimeMarkerOptions - // SkipRoles 列出不写进输出的 role；分区三个 role 都跳过时不做分区计算。 - // SkipRoles lists roles to leave out; skipping all three region roles skips the region computation. + // SkipRoles 列出不写进输出的 role；分区 role 全部跳过时不做分区计算。 + // SkipRoles lists roles to leave out; skipping every region role skips the region computation. SkipRoles []string } @@ -54,6 +55,12 @@ type planetTransitInput struct { separation float64 sunRadius float64 planetRad float64 + hasExternal bool + windowStart time.Time + windowEnd time.Time + gamma float64 + penumbra float64 + series basic.PlanetTransitSeries } // MarshalMercuryTransit 将水星凌日的全球可见分区、一触/四触地平线与可选站心等时线编码为 GeoJSON / encode a Mercury transit as GeoJSON. @@ -66,6 +73,8 @@ func (module *Generator) MarshalMercuryTransit(info mercury.TransitInfo, options start: info.Start, c2: info.InternalStart, greatest: info.Greatest, c3: info.InternalEnd, end: info.End, hasInternal: info.HasInternal, separation: info.MinimumSeparationArcsec, sunRadius: info.SunSemidiameterArcsec, planetRad: info.PlanetSemidiameterArcsec, + hasExternal: !info.Start.IsZero() && !info.End.IsZero(), windowStart: info.WindowStart, windowEnd: info.WindowEnd, + gamma: info.Gamma, penumbra: info.PenumbralRadius, series: info.Series, }, options, module) } @@ -79,10 +88,32 @@ func (module *Generator) MarshalVenusTransit(info venus.TransitInfo, options Pla start: info.Start, c2: info.InternalStart, greatest: info.Greatest, c3: info.InternalEnd, end: info.End, hasInternal: info.HasInternal, separation: info.MinimumSeparationArcsec, sunRadius: info.SunSemidiameterArcsec, planetRad: info.PlanetSemidiameterArcsec, + hasExternal: !info.Start.IsZero() && !info.End.IsZero(), windowStart: info.WindowStart, windowEnd: info.WindowEnd, + gamma: info.Gamma, penumbra: info.PenumbralRadius, series: info.Series, }, options, module) } +// planetTransitHasOnlyWindow 判断输入是否只给了站心求解窗口：地心无接触，接触时刻全为零。 +func planetTransitHasOnlyWindow(input planetTransitInput) bool { + if input.hasExternal || input.hasInternal { + return false + } + if !input.start.IsZero() || !input.end.IsZero() || !input.c2.IsZero() || !input.c3.IsZero() { + return false + } + if input.greatest.IsZero() || input.windowStart.IsZero() || input.windowEnd.IsZero() { + return false + } + if !input.windowStart.Before(input.greatest) || !input.greatest.Before(input.windowEnd) { + return false + } + return input.windowEnd.Sub(input.windowStart) <= 24*time.Hour +} + func validatePlanetTransitInput(input planetTransitInput) error { + if planetTransitHasOnlyWindow(input) { + return nil + } if input.start.IsZero() || input.greatest.IsZero() || input.end.IsZero() { return fmt.Errorf("geojson: %s transit contact times are required", input.planet) } @@ -140,7 +171,9 @@ func marshalPlanetTransit(input planetTransitInput, options PlanetTransitOptions longitudes[index] = -180 + 360*float64(index)/float64(longitudePoints) } refraction := basic.RefractionFromApparentAltitude(0, planetTransitRefractionPressureHPa, planetTransitRefractionTemperatureC) - base := map[string]interface{}{"planet": input.planet, "has_internal": input.hasInternal} + base := map[string]interface{}{"planet": input.planet, "has_internal": input.hasInternal, + "has_geocentric_contacts": input.hasExternal, "gamma": input.gamma, "penumbral_radius": input.penumbra, + "series": input.series.Number, "series_period_years": input.series.PeriodYears} withBase := func(extra map[string]interface{}) map[string]interface{} { properties := cloneProperties(base) for key, value := range extra { @@ -148,6 +181,9 @@ func marshalPlanetTransit(input planetTransitInput, options PlanetTransitOptions } return properties } + if !input.hasExternal { + return marshalPlanetTransitPartial(input, options, timeScale, withBase, module) + } features := make([]feature, 0, 12) columnsAt := func(at time.Time) [][]lunarLatitudeInterval { @@ -354,6 +390,13 @@ func planetTransitRegionFeatures( return features, nil } +func planetTransitResultOf(input planetTransitInput) basic.PlanetTransitResult { + if !input.hasExternal { + return basic.PlanetTransitResultFromWindow(input.planetIndex, input.windowStart, input.greatest, input.windowEnd) + } + return basic.PlanetTransitResultFromTimes(input.planetIndex, input.start, input.c2, input.greatest, input.c3, input.end) +} + // planetTransitIsochroneFeatures 把每个接触的站心等时线编码为一条带 times 的 MultiLineString： // 每段是一条支路，times 与坐标对齐，同一支路上所有点同刻。 func planetTransitIsochroneFeatures( @@ -364,7 +407,7 @@ func planetTransitIsochroneFeatures( withBase func(map[string]interface{}) map[string]interface{}, module *Generator, ) ([]feature, error) { - result := basic.PlanetTransitResultFromTimes(input.planetIndex, input.start, input.c2, input.greatest, input.c3, input.end) + result := planetTransitResultOf(input) features := make([]feature, 0, 4) for _, contact := range []struct { kind basic.PlanetTransitContactKind diff --git a/geojson/transit_partial.go b/geojson/transit_partial.go new file mode 100644 index 0000000..aefaf6c --- /dev/null +++ b/geojson/transit_partial.go @@ -0,0 +1,58 @@ +package geojson + +import ( + "fmt" + "time" + + "b612.me/astro/basic" +) + +func marshalPlanetTransitPartial(input planetTransitInput, options PlanetTransitOptions, scale basic.TimeScale, withBase func(map[string]interface{}) map[string]interface{}, module *Generator) ([]byte, error) { + result := planetTransitResultOf(input) + label := func(value time.Time) string { return formatTime(basic.LabelIn(scale, value.UTC())) } + features := make([]feature, 0, 4) + if !skippedLunarRole(options.SkipRoles, "transit-partial") { + columns := basic.PlanetTransitVisibilityColumnsWithInstance(result, module.inst) + longitudes := make([]float64, len(columns)) + intervals := make([][]lunarLatitudeInterval, len(columns)) + for index, column := range columns { + longitudes[index] = column.Longitude + for _, interval := range column.Intervals { + intervals[index] = append(intervals[index], lunarLatitudeInterval{low: interval[0], high: interval[1]}) + } + } + value, err := planetTransitRegionGeometry(longitudes, intervals) + if err != nil { + return nil, fmt.Errorf("geojson: partial transit: %w", err) + } + features = append(features, newFeature(planetTransitEvent, "transit-partial", value, + withBase(map[string]interface{}{"longitude_step_degrees": 2, "latitude_sample_degrees": 1}))) + } + if options.IsochroneStep > 0 { + lines, err := planetTransitIsochroneFeatures(input, options.IsochroneStep, scale, options.SkipRoles, withBase, module) + if err != nil { + return nil, err + } + features = append(features, lines...) + } + point := solarSubsolarPoint(input.greatest, module) + var err error + features, err = appendPointFeature(features, planetTransitEvent, "greatest", + pathSample{Time: basic.LabelIn(scale, input.greatest.UTC()), Longitude: point.Longitude, Latitude: point.Latitude}, + withBase(map[string]interface{}{"greatest": label(input.greatest), "window_start": label(input.windowStart), "window_end": label(input.windowEnd), + "minimum_separation_arcsec": input.separation, "sun_semidiameter_arcsec": input.sunRadius, "planet_semidiameter_arcsec": input.planetRad})) + if err != nil { + return nil, err + } + if options.TimeMarkers != nil { + markers, markerErr := planetTransitTimeMarkerSamples(input.windowStart, input.windowEnd, scale, *options.TimeMarkers, module) + if markerErr != nil { + return nil, markerErr + } + features, err = appendTimeMarkerPointFeatures(features, planetTransitEvent, "subsolar-track", markers, options.TimeMarkers.Location) + if err != nil { + return nil, err + } + } + return marshalFeatureCollectionWithTimeScale(dropFeaturesByRole(features, options.SkipRoles), scale) +} diff --git a/geojson/transit_partial_test.go b/geojson/transit_partial_test.go new file mode 100644 index 0000000..1751236 --- /dev/null +++ b/geojson/transit_partial_test.go @@ -0,0 +1,44 @@ +package geojson_test + +import ( + "encoding/json" + "testing" + "time" + + "b612.me/astro/geojson" + "b612.me/astro/venus" +) + +func TestMarshalVenusTransitGlobalPartial(tester *testing.T) { + info := venus.NextTransit(time.Date(2854, 1, 1, 0, 0, 0, 0, time.UTC)) + data, err := geojson.MarshalVenusTransit(info, geojson.PlanetTransitOptions{IsochroneStep: 30 * time.Second}) + if err != nil { + tester.Fatal(err) + } + collection, roles := decodeTransit(tester, data) + for _, role := range []string{"transit-partial", "greatest", "c1-isochrone", "c4-isochrone"} { + if roles[role] != 1 { + tester.Fatalf("missing %s: %v", role, roles) + } + } + for _, role := range []string{"c1-horizon", "c4-horizon", "visible-at-c1", "visible-at-c4", "c2-isochrone", "c3-isochrone", "transit-entire"} { + if roles[role] != 0 { + tester.Fatalf("fabricated %s", role) + } + } + for _, item := range collection.Features { + if item.Properties["role"] != "transit-partial" { + continue + } + var polygons [][][][]float64 + if err := json.Unmarshal(item.Geometry.Coordinates, &polygons); err != nil { + tester.Fatal(err) + } + if !geoJSONMultiPolygonContains(polygons, -135, -60) || geoJSONMultiPolygonContains(polygons, 0, 45) { + tester.Fatal("partial-transit footprint does not match the topocentric geometry") + } + if item.Properties["has_geocentric_contacts"] != false || item.Properties["series"] != float64(4) || item.Properties["gamma"] == nil { + tester.Fatalf("missing global-transit metadata: %v", item.Properties) + } + } +} diff --git a/kml/kml.go b/kml/kml.go index 9d7e1de..bb172e1 100644 --- a/kml/kml.go +++ b/kml/kml.go @@ -1254,6 +1254,7 @@ var lookAtPriority = map[string]bool{ "occultation-band": true, "total-footprint": true, "sublunar-track": true, + "transit-partial": true, } func lookAtElements(groups []roleGeometry) [][]element { @@ -1386,7 +1387,7 @@ func paletteStyle(event, role string) Style { return Style{LineColor: ColorOrange, LineWidth: 2} case "c1-horizon", "c4-horizon", "visible-at-c1", "visible-at-c4": return Style{LineColor: ColorOrange, LineWidth: 2} - case "transit-entire": + case "transit-entire", "transit-partial": return Style{LineColor: transitEntireColor, FillColor: transitEntireFill, LineWidth: 1.5} case "transit-at-sunrise": return Style{LineColor: transitSunriseColor, FillColor: transitSunriseFill, LineWidth: 1.5} @@ -1629,6 +1630,7 @@ var roleNames = map[string][2]string{ "visible-at-c1": {"一触可见区", "Visible at C1"}, "visible-at-c4": {"四触可见区", "Visible at C4"}, "transit-entire": {"全程可见", "Entire transit visible"}, + "transit-partial": {"部分凌日可见区", "Partial transit visible"}, "transit-at-sunrise": {"日出时已在凌日", "Transit in progress at sunrise"}, "transit-at-sunset": {"日落时仍在凌日", "Transit in progress at sunset"}, "c1-isochrone": {"一触等时线", "C1 isochrones"}, diff --git a/kml/transit_partial_test.go b/kml/transit_partial_test.go new file mode 100644 index 0000000..2a6c8f1 --- /dev/null +++ b/kml/transit_partial_test.go @@ -0,0 +1,35 @@ +package kml + +import ( + "strings" + "testing" + "time" + + "b612.me/astro/geojson" + "b612.me/astro/venus" +) + +func TestPlanetTransitGlobalPartialConversion(tester *testing.T) { + info := venus.NextTransit(time.Date(2854, 1, 1, 0, 0, 0, 0, time.UTC)) + data, err := geojson.MarshalVenusTransit(info, geojson.PlanetTransitOptions{IsochroneStep: time.Minute}) + if err != nil { + tester.Fatal(err) + } + root := convert(tester, string(data), Options{}) + if !strings.Contains(root.Document.Name, "2854-12-14") { + tester.Fatal(root.Document.Name) + } + if root.Document.LookAt == nil || root.Document.LookAt.Latitude > -50 { + tester.Fatal("KML must focus on the actual partial-transit footprint") + } + roles := map[string]bool{} + for _, item := range root.Document.Folders { + roles[item.Name] = true + } + if !roles["部分凌日可见区"] || !roles["一触等时线"] || !roles["四触等时线"] || roles["二触等时线"] || roles["一触视地平线"] { + tester.Fatalf("invalid partial-transit layers: %v", roles) + } + if paletteStyle(eventPlanetTransit, "transit-partial").FillColor == "" { + tester.Fatal("missing partial-transit fill") + } +} diff --git a/mercury/api_default.go b/mercury/api_default.go index 98cb10c..1c9efb4 100644 --- a/mercury/api_default.go +++ b/mercury/api_default.go @@ -78,7 +78,7 @@ func BrightLimbPositionAngleN(date time.Time, n int) float64 { return defaultModule.BrightLimbPositionAngleN(date, n) } -// ClosestTransit 最近一次地心水星凌日 / closest geocentric Mercury transit. +// ClosestTransit 最近一次全球水星凌日 / closest global Mercury transit. func ClosestTransit(date time.Time) TransitInfo { return defaultModule.ClosestTransit(date) } // CulminationTime 中天时刻，取 date 所在当地民用日 / culmination time, on the local civil day of date. @@ -178,10 +178,10 @@ func LastSuperiorConjunction(date time.Time) time.Time { return defaultModule.LastSuperiorConjunction(date) } -// LastTransit 上一次地心水星凌日 / previous geocentric Mercury transit. +// LastTransit 上一次全球水星凌日 / previous global Mercury transit. func LastTransit(date time.Time) TransitInfo { return defaultModule.LastTransit(date) } -// LocalTransit 地心凌日在观测地的站心接触与凌甚 / topocentric contacts of a transit seen from the site. +// LocalTransit 凌日在观测地的站心接触与凌甚 / topocentric contacts of a transit seen from the site. func LocalTransit(transit TransitInfo, lon, lat, height float64) LocalTransitInfo { return defaultModule.LocalTransit(transit, lon, lat, height) } @@ -227,7 +227,7 @@ func NextSuperiorConjunction(date time.Time) time.Time { return defaultModule.NextSuperiorConjunction(date) } -// NextTransit 下一次地心水星凌日 / next geocentric Mercury transit. +// NextTransit 下一次全球水星凌日 / next global Mercury transit. func NextTransit(date time.Time) TransitInfo { return defaultModule.NextTransit(date) } // ParallacticAngle 水星视差角（天顶方向角） / Mercury parallactic angle. diff --git a/mercury/transit.go b/mercury/transit.go index 259109d..db4917f 100644 --- a/mercury/transit.go +++ b/mercury/transit.go @@ -6,12 +6,12 @@ import ( "b612.me/astro/basic" ) -// TransitInfo 地心水星凌日信息 / geocentric Mercury transit information. +// TransitInfo 全球水星凌日与地心接触 / global Mercury transit and geocentric contacts. // // Start、Greatest、End、InternalStart、InternalEnd 都保持调用者输入的时区。 -// 内切接触不存在时 InternalStart / InternalEnd 为零值。 +// 地心接触不存在时对应时刻为零值；WindowStart、WindowEnd 只是站心求解窗口。 // Start, Greatest, End, InternalStart, and InternalEnd preserve the caller's timezone. -// InternalStart and InternalEnd are zero values when internal contacts do not exist. +// Absent geocentric contacts are zero; WindowStart and WindowEnd are only a local search window. type TransitInfo struct { Valid bool @@ -28,20 +28,26 @@ type TransitInfo struct { SunSemidiameterArcsec float64 PlanetSemidiameterArcsec float64 - HasInternal bool + HasInternal bool + HasExternal bool + WindowStart time.Time + WindowEnd time.Time + Gamma float64 + PenumbralRadius float64 + Series basic.PlanetTransitSeries } -// NextTransit 下一次地心水星凌日 / next geocentric Mercury transit. +// NextTransit 下一次全球水星凌日 / next global Mercury transit. func (module *Mercury) NextTransit(date time.Time) TransitInfo { return transitInfoFromBasic(basic.NextMercuryTransitWithInstance(basic.Date2JD(date.UTC()), module.inst), date.Location()) } -// LastTransit 上一次地心水星凌日 / previous geocentric Mercury transit. +// LastTransit 上一次全球水星凌日 / previous global Mercury transit. func (module *Mercury) LastTransit(date time.Time) TransitInfo { return transitInfoFromBasic(basic.LastMercuryTransitWithInstance(basic.Date2JD(date.UTC()), module.inst), date.Location()) } -// ClosestTransit 最近一次地心水星凌日 / closest geocentric Mercury transit. +// ClosestTransit 最近一次全球水星凌日 / closest global Mercury transit. func (module *Mercury) ClosestTransit(date time.Time) TransitInfo { return transitInfoFromBasic(basic.ClosestMercuryTransitWithInstance(basic.Date2JD(date.UTC()), module.inst), date.Location()) } @@ -50,19 +56,25 @@ func transitInfoFromBasic(result basic.PlanetTransitResult, loc *time.Location) if !result.Valid { return TransitInfo{} } - start := basic.JD2DateByZone(result.ExternalIngress, loc, false) greatest := basic.JD2DateByZone(result.Greatest, loc, false) - end := basic.JD2DateByZone(result.ExternalEgress, loc, false) info := TransitInfo{ Valid: true, - Start: start, Greatest: greatest, - End: end, - Duration: end.Sub(start), MinimumSeparationArcsec: result.MinimumSeparationArcsec, SunSemidiameterArcsec: result.SunSemidiameterArcsec, PlanetSemidiameterArcsec: result.PlanetSemidiameterArcsec, HasInternal: result.HasInternal, + HasExternal: result.HasExternal, + WindowStart: basic.JD2DateByZone(result.WindowStart, loc, false), + WindowEnd: basic.JD2DateByZone(result.WindowEnd, loc, false), + Gamma: result.Gamma, + PenumbralRadius: result.PenumbralRadius, + Series: result.Series, + } + if result.HasExternal { + info.Start = basic.JD2DateByZone(result.ExternalIngress, loc, false) + info.End = basic.JD2DateByZone(result.ExternalEgress, loc, false) + info.Duration = info.End.Sub(info.Start) } if result.HasInternal { info.InternalStart = basic.JD2DateByZone(result.InternalIngress, loc, false) @@ -81,7 +93,7 @@ type LocalTransitInfo = basic.PlanetTransitLocalInfo // LocalTransitContact 站心凌日接触 / topocentric transit contact. type LocalTransitContact = basic.PlanetTransitLocalContactInfo -// LocalTransit 地心凌日 transit 在观测地的站心接触与凌甚 / topocentric contacts of transit seen from the site. +// LocalTransit 凌日 transit 在观测地的站心接触与凌甚 / topocentric contacts of transit seen from the site. // // lon、lat 为度，东经、北纬为正；height 为海拔米。 // lon and lat are in degrees, east and north positive; height is in meters. @@ -89,8 +101,15 @@ func (module *Mercury) LocalTransit(transit TransitInfo, lon, lat, height float6 if !transit.Valid { return LocalTransitInfo{} } - result := basic.PlanetTransitResultFromTimes(1, - transit.Start, transit.InternalStart, transit.Greatest, transit.InternalEnd, transit.End) - local := basic.PlanetTransitLocalWithInstance(result, basic.Observer{Longitude: lon, Latitude: lat, Height: height}, module.inst) + local := basic.PlanetTransitLocalWithInstance(transitResultFromInfo(transit), basic.Observer{Longitude: lon, Latitude: lat, Height: height}, module.inst) return basic.PlanetTransitLocalInfoFromResult(local, transit.Greatest.Location()) } + +// transitResultFromInfo 按一触、四触是否为零选择地心接触或求解窗口还原 basic 结果。 +func transitResultFromInfo(transit TransitInfo) basic.PlanetTransitResult { + if transit.Start.IsZero() && transit.End.IsZero() { + return basic.PlanetTransitResultFromWindow(1, transit.WindowStart, transit.Greatest, transit.WindowEnd) + } + return basic.PlanetTransitResultFromTimes(1, + transit.Start, transit.InternalStart, transit.Greatest, transit.InternalEnd, transit.End) +} diff --git a/venus/api_default.go b/venus/api_default.go index 6f0c1dd..5ad80d4 100644 --- a/venus/api_default.go +++ b/venus/api_default.go @@ -78,7 +78,7 @@ func BrightLimbPositionAngleN(date time.Time, n int) float64 { return defaultModule.BrightLimbPositionAngleN(date, n) } -// ClosestTransit 最近一次地心金星凌日 / closest geocentric Venus transit. +// ClosestTransit 最近一次全球金星凌日 / closest global Venus transit. func ClosestTransit(date time.Time) TransitInfo { return defaultModule.ClosestTransit(date) } // CulminationTime 中天时刻，取 date 所在当地民用日 / culmination time, on the local civil day of date. @@ -178,10 +178,10 @@ func LastSuperiorConjunction(date time.Time) time.Time { return defaultModule.LastSuperiorConjunction(date) } -// LastTransit 上一次地心金星凌日 / previous geocentric Venus transit. +// LastTransit 上一次全球金星凌日 / previous global Venus transit. func LastTransit(date time.Time) TransitInfo { return defaultModule.LastTransit(date) } -// LocalTransit 地心凌日在观测地的站心接触与凌甚 / topocentric contacts of a transit seen from the site. +// LocalTransit 凌日在观测地的站心接触与凌甚 / topocentric contacts of a transit seen from the site. func LocalTransit(transit TransitInfo, lon, lat, height float64) LocalTransitInfo { return defaultModule.LocalTransit(transit, lon, lat, height) } @@ -227,7 +227,7 @@ func NextSuperiorConjunction(date time.Time) time.Time { return defaultModule.NextSuperiorConjunction(date) } -// NextTransit 下一次地心金星凌日 / next geocentric Venus transit. +// NextTransit 下一次全球金星凌日 / next global Venus transit. func NextTransit(date time.Time) TransitInfo { return defaultModule.NextTransit(date) } // ParallacticAngle 金星视差角（天顶方向角） / Venus parallactic angle. diff --git a/venus/transit.go b/venus/transit.go index c0bcd6f..2cbea24 100644 --- a/venus/transit.go +++ b/venus/transit.go @@ -6,12 +6,12 @@ import ( "b612.me/astro/basic" ) -// TransitInfo 地心金星凌日信息 / geocentric Venus transit information. +// TransitInfo 全球金星凌日与地心接触 / global Venus transit and geocentric contacts. // // Start、Greatest、End、InternalStart、InternalEnd 都保持调用者输入的时区。 -// 内切接触不存在时 InternalStart / InternalEnd 为零值。 +// 地心接触不存在时对应时刻为零值；WindowStart、WindowEnd 只是站心求解窗口。 // Start, Greatest, End, InternalStart, and InternalEnd preserve the caller's timezone. -// InternalStart and InternalEnd are zero values when internal contacts do not exist. +// Absent geocentric contacts are zero; WindowStart and WindowEnd are only a local search window. type TransitInfo struct { Valid bool @@ -28,20 +28,26 @@ type TransitInfo struct { SunSemidiameterArcsec float64 PlanetSemidiameterArcsec float64 - HasInternal bool + HasInternal bool + HasExternal bool + WindowStart time.Time + WindowEnd time.Time + Gamma float64 + PenumbralRadius float64 + Series basic.PlanetTransitSeries } -// NextTransit 下一次地心金星凌日 / next geocentric Venus transit. +// NextTransit 下一次全球金星凌日 / next global Venus transit. func (module *Venus) NextTransit(date time.Time) TransitInfo { return transitInfoFromBasic(basic.NextVenusTransitWithInstance(basic.Date2JD(date.UTC()), module.inst), date.Location()) } -// LastTransit 上一次地心金星凌日 / previous geocentric Venus transit. +// LastTransit 上一次全球金星凌日 / previous global Venus transit. func (module *Venus) LastTransit(date time.Time) TransitInfo { return transitInfoFromBasic(basic.LastVenusTransitWithInstance(basic.Date2JD(date.UTC()), module.inst), date.Location()) } -// ClosestTransit 最近一次地心金星凌日 / closest geocentric Venus transit. +// ClosestTransit 最近一次全球金星凌日 / closest global Venus transit. func (module *Venus) ClosestTransit(date time.Time) TransitInfo { return transitInfoFromBasic(basic.ClosestVenusTransitWithInstance(basic.Date2JD(date.UTC()), module.inst), date.Location()) } @@ -50,19 +56,25 @@ func transitInfoFromBasic(result basic.PlanetTransitResult, loc *time.Location) if !result.Valid { return TransitInfo{} } - start := basic.JD2DateByZone(result.ExternalIngress, loc, false) greatest := basic.JD2DateByZone(result.Greatest, loc, false) - end := basic.JD2DateByZone(result.ExternalEgress, loc, false) info := TransitInfo{ Valid: true, - Start: start, Greatest: greatest, - End: end, - Duration: end.Sub(start), MinimumSeparationArcsec: result.MinimumSeparationArcsec, SunSemidiameterArcsec: result.SunSemidiameterArcsec, PlanetSemidiameterArcsec: result.PlanetSemidiameterArcsec, HasInternal: result.HasInternal, + HasExternal: result.HasExternal, + WindowStart: basic.JD2DateByZone(result.WindowStart, loc, false), + WindowEnd: basic.JD2DateByZone(result.WindowEnd, loc, false), + Gamma: result.Gamma, + PenumbralRadius: result.PenumbralRadius, + Series: result.Series, + } + if result.HasExternal { + info.Start = basic.JD2DateByZone(result.ExternalIngress, loc, false) + info.End = basic.JD2DateByZone(result.ExternalEgress, loc, false) + info.Duration = info.End.Sub(info.Start) } if result.HasInternal { info.InternalStart = basic.JD2DateByZone(result.InternalIngress, loc, false) @@ -81,7 +93,7 @@ type LocalTransitInfo = basic.PlanetTransitLocalInfo // LocalTransitContact 站心凌日接触 / topocentric transit contact. type LocalTransitContact = basic.PlanetTransitLocalContactInfo -// LocalTransit 地心凌日 transit 在观测地的站心接触与凌甚 / topocentric contacts of transit seen from the site. +// LocalTransit 凌日 transit 在观测地的站心接触与凌甚 / topocentric contacts of transit seen from the site. // // lon、lat 为度，东经、北纬为正；height 为海拔米。 // lon and lat are in degrees, east and north positive; height is in meters. @@ -89,8 +101,15 @@ func (module *Venus) LocalTransit(transit TransitInfo, lon, lat, height float64) if !transit.Valid { return LocalTransitInfo{} } - result := basic.PlanetTransitResultFromTimes(2, - transit.Start, transit.InternalStart, transit.Greatest, transit.InternalEnd, transit.End) - local := basic.PlanetTransitLocalWithInstance(result, basic.Observer{Longitude: lon, Latitude: lat, Height: height}, module.inst) + local := basic.PlanetTransitLocalWithInstance(transitResultFromInfo(transit), basic.Observer{Longitude: lon, Latitude: lat, Height: height}, module.inst) return basic.PlanetTransitLocalInfoFromResult(local, transit.Greatest.Location()) } + +// transitResultFromInfo 按一触、四触是否为零选择地心接触或求解窗口还原 basic 结果。 +func transitResultFromInfo(transit TransitInfo) basic.PlanetTransitResult { + if transit.Start.IsZero() && transit.End.IsZero() { + return basic.PlanetTransitResultFromWindow(2, transit.WindowStart, transit.Greatest, transit.WindowEnd) + } + return basic.PlanetTransitResultFromTimes(2, + transit.Start, transit.InternalStart, transit.Greatest, transit.InternalEnd, transit.End) +} diff --git a/venus/transit_global_test.go b/venus/transit_global_test.go new file mode 100644 index 0000000..f2e4fae --- /dev/null +++ b/venus/transit_global_test.go @@ -0,0 +1,35 @@ +package venus + +import ( + "testing" + "time" +) + +func TestGlobalTransitWrapper(tester *testing.T) { + zone := time.FixedZone("UTC+8", 8*3600) + info := NextTransit(time.Date(2854, 1, 1, 0, 0, 0, 0, zone)) + if !info.Valid || info.Greatest.Year() != 2854 || info.HasExternal || !info.Start.IsZero() || !info.End.IsZero() || info.Duration != 0 || info.Series.Number != 4 { + tester.Fatalf("invalid global transit: %+v", info) + } + if info.Greatest.Location() != zone || info.WindowStart.Location() != zone || info.WindowEnd.Location() != zone { + tester.Fatal("global transit lost the caller's timezone") + } + local := LocalTransit(info, -135, -60, 0) + if !local.Valid || !local.Visible || local.HasInternal || local.Greatest.Time.Location() != zone { + tester.Fatalf("invalid local partial transit: %+v", local) + } + for _, result := range []TransitInfo{ClosestTransit(info.Greatest), LastTransit(info.Greatest.Add(time.Minute)), NextTransit(info.Greatest.Add(-time.Minute))} { + difference := result.Greatest.Sub(info.Greatest) + if !result.Valid || difference > time.Second || difference < -time.Second { + tester.Fatal("global event lost at a query boundary") + } + } +} + +func TestLocalTransitRequiresContactsOrWindow(tester *testing.T) { + info := NextTransit(time.Date(2854, 1, 1, 0, 0, 0, 0, time.UTC)) + info.WindowStart, info.WindowEnd = time.Time{}, time.Time{} + if local := LocalTransit(info, -135, -60, 0); local.Valid { + tester.Fatal("partial transit without a window must be rejected") + } +}