solar

@library("solar", "0.0.0");

Replace the placeholder with a published version, e.g. greycat install --bump=latest --branch=dev.

Solar position, sunrise and sunset, clear-sky irradiance, transposition, PV output, year-long PV simulations with an embedded NASA POWER cloud climatology, and a reverse model that infers installed power, tilt and azimuth from a production series. Native C, no dependency besides std.

Quick start

var lat = 48.85;
var lng = 2.35;

// Sun position on 2026-06-21 12:00 UTC (about 8 min after solar noon at this longitude)
var t = time::new(1782043200, DurationUnit::seconds);
var pos = SolarPosition::position(t, lat, lng);
println("altitude ${pos.altitude_deg} deg, azimuth ${pos.azimuth_deg} deg");

// Clear-sky irradiance and DC output of a 5 kW array, tilt 35, facing south
var cs = SolarIrradiance::clear_sky_at_point(t, lat, lng, 35.0, 2.0);
var pv = SolarPv::instant_output(5000.0, 0.85, 35.0, 180.0, cs.dni, cs.ghi, cs.dhi, pos.apparent_zenith_deg, pos.azimuth_deg, 0.25);
println("GHI ${cs.ghi} W/m2, power ${pv.power_w} W");

// Sunrise and sunset (UTC) of that day
var st = SolarPosition::sun_times_precise(time::new(1782000000, DurationUnit::seconds), lat, lng);
println("sunrise ${st.sunrise}, sunset ${st.sunset}");

// One year in monthly buckets, clouds from the embedded NASA data
var from = time::new(1767225600, DurationUnit::seconds); // 2026-01-01T00:00:00Z
var to = time::new(1798761600, DurationUnit::seconds); // 2027-01-01T00:00:00Z, exclusive
var months = SolarPv::simulate(lat, lng, 35.0, 5000.0, 35.0, 180.0, from, to, SimulationResolution::monthly, null);
for (_, m in months) {
    println("${m.at}: ${m.wh} Wh (clear sky ${m.wh_max}, std ${m.wh_std})");
}

Conventions

• Units: W, Wh (per bucket), W/m², Wh/m², degrees, Pa, °C, m/s, m.
• Azimuth (sun and panel): clockwise from North, 0 = N, 90 = E, 180 = S. Tilt: 0 = horizontal, 90 = vertical.
• Time: every time is UTC and results are UTC.
• Apparent vs geometric. Irradiance, airmass, clear-sky and tracker functions expect the apparent (refracted) zenith, as pvlib: pass pos.apparent_zenith_deg. altitude and is_daytime are geometric.
• Hour midpoints. simulate*, clear_sky_daily_wh, ideal_daily_wh and infer_system evaluate each hour or interval at its midpoint. to is exclusive; a trailing partial hour is weighted by its duration.
• day_start of sun_times* may be UTC or local midnight: the result is the solar day whose local noon falls in the 24 h after it. clear_sky_daily_wh covers exactly the 24 h from its day_start.
• Argument order. GCL has no named arguments, the order differs between functions, and swapped floats are never detected:

Function	Order
instant_output	capacity_w, efficiency, tilt, azimuth, dni, ghi, dhi, zenith, sun_azimuth, albedo (irradiance before sun angles)
poa_global, poa_global_perez	tilt, azimuth, zenith, sun_azimuth, dni, ghi, dhi, albedo [, dni_extra, airmass]
beam_component	dni, tilt, azimuth, zenith, sun_azimuth
klucher_ / hay_davies_ / perez_diffuse	tilt, azimuth, zenith, sun_azimuth, + dhi, ghi / dhi, dni, dni_extra / dhi, dni, dni_extra, airmass
simulate*, optimal_tilt*, infer_system	lat, lng, altitude_m, ...
clear_sky_daily_wh	capacity_w, efficiency, tilt, azimuth, lat, lng, altitude_m, linke_turbidity, albedo, day_start
horizon_shading	sun altitude (not zenith), sun azimuth, horizon arrays
cell_temperature	poa, ambient_c, wind_ms, a, b, delta_t (SAPM ΔT, not SolarAtmosphere::delta_t)

Function reference

SolarPosition

Function	Notes
position(t, lat, lng)	SunPositionResult, standard refraction (101325 Pa, 15 °C)
position_with_refraction(t, lat, lng, pressure_pa, temp_c)	Custom refraction
altitude(t, lat, lng), is_daytime(t, lat, lng)	Geometric altitude, altitude > 0
sun_times(day_start, lat, lng)	Fast (Spencer): about 1-4 min vs NREL SPA at mid-latitudes, poor within about 1° of the polar thresholds
sun_times_precise(day_start, lat, lng)	Full ephemeris, transit at hour angle 0, rise/set bisected: a few seconds vs NREL SPA, handles polar days
declination_spencer(doy), equation_of_time_spencer(doy), day_of_year(t)	Spencer 1971 (deg, min); 1-366 UTC

SunPositionResult: altitude_deg, apparent_altitude_deg, azimuth_deg, zenith_deg, apparent_zenith_deg, declination_deg, hour_angle_deg (negative in the morning), equation_of_time_min (from the ephemeris).

SunTimesResult: sunrise/sunset are the sun centre at -0.833° (refraction plus semi-diameter) for both functions. transit is always set; sunrise/sunset are null when the event does not happen (polar night, midnight sun); day_length_hours is sunset - sunrise (wrapped by 24 h when sunset precedes sunrise on a polar transition day), else 24 or 0.

SolarAtmosphere

Function	Notes
pressure_from_altitude(altitude_m)	ICAO standard atmosphere (Pa)
relative_airmass(apparent_zenith_deg), relative_airmass_kasten1966(zenith_deg)	Kasten & Young 1989, Kasten 1966; NaN above 90°
absolute_airmass(relative_airmass, pressure_pa)	Pressure-corrected
extraterrestrial_irradiance(doy)	1361 W/m² x Earth-Sun distance factor
delta_t(year)	TT - UT1 (s): Espenak & Meeus 2006 for 1800-2149, -20 + 32u² outside
linke_turbidity_estimate(lat, month)	Rough heuristic, typically 1-3 below climatologies (Remund 2003 / SoDa)
spectral_mismatch(airmass)	SAPM c-Si, absolute airmass

solar_constant (1361) is informational: Ineichen uses 1364 and DISC 1370, as pvlib.

SolarIrradiance

Function	Notes
clear_sky(apparent_zenith_deg, altitude_m, linke_turbidity)	Ineichen/Perez, fixed I0 1364 W/m², TL clamped to [1, 8], zeros at night
clear_sky_at_point(t, lat, lng, altitude_m, linke_turbidity)	Position + Ineichen with seasonal I0
ideal_ghi_at(...)	clear_sky_at_point(...).ghi
ideal_daily_wh(from, to, lat, lng, altitude_m, linke_turbidity)	Clear-sky GHI over [from, to) (Wh/m²)
aoi, aoi_projection, beam_component, isotropic_diffuse, ground_diffuse	Geometry and components
klucher_diffuse, hay_davies_diffuse, perez_diffuse	Klucher 1979, Hay-Davies 1980, Perez 1990 (relative airmass)
poa_global, poa_global_perez	PoaResult: beam + isotropic or Perez + ground
erbs_decomposition(ghi, zenith_deg, doy)	Erbs 1982; zenith > 86.27°: DNI 0, DHI = GHI
disc_decomposition(ghi, zenith_deg, doy, pressure_pa)	Maxwell 1987; zenith >= 87°: DNI 0, DHI = GHI
iam_ashrae(aoi_deg, b)	b typically 0.05
iam_physical(aoi_deg, n, K, L)	De Soto et al. 2006: n 1.526, K 4 1/m, L 0.002 m; only K·L matters
bifacial_ground(ghi, zenith_deg, albedo, gcr, height_m, pitch_m, bifaciality)	Simple empirical estimate ghi·albedo·min(1, (1-gcr)(0.5 + height/pitch))·bifaciality, not a view-factor model
horizon_shading(altitude_deg, azimuth_deg, azimuths, elevations)	1 above the interpolated horizon, else 0; azimuths in any order; same-size arrays

// Clear-sky GHI energy since local midnight
var midnight = time::now().calendar_floor(CalendarUnit::day, TimeZone::"Europe/Paris");
var wh_m2 = SolarIrradiance::ideal_daily_wh(midnight, time::now(), 48.85, 2.35, 35.0, 2.0);
println("${wh_m2} Wh/m2 so far today");

SolarTracking

single_axis(apparent_zenith, apparent_azimuth, axis_tilt, axis_azimuth, max_angle, backtrack, gcr) is an exact port of pvlib.tracking.singleaxis on flat ground, with backtracking (Anderson & Mikofski 2020); gcr = collector width / axis spacing. dual_axis(apparent_zenith, apparent_azimuth, max_zenith) faces the sun, tilt clamped to max_zenith. TrackerResult.tracker_theta is right-handed as in pvlib: positive turns toward axis_azimuth + 90 (west for a north-south axis). At night the panel is flat and aoi is 90.

var t = time::new(1782036000, DurationUnit::seconds); // 2026-06-21T10:00:00Z
var pos = SolarPosition::position(t, 48.85, 2.35);
var tr = SolarTracking::single_axis(pos.apparent_zenith_deg, pos.azimuth_deg, 0.0, 180.0, 60.0, true, 0.35);
var cs = SolarIrradiance::clear_sky_at_point(t, 48.85, 2.35, 35.0, 2.0);
var poa = SolarIrradiance::poa_global(tr.surface_tilt, tr.surface_azimuth, pos.apparent_zenith_deg, pos.azimuth_deg, cs.dni, cs.ghi, cs.dhi, 0.25);
println("theta ${tr.tracker_theta}, POA ${poa.poa_global} W/m2");

SolarPv

Function	Notes
instant_output(...)	poa / 1000 · capacity_w · efficiency, isotropic diffuse, not clipped
orientation_ratio(tilt, azimuth, zenith, sun_azimuth)	Beam factor Rb = cos(AOI)/cos(zenith), 0 behind the panel; not POA/GHI
clear_sky_daily_wh(...)	Equals the daily wh of simulate_clear_sky
cell_temperature(...)	SAPM; open rack glass/glass -3.47, -0.0594, 3
dc_to_ac(dc_w, ac_capacity_w, nominal_efficiency)	PVWatts v5 (pvlib.inverter.pvwatts), clipped at ac_capacity_w
snow_coverage_update(prev, ambient_c, tilt, snowfall_cm_per_hour, poa)	One 1-hour step of Marion et al. 2013 (pvlib.snow.coverage_nrel)
optimal_tilt(lat, lng, altitude_m, linke_turbidity, albedo)	Clear sky, equator-facing
optimal_tilt_real(..., azimuth_deg?)	NASA monthly cloud factors; azimuth defaults to 180 (N) / 0 (S)
optimal_azimuth(lat)	180 for lat >= 0, else 0
inter_row_shading(tilt, azimuth, zenith, sun_azimuth, gcr)	Unshaded row fraction, pvlib.shading.shaded_fraction1d

optimal_tilt* evaluate every whole degree from 0-90° on 12 mid-month days, including both endpoints, to select the global maximum of that objective.

Simulation

SolarPv::simulate(lat, lng, altitude_m, panel_power_w, tilt_deg, azimuth_deg, from, to, resolution, config?)
SolarPv::simulate_clear_sky(<same>)

Per hour midpoint: Ineichen clear sky (seasonal I0), Perez transposition, ASHRAE IAM on beam, SAPM spectral factor, low-irradiance factor below 200 W/m², optional cell temperature (-0.47 %/°C), soiling, efficiency, degradation, optional inverter. SimulationPoint: at, wh, wh_max (clear-sky energy, wh <= wh_max), wh_std.

Clouds. Precedence: measured_ghi > monthly_factors > embedded NASA data (for the local solar day). A cloud factor f is a mixture of clear-sky and overcast time: clear fraction p = clamp((f - 0.3) / 0.7, 0, 1), overcast light diffuse only at min(f, 0.3) x clear-sky GHI, wh = p·wh_max + (1 - p)·P(overcast). f = 1 gives wh_max. measured_ghi is split with Erbs instead. simulate_clear_sky ignores both arrays.

Optimistic by default. NASA factors are ratios to NASA’s own clear sky (real aerosols) but are applied to an Ineichen sky with linke_turbidity 2.0, very clean air. Set linke_turbidity from a climatology when absolute energy matters.

Buckets are UTC. hourly, daily and weekly are 1 h, 1 and 7 days from from (not calendar-aligned); monthly and yearly are calendar months and years. The first bucket is labelled from, the last may be partial: a one-year weekly run has 53 points. to <= from gives an empty array.

wh_std is filled only with the embedded NASA data: wh_max x the day-to-day std of the daily factor, summed within a local solar day and in quadrature across days. Days are treated as independent, so monthly and annual spreads are underestimated (roughly 2x). In the PV convention, P90 (exceeded with 90 % probability) is wh - 1.28·wh_std, P10 is wh + 1.28·wh_std.

SimulationConfig

All fields optional; invalid values raise a runtime error.

Field	Default	Meaning
efficiency	0.85	DC-side losses, excluding inverter and temperature, (0, 1.5]
albedo	0.25	[0, 1]
linke_turbidity	2.0	> 0, clamped to [1, 8]
ambient_temp_c	null	Constant ambient °C; enables temperature derating
wind_speed_ms	1.0	>= 0
ac_capacity_w	null	Inverter AC nameplate: PVWatts part-load curve + clipping, output is AC. Not an export limit (oversized gives near-zero energy). Null: DC energy
soiling_factors	null	12 transmission multipliers in [0, 1] (1 = clean, 0.98 = 2 % loss), UTC month
monthly_factors	null	12 cloud factors in [0, 1], UTC month
measured_ghi, measured_temp_c	null	One finite value per simulated hour (extra ignored)
degradation_rate, system_age_years	0, 0	Output x (1 - rate)^age; rate in [0, 1), age >= 0

var from = time::new(1767225600, DurationUnit::seconds);
var to = time::new(1798761600, DurationUnit::seconds);
var config = SimulationConfig {
    linke_turbidity: 3.5,
    ambient_temp_c: 15.0,
    ac_capacity_w: 4600.0,
    soiling_factors: [0.98, 0.98, 0.98, 0.98, 0.98, 0.98, 0.98, 0.98, 0.98, 0.98, 0.98, 0.98],
};
var year = SolarPv::simulate(48.85, 2.35, 35.0, 5000.0, 35.0, 180.0, from, to, SimulationResolution::yearly, config);
var p = year[0];
println("P50 ${p.wh} Wh, P90 ${p.wh - 1.28 * p.wh_std} Wh (too narrow), clear sky ${p.wh_max} Wh");

Aggregation

aggregate(points, from, to) buckets SimulationPoints by at like simulate: aggregate(simulate(.., hourly, c), hourly, X) equals simulate(.., X, c), including partial hours and polar days. Simulated points retain their interval end and local climate day. Their intervals are split at target bucket boundaries; overlapping stored intervals and duplicate at raise an error. Points are sorted. Manually created points use their stated resolution and the next timestamp to infer their duration. Energy is conserved. A target that is not coarser gives an empty array.

redistribute(values, from, to) does the same for bare Array<float>: values start on 1 January 00:00 UTC and cover at most one year (8784 hourly, 366 daily, 53 weekly or 12 monthly, more is an error); more than 8760 hourly or 365 daily values means a leap year. For multi-year or mid-year data use aggregate.

var from = time::new(1767225600, DurationUnit::seconds);
var to = time::new(1798761600, DurationUnit::seconds);
var daily = SolarPv::simulate(48.85, 2.35, 35.0, 5000.0, 35.0, 180.0, from, to, SimulationResolution::daily, null);
var monthly = SolarPv::aggregate(daily, SimulationResolution::daily, SimulationResolution::monthly);
var values = Array<float> {};
for (_, d in daily) {
    values.add(d.wh);
}
var per_month = SolarPv::redistribute(values, SimulationResolution::daily, SimulationResolution::monthly);
println("${monthly.size()} months, January ${per_month[0]} Wh");

Reverse model

SolarPv::infer_system(lat, lng, altitude_m, start, step, values: Array<float?>, kind, config?): InferredSystem
SolarPv::infer_system_series(lat, lng, altitude_m, samples: Array<Tuple<time, float>>, step?, kind, config?): InferredSystem
SolarPv::align_series(samples, step?, kind, start?): AlignedSeries

Estimates installed power, tilt and azimuth of one fixed plane, plus daily and monthly summaries. values[i] is the sample at start + i * step (UTC); null or non-finite is missing. step must be at most 1 hour (daily totals cannot determine the geometry).

MeasurementKind	Sample	Timestamp
energy_interval_start / _end	Wh over the step	interval start / end
power_average_start / _end	mean W over the step	interval start / end
power_instant	W at t	instant

align_series grids an irregular log. Input may be unsorted; for duplicate timestamps the last finite sample wins. The default step is the native spacing: exact regular cadences are kept (subsecond included), jittered ones are snapped to a nearby whole second or minute. power_instant samples are averaged into the nearest slot; interval kinds are spread over the slots they overlap (energy scaled by step / covered time, mean power time-averaged, slots less than half covered are null), for any step finer or coarser than the native one. The grid starts at start (default: the slot of the first sample). Non-finite values inform the cadence but add no production or coverage. infer_system_series does not return its grid: with weather, call align_series first and pass start, step, values to infer_system.

• Weather. Without it, the fit uses the clear-sky envelope of the data (clouds, curtailment and outages only reduce production). With InferConfig.weather_ghi (optionally weather_dhi, weather_dni, weather_temp_c, weather_wind_ms, same grid as values) every daylight sample is used. Weather that fits better shifted in time, or whose model correlates poorly with production (r < 0.7), is ignored with a warning.
• Detected: clock offset (time_offset_min: local time as UTC, or a wrong kind for steps of about 12.5 min and more), temporarily wrong clock, clipping (ac_limit_w), east/west split, morning/afternoon deficits, flat arrays, night production, low coverage. Thresholds are calibrated on about 140 real plants.
• Power. nameplate_power_w = effective_power_w / performance_ratio (0.85, same meaning as SimulationConfig.efficiency). Noise-free simulate_clear_sky data round-trips within 0.5 % and the exact geometry; simulate output has no clear days, so power and tilt come out low.
• Uncertainties. *_std_* are statistical lower bounds from the fit Jacobian, excluding model, clock and weather errors (often 10-30x smaller than real errors). fit_r2 is against the envelope without weather, the production with it.
• Summaries. daily/monthly UTC ProductionBuckets: measured_wh, coverage, clear_sky_wh, expected_wh (weather, else NASA typical weather with the simulate cloud model), performance.
• Limits. One plane; at least 6 months for a reliable tilt; power is only as good as the unthrottled hours.

InferConfig defaults are the InferConfig::*_default statics: performance_ratio 0.85 (]0, 1.5]), albedo 0.25 ([0, 1]), linke_turbidity 2.0, fit_time_offset true, detect_clipping true, angle_resolution_deg 1.0 (]0, 30], coarser is faster, results are rounded to it). weather_dhi and weather_dni require weather_ghi.

// A year of hourly energy, here the clear-sky output of a known system
var start = time::new(1767225600, DurationUnit::seconds);
var end = time::new(1798761600, DurationUnit::seconds);
var sim = SolarPv::simulate_clear_sky(48.85, 2.35, 35.0, 5000.0, 30.0, 170.0, start, end, SimulationResolution::hourly, null);
var hourly_wh = Array<float?> {};
for (_, p in sim) {
    hourly_wh.add(p.wh);
}
var r = SolarPv::infer_system(48.85, 2.35, 35.0, start, 1hour, hourly_wh, MeasurementKind::energy_interval_start, null);
println("${r.nameplate_power_w} W, tilt ${r.tilt_deg}, azimuth ${r.azimuth_deg}, warnings ${r.warnings}");

Climate data

SolarData reads an embedded NASA POWER climatology: the daily ratio ALLSKY_SFC_SW_DWN / CLRSKY_SFC_SW_DWN (SYN1DEG daily, 2001-2025) on a 1° grid in 52 week bins of 365/52 days, mean and day-to-day std. Interpolation is bilinear in lat/lng and Catmull-Rom in week. Any finite longitude works (0-360 and the antimeridian wrap); lat must be in [-90, 90], day_of_year in 1…366.

Function	Returns
estimate_monthly_factors(lat, lng), estimate_monthly_cloud	12 climatological month means
estimate_weekly_factors(lat, lng), estimate_weekly_cloud	The 52 stored week bins
estimate_factor(lat, lng, doy), estimate_cloud	One day of year
nasa_power_url(lat, lng)	POWER climatology API URL
parse_nasa_power_response(response)	12 factors from the parsed JSON; throws on missing data or -999

Factors are 0-1 (1 = clear sky); CloudFactor.std is the day-to-day std, not the uncertainty of the mean. The data is about 3.7 MB compressed in the library; each tensor (3.4 MB) is decompressed on first use. Regenerate with python3 script/generate_nasa_factors.py (dependencies in its docstring). Regeneration processes spatial tiles to limit memory and uses a validated atomic cache; pass --refresh to download again or --cache path to choose its location.

The live POWER API uses the same 1° grid over its own period. With @library("http", ...) in your project, pass the result as SimulationConfig.monthly_factors:

var response = Http<any>::get(SolarData::nasa_power_url(48.85, 2.35), null);
var factors = SolarData::parse_nasa_power_response(response);
println(factors);

From a body string: SolarData::parse_nasa_power_response(Json<any>::parse(body)).

Acknowledgement. NASA asks users of these data to state: “The data was obtained from National Aeronautics and Space Administration (NASA) Langley Research Center’s Prediction Of Worldwide Energy Resources (POWER) project funded through the NASA Earth Science Division”, with a data reference (POWER SYN1DEG daily, 2001-2025), and requests notification if POWER data is transmitted to other researchers (https://power.larc.nasa.gov/docs/referencing/).

Errors and thread safety

Runtime errors on invalid input:

• simulate*: lat, non-finite arguments, invalid SimulationConfig, time range too large.
• aggregate: duplicate at, overlapping stored intervals, non-finite wh, wh_max or wh_std. redistribute: over a year, non-finite values.
• clear_sky_daily_wh / ideal_daily_wh: time outside the supported range.
• horizon_shading: array sizes differ, non-numbers.
• SolarData: lat, non-finite lng, day of year outside 1…366. parse_nasa_power_response: missing object or month, non-finite or negative value (-999), zero clear-sky value.
• infer_system: step not in ]0, 1 h], lat/lng, altitude_m outside [-500, 9000], empty values, InferConfig out of range, weather length differs from values, weather_dhi/weather_dni without weather_ghi.
• align_series: empty samples, no finite value, no step derivable (a single sample without step), grid too large.

Everything else does not validate and returns NaN, 0 or a clamped value (position with lat 100 returns a result).

All functions are stateless and safe from concurrent tasks; the NASA tensors are published once and read-only afterwards. Result types are @volatile: copy values into your own types to persist them.

References

• Position: Sandia / Grover Hughes (1985) ephemeris as pvlib ephemeris, plus ΔT (Espenak & Meeus 2006), one-term nutation, parallax and Saemundsson refraction; about 0.01° vs NREL SPA (1950-2100). Fast paths: Spencer (1971).
• Clear sky: Ineichen & Perez (2002). Airmass: Kasten & Young (1989), Kasten (1966).
• Transposition: Klucher (1979), Hay & Davies (1980), Perez et al. (1990).
• Decomposition: Erbs et al. (1982), Maxwell (1987).
• IAM: ASHRAE, De Soto et al. (2006). Spectral, thermal: SAPM (King et al. 2004).
• Inverter: PVWatts v5 (Dobos 2014). Tracking: Anderson & Mikofski (2020).
• Snow: Marion et al. (2013). Shading: pvlib shaded_fraction1d.
• Climate: NASA POWER, https://power.larc.nasa.gov/