SPA update

.gitignore

@@ -16,3 +16,4 @@ node_modules
 **.old
 play/Project.toml
 docs/.CondaPkg/meta
+*.mem

docs/src/reference.md

@@ -24,6 +24,7 @@ SolarPosition.Positioning.solar_position!
 ## Observer and Position Types
 
 ```@docs
+SolarPosition.Positioning.AbstractObserver
 SolarPosition.Positioning.Observer
 SolarPosition.Positioning.SolPos
 SolarPosition.Positioning.ApparentSolPos

src/Positioning/Positioning.jl

@@ -14,9 +14,8 @@ using TimeZones: ZonedDateTime, UTC
 using StructArrays: StructArrays
 using Tables: Tables
 using DocStringExtensions: TYPEDFIELDS, TYPEDEF, TYPEDSIGNATURES
-using Reexport: @reexport
 import ..Refraction
-using ..Refraction: RefractionAlgorithm, NoRefraction, DefaultRefraction
+using ..Refraction: RefractionAlgorithm, NoRefraction, DefaultRefraction, SPARefraction
 
 """
     $(TYPEDEF)
@@ -32,6 +31,22 @@ struct MyAlgorithm <: SolarAlgorithm end
 """
 abstract type SolarAlgorithm end
 
+"""
+    $(TYPEDEF)
+
+Abstract base type for Observers.
+
+All concrete observer types must inherit from this type. Most algorithms will use the
+default `Observer` struct defined in this module. More specialized observer types can be
+defined by inheriting from this abstract type, see for example [`SPAObserver`](@ref).
+
+# Examples
+```julia
+struct MyObserver <: AbstractObserver{Float64} end
+```
+"""
+abstract type AbstractObserver{T<:AbstractFloat} end
+
 """
     $(TYPEDEF)
 
@@ -59,7 +74,7 @@ The `horizon` parameter can be specified as:
 - A number in degrees (e.g., `0.5667`)
 - A `degrees=>arcminutes` pair (e.g., `0=>34` for 34 arcminutes = 0.5667°)
 """
-struct Observer{T<:AbstractFloat}
+struct Observer{T<:AbstractFloat} <: AbstractObserver{T}
     "Geodetic latitude (+N)"
     latitude::T
     "Longitude (+E)"
@@ -172,7 +187,7 @@ an observer's geographic location and timestamp(s). It supports multiple input f
 automatically handles time zone conversions.
 
 # Arguments
-- `obs::Observer`: Observer location with latitude, longitude, and altitude
+- `obs::AbstractObserver`: Observer location with latitude, longitude, and altitude
 - `dt::DateTime` or `dt::ZonedDateTime`: Single timestamp
 - `dts::AbstractVector`: Vector of timestamps (DateTime or ZonedDateTime)
 - `alg::SolarAlgorithm`: Solar positioning algorithm (default: `PSA()`)
@@ -258,7 +273,7 @@ function _solar_position(obs, dt, alg::SolarAlgorithm, refraction::RefractionAlg
     # apply refraction correction
     refraction_correction_deg = Refraction.refraction(refraction, pos.elevation)
     apparent_elevation_deg = pos.elevation + refraction_correction_deg
-    apparent_zenith_deg = 90.0 - apparent_elevation_deg
+    apparent_zenith_deg = 90 - apparent_elevation_deg
 
     return ApparentSolPos(
         pos.azimuth,
@@ -270,26 +285,26 @@ function _solar_position(obs, dt, alg::SolarAlgorithm, refraction::RefractionAlg
 end
 
 function solar_position(
-    obs::Observer{T},
+    obs::AbstractObserver{T},
     dt::DateTime,
     alg::SolarAlgorithm = PSA(),
     refraction::RefractionAlgorithm = DefaultRefraction(),
 ) where {T<:AbstractFloat}
-    _solar_position(obs, dt, alg, refraction)
+    return _solar_position(obs, dt, alg, refraction)
 end
 
 function solar_position(
-    obs::Observer{T},
+    obs::AbstractObserver{T},
     dt::ZonedDateTime,
     alg::SolarAlgorithm = PSA(),
     refraction::RefractionAlgorithm = DefaultRefraction(),
 ) where {T<:AbstractFloat}
-    solar_position(obs, DateTime(dt, UTC), alg, refraction)
+    return solar_position(obs, DateTime(dt, UTC), alg, refraction)
 end
 
 function solar_position!(
     pos::StructArrays.StructVector{T},
-    obs::Observer,
+    obs::AbstractObserver,
     dts::AbstractVector{Union{DateTime,ZonedDateTime}},
     alg::SolarAlgorithm = PSA(),
     refraction::RefractionAlgorithm = DefaultRefraction(),
@@ -302,7 +317,7 @@ end
 
 function solar_position!(
     pos::StructArrays.StructVector{T},
-    obs::Observer,
+    obs::AbstractObserver,
     dts::AbstractVector{DateTime},
     alg::SolarAlgorithm = PSA(),
     refraction::RefractionAlgorithm = DefaultRefraction(),
@@ -315,7 +330,7 @@ end
 
 function solar_position!(
     pos::StructArrays.StructVector{T},
-    obs::Observer,
+    obs::AbstractObserver,
     dts::AbstractVector{ZonedDateTime},
     alg::SolarAlgorithm = PSA(),
     refraction::RefractionAlgorithm = DefaultRefraction(),
@@ -327,39 +342,39 @@ function solar_position!(
 end
 
 function solar_position(
-    obs::Observer{T},
+    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)
-    pos
+    return pos
 end
 
 function solar_position(
-    obs::Observer{T},
+    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)
-    pos
+    return pos
 end
 
 function solar_position(
-    obs::Observer{T},
+    obs::AbstractObserver{T},
     dts::AbstractVector{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)
-    pos
+    return pos
 end
 
 
@@ -370,7 +385,7 @@ Compute solar positions for all times in a table and add the results as new colu
 
 # Arguments
 - `table` : Table-like object with datetime column (must support Tables.jl interface).
-- `obs::Observer` : Observer location (latitude, longitude, altitude).
+- `obs::AbstractObserver` : Observer location (latitude, longitude, altitude).
 - `latitude, longitude, altitude` : Specify observer location directly.
 - `dt_col::Symbol` : Name of the datetime column (default: `:datetime`).
 - `alg::SolarAlgorithm` : Algorithm to use (default: `PSA()`).
@@ -386,7 +401,7 @@ The input table is modified **in-place** by adding new columns.
 """
 function solar_position!(
     table,
-    obs::Observer{T},
+    obs::AbstractObserver{T},
     alg::SolarAlgorithm = PSA(),
     refraction::RefractionAlgorithm = DefaultRefraction();
     dt_col::Symbol = :datetime,
@@ -414,7 +429,7 @@ See [`solar_position!`](@ref) for detailed documentation of arguments, examples,
 """
 function solar_position(
     table,
-    obs::Observer{T},
+    obs::AbstractObserver{T},
     alg::SolarAlgorithm = PSA(),
     refraction::RefractionAlgorithm = DefaultRefraction();
     kwargs...,

src/Positioning/spa.jl

@@ -65,7 +65,7 @@ multiple times at the same location.
 # Internal Fields
 $(TYPEDFIELDS)
 """
-struct SPAObserver{T<:AbstractFloat}
+struct SPAObserver{T<:AbstractFloat} <: AbstractObserver{T}
     "Geodetic latitude (+N)"
     latitude::T
     "Longitude (+E)"
@@ -448,71 +448,44 @@ function _solar_position(
     return SolPos{T}(az, e0, θz0)
 end
 
-# SPA with DefaultRefraction applies SPA's built-in atmospheric refraction model
-function _solar_position(
-    obs::Observer{T},
-    dt::DateTime,
-    alg::SPA,
-    ::Refraction.DefaultRefraction,
-) where {T<:AbstractFloat}
+function _solar_position(obs::Observer{T}, dt::DateTime, alg::SPA) where {T<:AbstractFloat}
     spa_obs = SPAObserver{T}(obs.latitude, obs.longitude, obs.altitude)
-
-    # Get base position without refraction
-    base_pos = _solar_position(spa_obs, dt, alg)
-
-    # Apply SPA's atmospheric refraction correction
-    Δe = atmospheric_refraction_correction(
-        alg.pressure,
-        alg.temperature,
-        base_pos.elevation,
-        alg.atmos_refract,
-    )
-
-    # Calculate apparent positions
-    e = base_pos.elevation + Δe  # apparent elevation
-    θz = 90.0 - e  # apparent zenith
-
-    return ApparentSolPos{T}(base_pos.azimuth, base_pos.elevation, base_pos.zenith, e, θz)
+    return _solar_position(spa_obs, dt, alg)
 end
 
-# SPA with NoRefraction returns SolPos (no refraction applied)
 function _solar_position(
-    obs::Observer{T},
-    dt::DateTime,
+    obs::AbstractObserver{T},
+    dt,
     alg::SPA,
-    ::Refraction.NoRefraction,
+    ::DefaultRefraction,
 ) where {T<:AbstractFloat}
-    spa_obs = SPAObserver{T}(obs.latitude, obs.longitude, obs.altitude)
-    return _solar_position(spa_obs, dt, alg)
+    return _solar_position(
+        obs,
+        dt,
+        alg,
+        SPARefraction{T}(
+            pressure = T(alg.pressure),
+            temperature = T(alg.temperature),
+            atmos_refract = T(alg.atmos_refract),
+        ),
+    )
 end
 
-# SPA with external refraction algorithm: apply the external refraction to base position
-function _solar_position(
-    obs::Observer{T},
-    dt::DateTime,
+function solar_position!(
+    pos::StructArrays.StructVector{S},
+    obs::AbstractObserver{T},
+    dts::AbstractVector{DateTime},
     alg::SPA,
-    refraction::Refraction.RefractionAlgorithm,
-) where {T<:AbstractFloat}
+    refraction::RefractionAlgorithm = DefaultRefraction(),
+) where {S<:AbstractSolPos,T<:AbstractFloat}
     spa_obs = SPAObserver{T}(obs.latitude, obs.longitude, obs.altitude)
-
-    # Get base position without refraction
-    base_pos = _solar_position(spa_obs, dt, alg)
-
-    # Apply external refraction correction
-    refraction_correction_deg = Refraction.refraction(refraction, base_pos.elevation)
-    apparent_elevation_deg = base_pos.elevation + refraction_correction_deg
-    apparent_zenith_deg = 90.0 - apparent_elevation_deg
-
-    return ApparentSolPos{T}(
-        base_pos.azimuth,
-        base_pos.elevation,
-        base_pos.zenith,
-        apparent_elevation_deg,
-        apparent_zenith_deg,
-    )
+    @inbounds for i in eachindex(dts, pos)
+        pos[i] = solar_position(spa_obs, dts[i], alg, refraction)
+    end
+    return pos
 end
 
 # SPA returns SolPos with NoRefraction, ApparentSolPos with any refraction
 result_type(::Type{SPA}, ::Type{NoRefraction}, ::Type{T}) where {T} = SolPos{T}
-result_type(::Type{SPA}, ::Type{<:Refraction.RefractionAlgorithm}, ::Type{T}) where {T} =
+result_type(::Type{SPA}, ::Type{<:RefractionAlgorithm}, ::Type{T}) where {T} =
     ApparentSolPos{T}

src/Positioning/utils.jl

@@ -1,7 +1,7 @@
 """Utility functions to be used across solar positioning algorithms."""
 
-# dt.instant.value = milliseconds since epoch
-@inline fractional_hour(dt::DateTime) = (Dates.value(dt) % 86_400_000) / 3_600_000
+# dt.instant.periods.value = milliseconds since epoch
+@inline fractional_hour(dt::DateTime) = (dt.instant.periods.value % 86_400_000) / 3_600_000
 
 # constants
 const EMR = 6371.01  # Earth Mean Radius in km

src/Refraction/Refraction.jl

@@ -9,7 +9,6 @@ for atmospheric refraction effects.
 module Refraction
 
 using DocStringExtensions: TYPEDFIELDS, TYPEDEF
-using Reexport: @reexport
 
 
 """

src/Refraction/spa.jl

@@ -13,9 +13,9 @@ Solar Position Algorithm (SPA).
 $(TYPEDFIELDS)
 
 # Constructor
-- `SPARefraction()`: Uses default parameters: pressure = 101325 Pa, temperature = 12 °C, refraction_limit = -0.5667°
-- `SPARefraction(pressure, temperature)`: Specify custom pressure [Pa] and temperature [°C], uses default refraction_limit
-- `SPARefraction(pressure, temperature, refraction_limit)`: Also specify refraction limit [degrees]
+- `SPARefraction()`: Uses default parameters: pressure = 101325 Pa, temperature = 12 °C, atmos_refract = -0.5667°
+- `SPARefraction(pressure, temperature)`: Specify custom pressure [Pa] and temperature [°C], uses default atmos_refract
+- `SPARefraction(pressure, temperature, atmos_refract)`: Also specify refraction limit [degrees]
 
 # Notes
 The equation to calculate the refraction correction is given by:
@@ -53,28 +53,25 @@ refraction_correction = refraction(spa, elevation)
 apparent_elevation = elevation + refraction_correction
 ```
 """
-struct SPARefraction{T} <: RefractionAlgorithm where {T<:AbstractFloat}
+@kwdef struct SPARefraction{T<:AbstractFloat} <: RefractionAlgorithm
     "Annual average atmospheric pressure [Pascal]"
-    pressure::T
+    pressure::T = 101325.0
     "Annual average temperature [°C]"
-    temperature::T
+    temperature::T = 12.0
     "Minimum elevation angle for refraction correction [degrees]"
-    refraction_limit::T
+    atmos_refract::T = -0.5667
 end
 
-SPARefraction() = SPARefraction{Float64}(101325.0, 12.0, -0.5667)
+# Positional constructor for 2 arguments (3-argument constructor is already created by @kwdef)
 SPARefraction(pressure::T, temperature::T) where {T<:AbstractFloat} =
     SPARefraction{T}(pressure, temperature, T(-0.5667))
 
 function _refraction(model::SPARefraction{T}, elevation_deg::T) where {T<:AbstractFloat}
     # Convert pressure from Pascal to hPa/mbar
     pressure_hPa = model.pressure / T(100.0)
 
-    # Check if sun is above horizon (elevation >= -0.26667 + refraction_limit)
-    # The sun diameter of 0.26667 degrees is added to the refraction limit
-    above_horizon = elevation_deg >= (T(-0.26667) + model.refraction_limit)
-
-    if !above_horizon
+    # Only apply correction when sun is above horizon accounting for refraction
+    if elevation_deg < model.atmos_refract - T(0.26667)
         return T(0.0)
     end
 

src/Utilities/Utilities.jl

@@ -1,11 +1,9 @@
 module Utilities
 
-using Reexport: @reexport
 using DocStringExtensions: TYPEDEF, TYPEDFIELDS, TYPEDSIGNATURES
 using ..Positioning: Observer, SPA, SolarAlgorithm, calculate_deltat
 import Dates: DateTime, Date, Day
 import TimeZones: ZonedDateTime, timezone, UTC
-using TimeZones: @tz_str
 
 include("spa.jl")
 include("srt.jl")

test/linting.jl

@@ -2,14 +2,14 @@ using SolarPosition
 using Test
 
 using Aqua: Aqua
-using JET: JET
 
 @testset "Aqua tests" begin
     @info "...with Aqua.jl"
     Aqua.test_all(SolarPosition)
 end
 
-if VERSION > v"1.11" # JET v0.10 requires Julia v1.12
+if VERSION == v"1.12" # JET compatibility
+    using JET: JET
     @testset "JET tests" begin
         @info "...with JET.jl"
         JET.test_package(SolarPosition; target_modules = (SolarPosition,))