214 lines
7.8 KiB
JavaScript
214 lines
7.8 KiB
JavaScript
const RAD = Math.PI / 180
|
|
const DEG = 180 / Math.PI
|
|
|
|
function jdn(date) {
|
|
return (date instanceof Date ? date : new Date(date)) / 86400000 + 2440587.5
|
|
}
|
|
|
|
function jdnToDate(jd) {
|
|
return new Date((jd - 2440587.5) * 86400000)
|
|
}
|
|
|
|
function sunPosition(lat, lon, date = new Date()) {
|
|
const d = jdn(date) - 2451545.0
|
|
const L = (280.46 + 0.9856474 * d) % 360
|
|
const g = (357.528 + 0.9856003 * d) % 360
|
|
const lam = L + 1.915 * Math.sin(g * RAD) + 0.02 * Math.sin(2 * g * RAD)
|
|
const eps = 23.439 - 0.0000004 * d
|
|
const sinL = Math.sin(lam * RAD)
|
|
const ra = Math.atan2(Math.cos(eps * RAD) * sinL, Math.cos(lam * RAD)) * DEG
|
|
const dec = Math.asin(Math.sin(eps * RAD) * sinL) * DEG
|
|
const UT = date.getUTCHours() + date.getUTCMinutes() / 60 + date.getUTCSeconds() / 3600
|
|
const GMST = (6.697375 + 0.0657098242 * d + UT) % 24
|
|
const LMST = (GMST + lon / 15) % 24
|
|
const ha = LMST * 15 - ra
|
|
const elev = Math.asin(
|
|
Math.sin(lat * RAD) * Math.sin(dec * RAD) +
|
|
Math.cos(lat * RAD) * Math.cos(dec * RAD) * Math.cos(ha * RAD)
|
|
) * DEG
|
|
const az = Math.atan2(
|
|
-Math.sin(ha * RAD),
|
|
Math.tan(dec * RAD) * Math.cos(lat * RAD) - Math.sin(lat * RAD) * Math.cos(ha * RAD)
|
|
) * DEG
|
|
return {
|
|
sun_elevation: Math.round(elev * 10) / 10,
|
|
sun_azimuth: Math.round(((az + 360) % 360) * 10) / 10,
|
|
}
|
|
}
|
|
|
|
function sunriseSunset(lat, lon, date = new Date()) {
|
|
const d = Math.floor(jdn(date)) - 2451545
|
|
const noon = 2451545 + 0.0009 + ((-lon) / 360) + Math.round(d - (-lon) / 360)
|
|
const M = (357.5291 + 0.98560028 * (noon - 2451545)) % 360
|
|
const C = 1.9148 * Math.sin(M * RAD) + 0.02 * Math.sin(2 * M * RAD) + 0.0003 * Math.sin(3 * M * RAD)
|
|
const lam = (M + C + 180 + 102.9372) % 360
|
|
const jnoon = noon + 0.0053 * Math.sin(M * RAD) - 0.0069 * Math.sin(2 * lam * RAD)
|
|
const dec = Math.asin(Math.sin(23.4397 * RAD) * Math.sin(lam * RAD)) * DEG
|
|
const cosH = (Math.sin(-0.8333 * RAD) - Math.sin(lat * RAD) * Math.sin(dec * RAD)) / (Math.cos(lat * RAD) * Math.cos(dec * RAD))
|
|
|
|
if (Math.abs(cosH) > 1) {
|
|
return { sun_sunrise: null, sun_sunset: null, sun_solar_noon: jdnToDate(jnoon).toISOString(), sun_polar: cosH < -1 ? 'day' : 'night' }
|
|
}
|
|
|
|
const H = Math.acos(cosH) * DEG
|
|
const rise = jdnToDate(jnoon - H / 360)
|
|
const set = jdnToDate(jnoon + H / 360)
|
|
const noon_d = jdnToDate(jnoon)
|
|
|
|
// Round to nearest second, drop milliseconds
|
|
const fmt = d => new Date(Math.round(d.getTime() / 1000) * 1000).toISOString()
|
|
|
|
return {
|
|
sun_sunrise: fmt(rise),
|
|
sun_sunset: fmt(set),
|
|
sun_solar_noon: fmt(noon_d),
|
|
sun_day_length_hours: Math.round((set - rise) / 36000) / 100,
|
|
sun_golden_hour_morning_start: fmt(new Date(rise - 30 * 60000)),
|
|
sun_golden_hour_morning_end: fmt(new Date(rise + 30 * 60000)),
|
|
sun_golden_hour_evening_start: fmt(new Date(set - 30 * 60000)),
|
|
sun_golden_hour_evening_end: fmt(new Date(set + 30 * 60000)),
|
|
}
|
|
}
|
|
|
|
// Smooth sunspot number approximation from solar flux F10.7
|
|
// SSN ≈ 1.61 * F10.7 - 63.7 (linear regression, valid for F10.7 > 70)
|
|
export function ssnFromFlux(f107) {
|
|
if (!f107) return null
|
|
return Math.max(0, Math.round(1.61 * f107 - 63.7))
|
|
}
|
|
|
|
export function sunActivityLabel(f107) {
|
|
if (!f107) return null
|
|
if (f107 < 80) return 'Very Low'
|
|
if (f107 < 100) return 'Low'
|
|
if (f107 < 150) return 'Moderate'
|
|
if (f107 < 200) return 'High'
|
|
return 'Very High'
|
|
}
|
|
|
|
function moonPhase(date = new Date()) {
|
|
const jd = jdn(date)
|
|
const cycle = 29.53058867
|
|
const known = 2451550.1
|
|
const phase = ((jd - known) % cycle + cycle) % cycle
|
|
const illum = Math.round((1 - Math.cos(phase / cycle * 2 * Math.PI)) / 2 * 100)
|
|
let name
|
|
if (phase < 1.85) name = 'New Moon'
|
|
else if (phase < 7.38) name = 'Waxing Crescent'
|
|
else if (phase < 9.22) name = 'First Quarter'
|
|
else if (phase < 14.76) name = 'Waxing Gibbous'
|
|
else if (phase < 16.61) name = 'Full Moon'
|
|
else if (phase < 22.15) name = 'Waning Gibbous'
|
|
else if (phase < 23.99) name = 'Last Quarter'
|
|
else name = 'Waning Crescent'
|
|
return { moon_phase: name, moon_illumination_pct: illum }
|
|
}
|
|
|
|
function nextMoonEvents(date = new Date()) {
|
|
const cycle = 29.53058867
|
|
const jd = jdn(date)
|
|
const known = 2451550.1
|
|
const phase = ((jd - known) % cycle + cycle) % cycle
|
|
const toNew = (cycle - phase) % cycle
|
|
const toFull = phase < 14.76 ? 14.76 - phase : cycle - phase + 14.76
|
|
return {
|
|
moon_next_new: jdnToDate(jd + toNew).toISOString(),
|
|
moon_next_full: jdnToDate(jd + toFull).toISOString(),
|
|
}
|
|
}
|
|
|
|
function moonriseMoonset(lat, lon, date = new Date()) {
|
|
const base = new Date(date)
|
|
base.setUTCHours(0, 0, 0, 0)
|
|
let prev = null
|
|
let moonrise = null
|
|
let moonset = null
|
|
|
|
// Scan 48h in 10min steps — covers edge cases where moon rises/sets next UTC day
|
|
for (let m = 0; m <= 2880; m += 10) {
|
|
const t = new Date(base.getTime() + m * 60000)
|
|
const d = jdn(t) - 2451545
|
|
const L = (218.316 + 13.176396 * d) % 360
|
|
const M = (134.963 + 13.064993 * d) % 360
|
|
const F = (93.272 + 13.229350 * d) % 360
|
|
const lon_ = L + 6.289 * Math.sin(M * RAD)
|
|
const b = 5.128 * Math.sin(F * RAD)
|
|
const dec = Math.asin(
|
|
Math.sin(b * RAD) * Math.cos(23.4397 * RAD) +
|
|
Math.cos(b * RAD) * Math.sin(23.4397 * RAD) * Math.sin(lon_ * RAD)
|
|
) * DEG
|
|
const GMST = (6.697375 + 0.0657098242 * d + (t.getUTCHours() + t.getUTCMinutes() / 60)) % 24
|
|
const LMST = (GMST + lon / 15) % 24
|
|
const ra = Math.atan2(
|
|
Math.sin(lon_ * RAD) * Math.cos(23.4397 * RAD) - Math.tan(b * RAD) * Math.sin(23.4397 * RAD),
|
|
Math.cos(lon_ * RAD)
|
|
) * DEG
|
|
const ha = LMST * 15 - ra
|
|
const elev = Math.asin(
|
|
Math.sin(lat * RAD) * Math.sin(dec * RAD) +
|
|
Math.cos(lat * RAD) * Math.cos(dec * RAD) * Math.cos(ha * RAD)
|
|
) * DEG
|
|
|
|
if (prev !== null) {
|
|
if (prev < 0 && elev >= 0 && !moonrise) moonrise = new Date(t.getTime() - 5 * 60000).toISOString()
|
|
if (prev >= 0 && elev < 0 && !moonset) moonset = new Date(t.getTime() - 5 * 60000).toISOString()
|
|
}
|
|
prev = elev
|
|
if (moonrise && moonset) break
|
|
}
|
|
|
|
return { moon_moonrise: moonrise, moon_moonset: moonset }
|
|
}
|
|
|
|
function nextSolsticeEquinox(date = new Date()) {
|
|
const y = date.getFullYear()
|
|
const events = [
|
|
{ name: 'March Equinox', date: new Date(Date.UTC(y, 2, 20)) },
|
|
{ name: 'June Solstice', date: new Date(Date.UTC(y, 5, 21)) },
|
|
{ name: 'September Equinox', date: new Date(Date.UTC(y, 8, 22)) },
|
|
{ name: 'December Solstice', date: new Date(Date.UTC(y, 11, 21)) },
|
|
{ name: 'March Equinox', date: new Date(Date.UTC(y+1, 2, 20)) },
|
|
]
|
|
const next = events.find(e => e.date > date)
|
|
const prev = [...events].reverse().find(e => e.date <= date)
|
|
return {
|
|
season_next_event: next.name,
|
|
season_next_event_date: next.date.toISOString(),
|
|
season_last_event: prev?.name || null,
|
|
season_last_event_date: prev?.date.toISOString() || null,
|
|
}
|
|
}
|
|
|
|
export function getCelestialCurrent(lat, lon, date = new Date()) {
|
|
return {
|
|
...sunPosition(lat, lon, date),
|
|
...sunriseSunset(lat, lon, date),
|
|
...moonPhase(date),
|
|
...nextMoonEvents(date),
|
|
...moonriseMoonset(lat, lon, date),
|
|
...nextSolsticeEquinox(date),
|
|
}
|
|
}
|
|
|
|
export function getCelestialHourly(lat, lon, hours) {
|
|
return hours.map(({ time }) => {
|
|
const d = new Date(time)
|
|
const pos = sunPosition(lat, lon, d)
|
|
const phase = moonPhase(d)
|
|
return { time, ...pos, moon_illumination_pct: phase.moon_illumination_pct, moon_phase: phase.moon_phase }
|
|
})
|
|
}
|
|
|
|
export function getCelestialDaily(lat, lon, days) {
|
|
return days.map(({ time }) => {
|
|
const d = new Date(time)
|
|
return {
|
|
time,
|
|
...sunriseSunset(lat, lon, d),
|
|
...moonPhase(d),
|
|
...nextMoonEvents(d),
|
|
...moonriseMoonset(lat, lon, d),
|
|
}
|
|
})
|
|
}
|