Replace single DateTime interface with Vector{DateTime}

.gitignore

@@ -15,6 +15,6 @@ node_modules
 .vscode/settings.json
 **.old
 play/Project.toml
-docs/.CondaPkg/meta
+.CondaPkg
 *.mem
 **/Manifest-*.toml

CLAUDE.md

@@ -40,30 +40,35 @@ The package provides a unified interface to multiple solar position algorithms.
 
 ```text
 src/
-  SolarPosition.jl       # Main module, re-exports all submodules
+  SolarPosition.jl       # Main module, re-exports all submodules via @reexport
   Positioning/           # Solar position algorithms
-    Positioning.jl       # Observer struct, SolPos/ApparentSolPos types, solar_position() API
+    Positioning.jl       # Observer, SolPos/ApparentSolPos types, solar_position() dispatch
     psa.jl               # PSA algorithm (default, ±0.0083°)
     noaa.jl, spa.jl,     # Other algorithms (NOAA, SPA, Walraven, USNO)
     walraven.jl, usno.jl
     deltat.jl            # Delta T / leap seconds
   Refraction/            # Atmospheric refraction correction models
-    Refraction.jl        # Abstract base + interface
-    hughes.jl, bennett.jl, sg2.jl, spa.jl, ...
+    Refraction.jl        # RefractionAlgorithm abstract type + refraction() interface
+    hughes.jl, archer.jl, bennett.jl, michalsky.jl, sg2.jl, spa.jl
   Utilities/             # Sunrise/sunset/transit calculations
-    srt.jl, spa.jl
-ext/                     # Weak dependency extensions
-  SolarPositionMakieExt.jl          # Sun path plotting
-  SolarPositionOhMyThreadsExt.jl    # Parallel solar position computation
-  SolarPositionModelingToolkitExt.jl # Symbolic models for MTK
+    srt.jl               # transit_sunrise_sunset(), next_sunrise/sunset/solar_noon, etc.
+    spa.jl               # SPA-based utility helpers
+ext/                     # Weak dependency extensions (loaded via [weakdeps] in Project.toml)
+  SolarPositionMakieExt.jl          # analemmas!() sun path plotting
+  SolarPositionOhMyThreadsExt.jl    # Parallel solar_position via OhMyThreads
+  SolarPositionModelingToolkitExt.jl # SolarPositionBlock for MTK symbolic models
 ```
 
-**Core API pattern:** `solar_position(algorithm, observer, datetime)` returns a `SolPos` or `ApparentSolPos` (if refraction is included). The `Observer` struct holds location (lat/lon/altitude) and optional atmospheric parameters (pressure, temperature).
+**Core API pattern:** `solar_position(observer, datetime, algorithm, refraction)` returns a `SolPos` or `ApparentSolPos` (if refraction model is provided). Also accepts vectors of datetimes (returns `StructVector`) and Tables.jl-compatible tables (adds columns in-place).
 
-**Test discovery:** Test files matching `test-*.jl` under `test/` are automatically discovered and wrapped in `@testset`. Reference values live in `expected-values.jl` files alongside algorithm tests.
+**Adding a new algorithm:** Subtype `SolarAlgorithm`, implement `_solar_position(obs, dt, alg::YourAlg)::SolPos`, and the dispatch in `Positioning.jl` handles refraction wrapping automatically. Same pattern for refraction: subtype `RefractionAlgorithm` and implement `_refraction(model, elevation)`.
+
+**Test discovery:** Test files matching `test-*.jl` under `test/` are automatically discovered and wrapped in `@testset`. Don't add tests to `runtests.jl` directly. Reference values live in `expected-values.jl` files alongside algorithm tests.
 
 ## Code Style
 
-- Formatter: **Runic.jl** (enforced via pre-commit). Run before committing.
-- Imports: All used symbols must be explicitly imported (checked by ExplicitImports.jl).
+- Formatter: **Runic.jl** (enforced via pre-commit). Run `pre-commit run runic --all-files` before committing.
+- Imports: All used symbols must be explicitly imported (checked by ExplicitImports.jl pre-commit hook).
+- Spelling: **typos** is enforced via pre-commit. If a false positive occurs, add it to `_typos.toml`.
 - Package quality is checked with **Aqua.jl** and type inference with **JET.jl** (Julia 1.12 only).
+- Pre-commit `no-commit-to-branch` blocks direct commits to `main`; work on feature branches.

src/Positioning/Positioning.jl

@@ -263,18 +263,77 @@ See also: [`solar_position!`](@ref), [`Observer`](@ref), [`PSA`](@ref), [`NOAA`]
 """
 function solar_position end
 
+# Resolve DefaultRefraction to a concrete refraction algorithm.
+# Each algorithm file overrides this for its own DefaultRefraction mapping.
+# The observer type T is passed to allow type-parameterized refraction models.
+_resolve_refraction(::SolarAlgorithm, r::RefractionAlgorithm, ::Type{T}) where {T} = r
+
+# --- Vectorized internal interface ---
+
+# Generic fallback: loops over scalar _solar_position for any algorithm that only
+# defines the scalar method. Built-in algorithms override with optimized versions.
+function _solar_position!(
+        pos::StructArrays.StructVector{SolPos{T}},
+        obs::AbstractObserver{T},
+        dts::AbstractVector{DateTime},
+        alg::SolarAlgorithm,
+    ) where {T}
+    @inbounds for i in eachindex(dts, pos)
+        pos[i] = _solar_position(obs, dts[i], alg)
+    end
+    return pos
+end
+
+# NoRefraction: just compute raw positions
+function _solar_position!(
+        pos::StructArrays.StructVector{SolPos{T}},
+        obs::AbstractObserver{T},
+        dts::AbstractVector{DateTime},
+        alg::SolarAlgorithm,
+        ::NoRefraction,
+    ) where {T}
+    return _solar_position!(pos, obs, dts, alg)
+end
+
+# With refraction: compute raw positions into shared columns, then apply refraction
+function _solar_position!(
+        pos::StructArrays.StructVector{ApparentSolPos{T}},
+        obs::AbstractObserver{T},
+        dts::AbstractVector{DateTime},
+        alg::SolarAlgorithm,
+        refraction::RefractionAlgorithm,
+    ) where {T}
+    # Create SolPos view backed by the same column arrays (no heap allocation)
+    raw = StructArrays.StructVector{SolPos{T}}(
+        (
+            azimuth = pos.azimuth,
+            elevation = pos.elevation,
+            zenith = pos.zenith,
+        )
+    )
+    _solar_position!(raw, obs, dts, alg)
+
+    # Apply refraction in second pass
+    @inbounds for i in eachindex(pos)
+        correction = Refraction.refraction(refraction, pos.elevation[i])
+        pos.apparent_elevation[i] = pos.elevation[i] + correction
+        pos.apparent_zenith[i] = T(90) - pos.apparent_elevation[i]
+    end
+    return pos
+end
+
+# --- Scalar internal interface ---
+
+# Scalar refraction wrappers (zero-allocation path for single datetimes)
 function _solar_position(obs, dt, alg::SolarAlgorithm, ::NoRefraction)
     return _solar_position(obs, dt, alg)
 end
 
 function _solar_position(obs, dt, alg::SolarAlgorithm, refraction::RefractionAlgorithm)
     pos = _solar_position(obs, dt, alg)
-
-    # apply refraction correction
-    refraction_correction_deg = Refraction.refraction(refraction, pos.elevation)
-    apparent_elevation_deg = pos.elevation + refraction_correction_deg
+    correction = Refraction.refraction(refraction, pos.elevation)
+    apparent_elevation_deg = pos.elevation + correction
     apparent_zenith_deg = 90 - apparent_elevation_deg
-
     return ApparentSolPos(
         pos.azimuth,
         pos.elevation,
@@ -284,13 +343,15 @@ function _solar_position(obs, dt, alg::SolarAlgorithm, refraction::RefractionAlg
     )
 end
 
+# Scalar public API (zero-allocation)
 function solar_position(
         obs::AbstractObserver{T},
         dt::DateTime,
         alg::SolarAlgorithm = PSA(),
         refraction::RefractionAlgorithm = DefaultRefraction(),
     ) where {T <: AbstractFloat}
-    return _solar_position(obs, dt, alg, refraction)
+    resolved = _resolve_refraction(alg, refraction, T)
+    return _solar_position(obs, dt, alg, resolved)
 end
 
 function solar_position(
@@ -303,77 +364,49 @@ function solar_position(
 end
 
 function solar_position!(
-        pos::StructArrays.StructVector{T},
-        obs::AbstractObserver,
-        dts::AbstractVector{Union{DateTime, ZonedDateTime}},
-        alg::SolarAlgorithm = PSA(),
-        refraction::RefractionAlgorithm = DefaultRefraction(),
-    ) where {T <: AbstractSolPos}
-    @inbounds for i in eachindex(dts, pos)
-        pos[i] = solar_position(obs, dts[i], alg, refraction)
-    end
-    return pos
-end
-
-function solar_position!(
-        pos::StructArrays.StructVector{T},
-        obs::AbstractObserver,
+        pos::StructArrays.StructVector{<:AbstractSolPos},
+        obs::AbstractObserver{T},
         dts::AbstractVector{DateTime},
         alg::SolarAlgorithm = PSA(),
         refraction::RefractionAlgorithm = DefaultRefraction(),
-    ) where {T <: AbstractSolPos}
-    @inbounds for i in eachindex(dts, pos)
-        pos[i] = solar_position(obs, dts[i], alg, refraction)
-    end
+    ) where {T}
+    resolved = _resolve_refraction(alg, refraction, T)
+    _solar_position!(pos, obs, dts, alg, resolved)
     return pos
 end
 
 function solar_position!(
-        pos::StructArrays.StructVector{T},
-        obs::AbstractObserver,
+        pos::StructArrays.StructVector{<:AbstractSolPos},
+        obs::AbstractObserver{T},
         dts::AbstractVector{ZonedDateTime},
         alg::SolarAlgorithm = PSA(),
         refraction::RefractionAlgorithm = DefaultRefraction(),
-    ) where {T <: AbstractSolPos}
-    @inbounds for i in eachindex(dts, pos)
-        pos[i] = solar_position(obs, dts[i], alg, refraction)
-    end
-    return pos
+    ) where {T}
+    utc_dts = DateTime.(dts, Ref(UTC))
+    return solar_position!(pos, obs, utc_dts, alg, refraction)
 end
 
-function solar_position(
+function solar_position!(
+        pos::StructArrays.StructVector{<:AbstractSolPos},
         obs::AbstractObserver{T},
         dts::AbstractVector{Union{DateTime, ZonedDateTime}},
         alg::SolarAlgorithm = PSA(),
         refraction::RefractionAlgorithm = DefaultRefraction(),
-    ) where {T <: AbstractFloat}
-    RetType = result_type(typeof(alg), typeof(refraction), T)
-    pos = StructArrays.StructVector{RetType}(undef, length(dts))
-    solar_position!(pos, obs, dts, alg, refraction)
-    return pos
-end
-
-function solar_position(
-        obs::AbstractObserver{T},
-        dts::AbstractVector{DateTime},
-        alg::SolarAlgorithm = PSA(),
-        refraction::RefractionAlgorithm = DefaultRefraction(),
-    ) where {T <: AbstractFloat}
-    RetType = result_type(typeof(alg), typeof(refraction), T)
-    pos = StructArrays.StructVector{RetType}(undef, length(dts))
-    solar_position!(pos, obs, dts, alg, refraction)
-    return pos
+    ) where {T}
+    utc_dts = [dt isa ZonedDateTime ? DateTime(dt, UTC) : dt for dt in dts]
+    return solar_position!(pos, obs, utc_dts, alg, refraction)
 end
 
 function solar_position(
         obs::AbstractObserver{T},
-        dts::AbstractVector{ZonedDateTime},
+        dts::AbstractVector{<:Union{DateTime, ZonedDateTime}},
         alg::SolarAlgorithm = PSA(),
         refraction::RefractionAlgorithm = DefaultRefraction(),
     ) where {T <: AbstractFloat}
-    RetType = result_type(typeof(alg), typeof(refraction), T)
+    resolved = _resolve_refraction(alg, refraction, T)
+    RetType = result_type(typeof(alg), typeof(resolved), T)
     pos = StructArrays.StructVector{RetType}(undef, length(dts))
-    solar_position!(pos, obs, dts, alg, refraction)
+    solar_position!(pos, obs, dts, alg, resolved)
     return pos
 end
 

src/Positioning/noaa.jl

@@ -33,14 +33,10 @@ end
 
 NOAA() = NOAA(67.0)  # default delta_t value (2020 default from pvlib)
 
-
-function _solar_position(obs::Observer{T}, dt::DateTime, alg::NOAA) where {T}
-    δt = if alg.delta_t === nothing
-        calculate_deltat(dt)
-    else
-        alg.delta_t
-    end
-
+# Core computation shared by scalar and vectorized paths
+@inline function _noaa_kernel(
+        dt::DateTime, sin_lat, cos_lat, longitude, ::Type{T},
+    ) where {T}
     # convert to Julian date and Julian century
     jd = datetime2julian(dt)
     jc = (jd - 2451545.0) / 36525.0
@@ -51,7 +47,7 @@ function _solar_position(obs::Observer{T}, dt::DateTime, alg::NOAA) where {T}
     # mean anomaly [degrees]
     mean_anom = 357.52911 + jc * (35999.05029 - 0.0001537 * jc)
 
-    # cccentricity of Earth's orbit
+    # eccentricity of Earth's orbit
     eccent = 0.016708634 - jc * (0.000042037 + 0.0000001267 * jc)
 
     # sun equation of center [degrees]
@@ -68,7 +64,8 @@ function _solar_position(obs::Observer{T}, dt::DateTime, alg::NOAA) where {T}
     # mean obliquity of ecliptic [degrees]
     mean_obliq =
         23.0 +
-        (26.0 + (21.448 - jc * (46.815 + jc * (0.00059 - jc * 0.001813))) / 60.0) / 60.0
+        (26.0 + (21.448 - jc * (46.815 + jc * (0.00059 - jc * 0.001813))) / 60.0) /
+        60.0
 
     # obliquity correction [degrees]
     obliq_corr = mean_obliq + 0.00256 * cosd(125.04 - 1934.136 * jc)
@@ -80,27 +77,26 @@ function _solar_position(obs::Observer{T}, dt::DateTime, alg::NOAA) where {T}
         4.0 * rad2deg(
         var_y * sind(2.0 * mean_long) - 2.0 * eccent * sind(mean_anom) +
             4.0 * eccent * var_y * sind(mean_anom) * cosd(2.0 * mean_long) -
-            0.5 * var_y^2 * sind(4.0 * mean_long) - 1.25 * eccent^2 * sind(2.0 * mean_anom),
+            0.5 * var_y^2 * sind(4.0 * mean_long) -
+            1.25 * eccent^2 * sind(2.0 * mean_anom),
     )
 
     # true solar time [minutes]
     hour_frac = fractional_hour(dt)
     minutes = hour_frac * 60.0
-    true_solar_time = mod(minutes + eot + 4.0 * obs.longitude, 1440.0)
+    true_solar_time = mod(minutes + eot + 4.0 * longitude, 1440.0)
 
     # hour angle [degrees]
-    # true_solar_time is in [0, 1440) minutes, so true_solar_time/4 is in [0, 360) degrees
-    # Convert to standard hour angle range (-180, 180] where 0 is solar noon
     hour_angle = true_solar_time / 4.0 - 180.0
 
     # zenith angle [degrees]
     zenith = acosd(
-        obs.sin_lat * sind(sun_declin) + obs.cos_lat * cosd(sun_declin) * cosd(hour_angle),
+        sin_lat * sind(sun_declin) + cos_lat * cosd(sun_declin) * cosd(hour_angle),
     )
 
     # azimuth angle [degrees]
-    azimuth_numerator = obs.sin_lat * cosd(zenith) - sind(sun_declin)
-    azimuth_denominator = obs.cos_lat * sind(zenith)
+    azimuth_numerator = sin_lat * cosd(zenith) - sind(sun_declin)
+    azimuth_denominator = cos_lat * sind(zenith)
 
     azimuth = if hour_angle > 0.0
         mod(acosd(azimuth_numerator / azimuth_denominator) + 180.0, 360.0)
@@ -111,10 +107,28 @@ function _solar_position(obs::Observer{T}, dt::DateTime, alg::NOAA) where {T}
     return SolPos{T}(azimuth, 90.0 - zenith, zenith)
 end
 
-function _solar_position(obs, dt, alg::NOAA, ::DefaultRefraction)
-    return _solar_position(obs, dt, alg, HUGHES())
+function _solar_position(obs::Observer{T}, dt::DateTime, alg::NOAA) where {T}
+    return _noaa_kernel(dt, obs.sin_lat, obs.cos_lat, obs.longitude, T)
 end
 
+function _solar_position!(
+        pos::StructArrays.StructVector{SolPos{T}},
+        obs::Observer{T},
+        dts::AbstractVector{DateTime},
+        alg::NOAA,
+    ) where {T}
+    sin_lat = obs.sin_lat
+    cos_lat = obs.cos_lat
+    longitude = obs.longitude
+
+    @inbounds for i in eachindex(dts, pos)
+        pos[i] = _noaa_kernel(dts[i], sin_lat, cos_lat, longitude, T)
+    end
+    return pos
+end
+
+_resolve_refraction(::NOAA, ::DefaultRefraction, ::Type{T}) where {T} = HUGHES()
+
 # NOAA with DefaultRefraction returns ApparentSolPos (uses HUGHES refraction)
 result_type(::Type{NOAA}, ::Type{DefaultRefraction}, ::Type{T}) where {T} =
     ApparentSolPos{T}