diff --git a/SAM/SAM.Tests/Resources/Weather/AUS_NSW_Sydney.Intl.AP.947670_TMYx.2004-2018.epw.gz b/SAM/SAM.Tests/Resources/Weather/AUS_NSW_Sydney.Intl.AP.947670_TMYx.2004-2018.epw.gz new file mode 100644 index 000000000..266b0547e Binary files /dev/null and b/SAM/SAM.Tests/Resources/Weather/AUS_NSW_Sydney.Intl.AP.947670_TMYx.2004-2018.epw.gz differ diff --git a/SAM/SAM.Tests/Resources/Weather/HUN_BU_Budapest.Met.Center.128400_TMYx.epw.gz b/SAM/SAM.Tests/Resources/Weather/HUN_BU_Budapest.Met.Center.128400_TMYx.epw.gz new file mode 100644 index 000000000..22b51505e Binary files /dev/null and b/SAM/SAM.Tests/Resources/Weather/HUN_BU_Budapest.Met.Center.128400_TMYx.epw.gz differ diff --git a/SAM/SAM.Tests/Resources/Weather/README.md b/SAM/SAM.Tests/Resources/Weather/README.md new file mode 100644 index 000000000..5b8e6c279 --- /dev/null +++ b/SAM/SAM.Tests/Resources/Weather/README.md @@ -0,0 +1,36 @@ +# Reference EPW fixtures (EPW writer tests) + +Five TMYx EPW files, gzip-compressed (`gzip -9 -n`) without modification. `SHA256SUMS.txt` holds the +SHA-256 of each **uncompressed** original. `WeatherEPWExportTests.ReferenceFixtures_MatchRecordedHashes` +fails if a fixture changes. + +| File | Location | Lat / lon / time zone | Source | Why | +|---|---|---|---|---| +| `HUN_BU_Budapest.Met.Center.128400_TMYx.epw` | Budapest, HU | 47.517 / 19.033 / +1 | SRC-TMYx (ERA5-based) | The F1 incident file: the embedded weather of the SAM_UI PR #239 box model came from it | +| `USA_AK_Fairbanks.Intl.AP.702610_TMYx.2004-2018.epw` | Fairbanks, AK | 64.804 / -147.876 / -9 | ISD-TMYx | High latitude and low sun; 13 degrees west of the zone meridian | +| `AUS_NSW_Sydney.Intl.AP.947670_TMYx.2004-2018.epw` | Sydney, AU | -33.946 / 151.177 / +10 | ISD-TMYx | Southern hemisphere | +| `SGP_SG_Changi.Intl.AP.486980_TMYx.2004-2018.epw` | Singapore | 1.350 / 103.994 / +8 | ISD-TMYx | Tropics; 14 degrees east of the zone meridian | +| `USA_AZ_Phoenix-Sky.Harbor.Intl.AP.722780_TMYx.2004-2018.epw` | Phoenix, AZ | 33.433 / -112.017 / -7 | ISD-TMYx | Clear, high direct normal | + +## Provenance + +- Budapest: copied 2026-10-19 from the SAM_UI PR #239 validation inputs + (`C:\TasOut\licensed-e2e-20261019\inputs-box\Budapest.epw`, machine-local, byte-identical to the + `-EpwPath` file of that validation). Its `COMMENTS 2` header names climate.onebuilding.org as the source. +- The other four: copied 2026-10-19 from the local Ladybug Tools bundle + (`C:\Program Files\ladybug_tools\resources\weather\\.epw`), which redistributes the TMYx + files published at https://climate.onebuilding.org/. + +## Licence + +TMYx files from climate.onebuilding.org are free to use with attribution: "Climatic data provided by +climate.onebuilding.org". + +## Known properties used by the tests + +On the EPW time base (hour ending, local standard time), the Fairbanks, Sydney and Singapore direct +normal closes with its own global - diffuse horizontal within about 1.3 %. The Budapest and Phoenix +files do not (about 9 %). Their direct normal closes with global - diffuse when the solar geometry is +taken about 30 minutes earlier, so within each file the direct normal and the horizontal components +are on time bases about half an hour apart. Budapest's own extraterrestrial column is a further +30 minutes later than standard-time geometry. The tests measure this rather than assume it. diff --git a/SAM/SAM.Tests/Resources/Weather/SGP_SG_Changi.Intl.AP.486980_TMYx.2004-2018.epw.gz b/SAM/SAM.Tests/Resources/Weather/SGP_SG_Changi.Intl.AP.486980_TMYx.2004-2018.epw.gz new file mode 100644 index 000000000..519cd042a Binary files /dev/null and b/SAM/SAM.Tests/Resources/Weather/SGP_SG_Changi.Intl.AP.486980_TMYx.2004-2018.epw.gz differ diff --git a/SAM/SAM.Tests/Resources/Weather/SHA256SUMS.txt b/SAM/SAM.Tests/Resources/Weather/SHA256SUMS.txt new file mode 100644 index 000000000..e4fa46b39 --- /dev/null +++ b/SAM/SAM.Tests/Resources/Weather/SHA256SUMS.txt @@ -0,0 +1,5 @@ +7641aaa476462f76c8598d0836205ccf3df8f8b79aa6343910349303a04d180d *AUS_NSW_Sydney.Intl.AP.947670_TMYx.2004-2018.epw +472ffb4e2f10867600949a392ce46891fc1bd05f2fc8d8917a99bca7d26209c3 *HUN_BU_Budapest.Met.Center.128400_TMYx.epw +d47089f83c4c35f74407050c45e54661f6148c6a35522c9a32af26e4b51144a9 *SGP_SG_Changi.Intl.AP.486980_TMYx.2004-2018.epw +d1591f5ca40e07c229aa781d0514e84722a501678ff276f020698700a5da16b0 *USA_AK_Fairbanks.Intl.AP.702610_TMYx.2004-2018.epw +dec64d7722ed677992f9023160c9c98a1043bb74cbcce696736cbf421a513a80 *USA_AZ_Phoenix-Sky.Harbor.Intl.AP.722780_TMYx.2004-2018.epw diff --git a/SAM/SAM.Tests/Resources/Weather/USA_AK_Fairbanks.Intl.AP.702610_TMYx.2004-2018.epw.gz b/SAM/SAM.Tests/Resources/Weather/USA_AK_Fairbanks.Intl.AP.702610_TMYx.2004-2018.epw.gz new file mode 100644 index 000000000..4921e770c Binary files /dev/null and b/SAM/SAM.Tests/Resources/Weather/USA_AK_Fairbanks.Intl.AP.702610_TMYx.2004-2018.epw.gz differ diff --git a/SAM/SAM.Tests/Resources/Weather/USA_AZ_Phoenix-Sky.Harbor.Intl.AP.722780_TMYx.2004-2018.epw.gz b/SAM/SAM.Tests/Resources/Weather/USA_AZ_Phoenix-Sky.Harbor.Intl.AP.722780_TMYx.2004-2018.epw.gz new file mode 100644 index 000000000..4bbd785b5 Binary files /dev/null and b/SAM/SAM.Tests/Resources/Weather/USA_AZ_Phoenix-Sky.Harbor.Intl.AP.722780_TMYx.2004-2018.epw.gz differ diff --git a/SAM/SAM.Tests/SAM.Tests.csproj b/SAM/SAM.Tests/SAM.Tests.csproj index aaa985fbd..a59f66cc2 100644 --- a/SAM/SAM.Tests/SAM.Tests.csproj +++ b/SAM/SAM.Tests/SAM.Tests.csproj @@ -82,6 +82,10 @@ PreserveNewest + + + PreserveNewest + diff --git a/SAM/SAM.Tests/WeatherEPWExportTests.cs b/SAM/SAM.Tests/WeatherEPWExportTests.cs new file mode 100644 index 000000000..bf85d6090 --- /dev/null +++ b/SAM/SAM.Tests/WeatherEPWExportTests.cs @@ -0,0 +1,794 @@ +// SPDX-License-Identifier: LGPL-3.0-or-later +// Copyright (c) 2020–2026 Michal Dengusiak & Jakub Ziolkowski and contributors + +using SAM.Weather; +using System; +using System.Collections.Generic; +using System.Globalization; +using System.IO; +using System.IO.Compression; +using System.Linq; +using System.Security.Cryptography; +using Xunit; +using Xunit.Abstractions; + +namespace SAM.Tests +{ + /// + /// EPW export of SAM weather (Query.DataString / Convert.ToEPW): radiation field mapping, the + /// horizontal-beam to direct-normal conversion, sky fields and missing values. + /// + /// References: five unmodified TMYx EPWs (Resources/Weather, climate.onebuilding.org) whose + /// direct-normal field is an independent dataset value. SAM imports only global and diffuse + /// horizontal from them, so direct normal has to be reconstructed by the writer and can be + /// compared with the reference. + /// + /// + /// The geometry check uses , the Michalsky (1988) Astronomical Almanac + /// algorithm (sidereal-time route), deliberately not the NOAA/Spencer equation-of-time route the + /// writer uses, so an agreement is not one implementation agreeing with itself. + /// + /// + public class WeatherEPWExportTests + { + private readonly ITestOutputHelper output; + + public WeatherEPWExportTests(ITestOutputHelper output) + { + this.output = output; + } + + public static IEnumerable References() + { + yield return new object[] { "HUN_BU_Budapest.Met.Center.128400_TMYx.epw.gz" }; + yield return new object[] { "USA_AK_Fairbanks.Intl.AP.702610_TMYx.2004-2018.epw.gz" }; + yield return new object[] { "AUS_NSW_Sydney.Intl.AP.947670_TMYx.2004-2018.epw.gz" }; + yield return new object[] { "SGP_SG_Changi.Intl.AP.486980_TMYx.2004-2018.epw.gz" }; + yield return new object[] { "USA_AZ_Phoenix-Sky.Harbor.Intl.AP.722780_TMYx.2004-2018.epw.gz" }; + } + + #region Reference EPWs + + [Fact] + public void ReferenceFixtures_MatchRecordedHashes() + { + string directory = WeatherDirectory(); + Dictionary hashes = File.ReadAllLines(Path.Combine(directory, "SHA256SUMS.txt")) + .Where(x => !string.IsNullOrWhiteSpace(x)) + .Select(x => x.Split(new[] { ' ' }, 2)) + .ToDictionary(x => x[1].TrimStart('*').Trim(), x => x[0]); + + foreach (object[] reference in References()) + { + string name = (string)reference[0]; + using SHA256 sHA256 = SHA256.Create(); + string hash = string.Concat(sHA256.ComputeHash(Decompress(Path.Combine(directory, name))).Select(x => x.ToString("x2"))); + + // The hash is of the original, uncompressed EPW. + Assert.Equal(hashes[name.Substring(0, name.Length - ".gz".Length)], hash); + } + } + + [Theory] + [MemberData(nameof(References))] + public void ReferenceEpw_Export_PreservesHorizontalAndSkyFields(string name) + { + List lines = ReadReference(name).Split('\n').Select(x => x.TrimEnd('\r')).ToList(); + List reference = Rows(lines); + List exported = Rows(Weather.Convert.ToEPW(Weather.Convert.ToSAM(lines)).Split('\n')); + + Assert.Equal(8760, reference.Count); + Assert.Equal(8760, exported.Count); + + for (int i = 0; i < reference.Count; i++) + { + EpwRow r = reference[i]; + EpwRow e = exported[i]; + + Assert.Equal((r.Month, r.Day, r.Hour), (e.Month, e.Day, e.Hour)); + Assert.Equal(32, e.FieldCount); + + // Stated values are written through unchanged. + Assert.Equal(r.DryBulb, e.DryBulb, 6); + Assert.Equal(r.DewPoint, e.DewPoint, 6); + Assert.Equal(r.GlobalHorizontal, e.GlobalHorizontal, 6); + Assert.Equal(r.DiffuseHorizontal, e.DiffuseHorizontal, 6); + Assert.Equal(r.TotalSkyCover, e.TotalSkyCover, 6); + + // SAM holds one cloud amount: it is written to opaque sky cover too. + Assert.Equal(r.TotalSkyCover, e.OpaqueSkyCover, 6); + + // Not held by SAM: missing, so EnergyPlus derives it (never a 0 W/m2 sky). + Assert.Equal(9999, e.HorizontalInfrared); + + Assert.InRange(e.DirectNormal, 0, 1500); + Assert.InRange(e.ExtraterrestrialHorizontal, 0, 1420); + } + } + + /// + /// References whose own direct normal closes with their own global - diffuse on the EPW + /// time base (hour ending, local standard time): reconstructing direct normal from the + /// horizontal components must reproduce the dataset value - annually, hourly with the sun up, + /// and as beam on the horizontal and on vertical facades facing N, E, S and W. The statistics + /// are written to the test output (they are the evidence quoted in the PR record). + /// + [Theory] + [InlineData("USA_AK_Fairbanks.Intl.AP.702610_TMYx.2004-2018.epw.gz")] + [InlineData("AUS_NSW_Sydney.Intl.AP.947670_TMYx.2004-2018.epw.gz")] + [InlineData("SGP_SG_Changi.Intl.AP.486980_TMYx.2004-2018.epw.gz")] + public void ConsistentReferenceEpw_Export_ReproducesDirectNormalAndFacadeBeam(string name) + { + Comparison comparison = Compare(name); + + Assert.InRange(comparison.ReferenceClosureError, 0, 0.02); + + Assert.InRange(comparison.ExportedBeamOnHorizontal / comparison.HorizontalBeam, 0.995, 1.005); + Assert.InRange(comparison.ExportedDirectNormal / comparison.ReferenceDirectNormal, 0.99, 1.01); + Assert.InRange(System.Math.Abs(comparison.MedianRelativeError), 0, 0.005); + Assert.InRange(comparison.Percentile90AbsoluteRelativeError, 0, 0.02); + + for (int k = 0; k < 4; k++) + { + // Facades that receive almost no beam (south in Sydney) are compared absolutely. + Assert.True(System.Math.Abs(comparison.ExportedFacadeBeam[k] - comparison.ReferenceFacadeBeam[k]) <= System.Math.Max(0.015 * comparison.ReferenceFacadeBeam[k], 1000), $"facade {AlmanacSun.FacadeNames[k]}"); + } + } + + /// + /// References whose own direct normal does NOT close with their own global - diffuse on the + /// EPW time base: the closure is best with the geometry about 30 minutes earlier, i.e. the + /// dataset's direct normal and its horizontal components are on time bases half an hour apart + /// (and Budapest's own extraterrestrial column is 30 minutes later still). No conversion that + /// follows the EPW convention can reproduce such a direct normal from global - diffuse. What + /// is required is that the export keeps the stated horizontal beam (energy closure) and that + /// the inconsistency is the dataset's, which this test measures rather than assumes. + /// + [Theory] + [InlineData("HUN_BU_Budapest.Met.Center.128400_TMYx.epw.gz")] + [InlineData("USA_AZ_Phoenix-Sky.Harbor.Intl.AP.722780_TMYx.2004-2018.epw.gz")] + public void InconsistentReferenceEpw_Export_KeepsHorizontalBeam(string name) + { + Comparison comparison = Compare(name); + + Assert.True(comparison.ReferenceClosureError > 0.05, "the reference is internally inconsistent on the EPW time base"); + Assert.True(comparison.ReferenceClosureErrorShifted < comparison.ReferenceClosureError / 2, "and consistent with the geometry 30 minutes earlier"); + + // All of the stated horizontal beam is placed except what exceeds the Rayleigh bound at grazing sun. + Assert.InRange(comparison.ExportedBeamOnHorizontal / comparison.HorizontalBeam, 0.99, 1.0005); + + // Where the sun is high the half-hour inconsistency is small. + Assert.InRange(System.Math.Abs(comparison.MedianRelativeError), 0, 0.03); + } + + /// + /// How often the Rayleigh-atmosphere bound on direct normal binds, and what it removes, against + /// the same conversion without the bound. The writer's own hour geometry is read back from the + /// export (extraterrestrial horizontal / extraterrestrial normal = hour-mean cos zenith), so + /// "bound" means exactly the hours the bound changed (more than 1 %, above the 2-decimal + /// rounding of the fields). Beam stated in hours where the sun is below the horizon all hour is + /// a separate rule (direct normal 0) and is counted separately. Facade beam uses the independent + /// incidence geometry for both variants. The statistics are written to the test output (PR + /// record evidence). + /// + [Theory] + [MemberData(nameof(References))] + public void RayleighBound_BindsOnlyAtGrazingSun_AndRemovesLittleEnergy(string name) + { + List lines = ReadReference(name).Split('\n').Select(x => x.TrimEnd('\r')).ToList(); + List reference = Rows(lines); + WeatherData weatherData = Weather.Convert.ToSAM(lines); + List exported = Rows(Weather.Convert.ToEPW(weatherData).Split('\n')); + double timeZone = Core.Query.Double(Core.Query.UTC(weatherData.GetValue(WeatherDataParameter.TimeZone))); + + int beamHours = 0, boundHours = 0, sunDownBeamHours = 0, momentaryHours = 0; + double sunDownBeam = 0, maximumCosZenithBound = 0, maximumUnbounded = 0; + double directNormalBounded = 0, directNormalUnbounded = 0; + double beamBounded = 0, beamUnbounded = 0; + double[] facadeBounded = new double[4], facadeUnbounded = new double[4]; + for (int i = 0; i < reference.Count; i++) + { + EpwRow r = reference[i]; + EpwRow e = exported[i]; + + double beam = System.Math.Max(0, r.GlobalHorizontal - r.DiffuseHorizontal); + if (beam <= 0) + continue; + + beamHours++; + if (e.ExtraterrestrialDirectNormal > 0 && e.ExtraterrestrialHorizontal <= 0) + { + // Sun up for about a minute: the extraterrestrial horizontal rounds to 0.00, so the + // geometry cannot be read back; the beam is far beyond any atmosphere and is bounded. + momentaryHours++; + boundHours++; + continue; + } + + double meanCosZenith = e.ExtraterrestrialDirectNormal > 0 ? e.ExtraterrestrialHorizontal / e.ExtraterrestrialDirectNormal : 0; + if (meanCosZenith <= 0) + { + sunDownBeamHours++; + sunDownBeam += beam; + Assert.Equal(0, e.DirectNormal); + continue; + } + + double unbounded = beam / meanCosZenith; + if (e.DirectNormal < 0.99 * unbounded) + { + boundHours++; + maximumCosZenithBound = System.Math.Max(maximumCosZenithBound, meanCosZenith); + maximumUnbounded = System.Math.Max(maximumUnbounded, unbounded); + } + + AlmanacSun.HourGeometry geometry = AlmanacSun.Hour(2018, r.Month, r.Day, r.Hour, weatherData.Latitude, weatherData.Longitude, timeZone, 0, 120); + directNormalBounded += e.DirectNormal; + directNormalUnbounded += unbounded; + beamBounded += e.DirectNormal * meanCosZenith; + beamUnbounded += beam; + for (int k = 0; k < 4; k++) + { + facadeBounded[k] += e.DirectNormal * geometry.MeanCosIncidence[k]; + facadeUnbounded[k] += unbounded * geometry.MeanCosIncidence[k]; + } + } + + static string Change(double bounded, double unbounded) => unbounded > 0 ? $"{(bounded / unbounded - 1) * 100:+0.00;-0.00;0.00} %" : "n/a"; + + output.WriteLine($"{name}: bound binds in {boundHours} of {beamHours} beam hours ({100.0 * boundHours / beamHours:F2} %), all with hour-mean cos zenith <= {maximumCosZenithBound:F3} (sun below ~{System.Math.Asin(maximumCosZenithBound) * 180 / System.Math.PI:F1} deg); largest unbounded DNI {maximumUnbounded:F0} W/m2"); + output.WriteLine($" beam with the sun below the horizon all hour (DNI 0 by rule): {sunDownBeamHours} hours, {sunDownBeam / 1000:F2} kWh/m2"); + output.WriteLine($" of the bound hours, {momentaryHours} have the sun up for only about a minute (geometry not read back; excluded from the energy sums)"); + output.WriteLine($" annual DNI: unbounded {directNormalUnbounded / 1000:F1}, bounded {directNormalBounded / 1000:F1} kWh/m2 ({Change(directNormalBounded, directNormalUnbounded)})"); + output.WriteLine($" annual horizontal beam: stated {beamUnbounded / 1000:F1}, exported {beamBounded / 1000:F1} kWh/m2 ({Change(beamBounded, beamUnbounded)})"); + for (int k = 0; k < 4; k++) + output.WriteLine($" facade {AlmanacSun.FacadeNames[k]} annual beam: unbounded {facadeUnbounded[k] / 1000:F1}, bounded {facadeBounded[k] / 1000:F1} kWh/m2 ({Change(facadeBounded[k], facadeUnbounded[k])})"); + + // Never for a self-consistent dataset; otherwise rare, only with the sun low, and a + // negligible share of the beam on the horizontal. + if (!name.StartsWith("HUN_") && !name.StartsWith("USA_AZ_")) + Assert.Equal(0, boundHours); + + Assert.InRange(boundHours, 0, beamHours / 10); + Assert.InRange(maximumCosZenithBound, 0, 0.4); + Assert.InRange(beamBounded / beamUnbounded, 0.99, 1.0005); + } + + /// + /// The F1 incident hours (SAM_UI PR #239): the embedded-weather export wrote global horizontal 0 + /// and the horizontal beam (global - diffuse) into the direct-normal field. + /// + [Theory] + [InlineData(6, 16, 8, 403, 569, 120)] + [InlineData(6, 16, 10, 623, 600, 174)] + public void BudapestIncidentHours_ExportGlobalAndDirectNormal(int month, int day, int hour, double globalHorizontal, double directNormal, double diffuseHorizontal) + { + List lines = ReadReference("HUN_BU_Budapest.Met.Center.128400_TMYx.epw.gz").Split('\n').Select(x => x.TrimEnd('\r')).ToList(); + EpwRow e = Rows(Weather.Convert.ToEPW(Weather.Convert.ToSAM(lines)).Split('\n')).Single(x => x.Month == month && x.Day == day && x.Hour == hour); + + Assert.Equal(globalHorizontal, e.GlobalHorizontal, 6); + Assert.Equal(diffuseHorizontal, e.DiffuseHorizontal, 6); + + // The direct normal puts exactly the stated beam on the horizontal... + double meanCosZenith = AlmanacSun.Hour(2018, month, day, hour, 47.517, 19.033, 1).MeanCosZenith; + Assert.InRange(e.DirectNormal * meanCosZenith, (globalHorizontal - diffuseHorizontal) * 0.995, (globalHorizontal - diffuseHorizontal) * 1.005); + + // ...is not the horizontal beam any more, and is near the dataset's own value, which sits + // on a time base half an hour earlier (see InconsistentReferenceEpw_Export_KeepsHorizontalBeam). + Assert.True(e.DirectNormal - (globalHorizontal - diffuseHorizontal) > 100); + Assert.InRange(e.DirectNormal, directNormal * 0.8, directNormal * 1.05); + } + + [Theory] + [MemberData(nameof(References))] + public void ReferenceEpw_ExportReimportExport_IsStable(string name) + { + List lines = ReadReference(name).Split('\n').Select(x => x.TrimEnd('\r')).ToList(); + string first = Weather.Convert.ToEPW(Weather.Convert.ToSAM(lines)); + string second = Weather.Convert.ToEPW(Weather.Convert.ToSAM(first.Split('\n').ToList())); + + Assert.Equal(first, second); + } + + [Fact] + public void ReferenceEpw_Import_IsUnchanged() + { + // The explicit-EPW route reads the file itself; the SAM reader keeps its mapping. + List lines = ReadReference("HUN_BU_Budapest.Met.Center.128400_TMYx.epw.gz").Split('\n').Select(x => x.TrimEnd('\r')).ToList(); + WeatherData weatherData = Weather.Convert.ToSAM(lines); + WeatherHour weatherHour = weatherData.GetWeatherHour(new DateTime(2018, 6, 16, 9, 0, 0)); + + Assert.Equal(623, weatherHour[WeatherDataType.GlobalSolarRadiation]); + Assert.Equal(174, weatherHour[WeatherDataType.DiffuseSolarRadiation]); + Assert.True(double.IsNaN(weatherHour[WeatherDataType.DirectSolarRadiation])); + Assert.Equal(17.3, weatherHour[WeatherDataType.DryBulbTemperature]); + Assert.Equal(5.9, weatherHour[WeatherDataType.WetBulbTemperature]); + Assert.Equal(0.2, weatherHour[WeatherDataType.CloudCover], 12); + Assert.Equal("UTC+01:00", weatherData.GetValue(WeatherDataParameter.TimeZone)); + } + + #endregion + + #region Synthetic hours + + /// + /// A known direct normal at solar noon, low sun, winter and summer, both hemispheres and + /// longitudes far from the zone meridian: the beam put on the horizontal with the independent + /// geometry comes back as that direct normal. Within 1 % with the sun up; within 3 % in hours + /// where the sun stays below ~5 degrees, where the two solar-position algorithms' few + /// arc-minutes of difference are a percent of the altitude. + /// + [Theory] + [InlineData(0.0, 0.0, 0.0, 3, 20, 11, 900.0, 0.01)] // equator, equinox, hour ending 12:00 (solar noon at 12:07) + [InlineData(51.5, -0.1, 0.0, 6, 21, 12, 850.0, 0.01)] // London, summer, hour 12-13 + [InlineData(51.5, -0.1, 0.0, 12, 21, 11, 600.0, 0.01)] // London, winter noon, altitude ~15 deg + [InlineData(51.5, -0.1, 0.0, 12, 21, 8, 150.0, 0.03)] // London, winter, first hour after sunrise + [InlineData(-33.9, 151.2, 10.0, 12, 21, 12, 900.0, 0.01)] // Sydney, southern summer + [InlineData(-33.9, 151.2, 10.0, 6, 21, 7, 200.0, 0.03)] // Sydney, southern winter, sunrise hour + [InlineData(64.8, -147.9, -9.0, 12, 21, 13, 300.0, 0.03)] // Fairbanks, ~2 deg sun, 13 deg west of the zone meridian + [InlineData(1.35, 104.0, 8.0, 9, 1, 13, 800.0, 0.01)] // Singapore, 14 deg east of the zone meridian + [InlineData(28.6, 77.1, 5.5, 4, 15, 10, 850.0, 0.01)] // New Delhi, half-hour time zone + public void KnownDirectNormal_PutOnHorizontal_IsRecovered(double latitude, double longitude, double timeZone, int month, int day, int hour, double directNormal, double tolerance) + { + double meanCosZenith = AlmanacSun.MeanCosZenith(2025, month, day, hour + 1, latitude, longitude, timeZone); + Assert.True(meanCosZenith > 0.01, "the case must have the sun up"); + + double diffuse = 80; + EpwRow e = SingleHour(latitude, longitude, timeZone, month, day, hour, global: directNormal * meanCosZenith + diffuse, diffuse: diffuse); + + output.WriteLine($"mean cos zenith {meanCosZenith:F4}, exported DNI {e.DirectNormal}"); + Assert.InRange(e.DirectNormal, directNormal * (1 - tolerance), directNormal * (1 + tolerance)); + Assert.Equal(System.Math.Round(directNormal * meanCosZenith + diffuse, 2), e.GlobalHorizontal, 6); + Assert.Equal(diffuse, e.DiffuseHorizontal, 6); + + // Extraterrestrial horizontal is the solar constant (eccentricity corrected) on the same geometry. + Assert.InRange(e.ExtraterrestrialHorizontal / (e.ExtraterrestrialDirectNormal * meanCosZenith), 1 - tolerance, 1 + tolerance); + Assert.InRange(e.ExtraterrestrialDirectNormal, 1320, 1415); + } + + [Fact] + public void SolarNoon_Equinox_Equator_GivesNearOverheadGeometry() + { + // Hours ending 12:00 and 13:00 straddle solar noon (12:07 at lon 0 on 20 March). + EpwRow before = SingleHour(0, 0, 0, 3, 20, 11, global: 900, diffuse: 100); + EpwRow after = SingleHour(0, 0, 0, 3, 20, 12, global: 900, diffuse: 100); + + // Sun within ~16 deg of the zenith all the time: DNI only slightly above the horizontal beam. + Assert.InRange(before.DirectNormal, 800, 830); + Assert.InRange(after.DirectNormal, 800, 830); + Assert.InRange(before.ExtraterrestrialHorizontal, 1340, 1370); + } + + [Fact] + public void Night_WritesZeroRadiation_AndNoExtraterrestrial() + { + EpwRow e = SingleHour(51.5, -0.1, 0, 1, 15, 1, global: 0, diffuse: 0); + + Assert.Equal(0, e.GlobalHorizontal); + Assert.Equal(0, e.DirectNormal); + Assert.Equal(0, e.DiffuseHorizontal); + Assert.Equal(0, e.ExtraterrestrialHorizontal); + Assert.Equal(0, e.ExtraterrestrialDirectNormal); + } + + [Fact] + public void BeamWithSunBelowHorizon_IsNotInventedAsDirectNormal_AndGlobalIsKept() + { + // 02:00-03:00 in January in London: a horizontal beam here is a timing or data error; + // it cannot be a direct-normal value. The stated global and diffuse stay as stated. + EpwRow e = SingleHour(51.5, -0.1, 0, 1, 15, 2, global: 30, diffuse: 10); + + Assert.Equal(0, e.DirectNormal); + Assert.Equal(30, e.GlobalHorizontal); + Assert.Equal(10, e.DiffuseHorizontal); + } + + [Fact] + public void SunriseHour_DirectNormalIsLimitedToRayleighAtmosphere() + { + // London, 21 December, 08:00-09:00: the sun rises at ~08:04 and stays below ~4 deg. + // 180 W/m2 of horizontal beam would need ~4500 W/m2 direct normal - impossible. The + // writer keeps the stated global and diffuse and bounds the direct normal by what a pure + // Rayleigh atmosphere transmits at that air mass (here roughly half the extraterrestrial). + double meanCosZenith = AlmanacSun.MeanCosZenith(2025, 12, 21, 9, 51.5, -0.1, 0); + Assert.InRange(meanCosZenith, 0.02, 0.06); + + EpwRow e = SingleHour(51.5, -0.1, 0, 12, 21, 8, global: 200, diffuse: 20); + + Assert.InRange(e.DirectNormal, 0.3 * e.ExtraterrestrialDirectNormal, 0.7 * e.ExtraterrestrialDirectNormal); + Assert.True(e.DirectNormal * meanCosZenith < 180); + Assert.Equal(200, e.GlobalHorizontal); + Assert.Equal(20, e.DiffuseHorizontal); + + // A physically possible beam at the same hour is not touched by the bound. + EpwRow possible = SingleHour(51.5, -0.1, 0, 12, 21, 8, global: 20 + 100 * meanCosZenith, diffuse: 20); + Assert.InRange(possible.DirectNormal, 97, 103); + } + + [Fact] + public void HighSun_RayleighBoundIsNotReached_ByARealClearSky() + { + // Phoenix-like clear noon, 1050 W/m2 direct normal (above typical sea-level clear-sky maxima). + double meanCosZenith = AlmanacSun.MeanCosZenith(2025, 6, 21, 13, 33.4, -112.0, -7); + EpwRow e = SingleHour(33.4, -112.0, -7, 6, 21, 12, global: 1050 * meanCosZenith + 90, diffuse: 90); + + Assert.InRange(e.DirectNormal, 1040, 1060); + } + + [Fact] + public void MissingLocation_WritesDirectNormalMissing_NotTheHorizontalBeam() + { + WeatherDay weatherDay = Day(global: 500, diffuse: 200); + WeatherYear weatherYear = new WeatherYear(2025); + weatherYear[170] = weatherDay; + + // No time zone stated: the solar geometry is undefined. + WeatherData weatherData = new WeatherData("Station", "Description", 51.5, -0.1, 10, weatherYear); + EpwRow e = Rows(weatherData.DataString().Split('\n'))[12]; + + Assert.Equal(9999, e.DirectNormal); + Assert.Equal(9999, e.ExtraterrestrialHorizontal); + Assert.Equal(500, e.GlobalHorizontal); + Assert.Equal(200, e.DiffuseHorizontal); + + // The location-less public overloads behave the same. + EpwRow overload = Rows(Query.DataString(2025, 6, 20, 12, weatherDay).Split('\n'))[0]; + Assert.Equal(9999, overload.DirectNormal); + Assert.Equal(500, overload.GlobalHorizontal); + } + + [Fact] + public void MissingRadiationAndSky_AreWrittenAsEpwMissingValues_NeverZero() + { + WeatherDay weatherDay = new WeatherDay(); + weatherDay[WeatherDataType.DryBulbTemperature] = Enumerable.Repeat(10.0, 24).ToArray(); + + EpwRow e = SingleHour(weatherDay, 51.5, -0.1, 0, 6, 20, 12); + + Assert.Equal(9999, e.GlobalHorizontal); + Assert.Equal(9999, e.DirectNormal); + Assert.Equal(9999, e.DiffuseHorizontal); + Assert.Equal(9999, e.HorizontalInfrared); + Assert.Equal(99, e.TotalSkyCover); + Assert.Equal(99, e.OpaqueSkyCover); + Assert.Equal(99.9, e.DewPoint); + } + + [Fact] + public void PartialRadiation_IsCompletedFromGlobalEqualsDirectPlusDiffuse() + { + // Direct (horizontal) and diffuse stated, global not: global is their sum. + WeatherDay weatherDay = new WeatherDay(); + weatherDay[WeatherDataType.DirectSolarRadiation] = Enumerable.Repeat(300.0, 24).ToArray(); + weatherDay[WeatherDataType.DiffuseSolarRadiation] = Enumerable.Repeat(100.0, 24).ToArray(); + + EpwRow e = SingleHour(weatherDay, 51.5, -0.1, 0, 6, 21, 12); + double meanCosZenith = AlmanacSun.MeanCosZenith(2025, 6, 21, 13, 51.5, -0.1, 0); + + Assert.Equal(400, e.GlobalHorizontal); + Assert.Equal(100, e.DiffuseHorizontal); + Assert.InRange(e.DirectNormal, 300 / meanCosZenith * 0.99, 300 / meanCosZenith * 1.01); + } + + [Fact] + public void TasStyleWeather_DerivesDewPoint_AndWritesCloudInTenths() + { + // TAS weather has no dew point and no pressure; cloud cover is 0-1. + WeatherDay weatherDay = new WeatherDay(); + weatherDay[WeatherDataType.DryBulbTemperature] = Enumerable.Repeat(20.0, 24).ToArray(); + weatherDay[WeatherDataType.RelativeHumidity] = Enumerable.Repeat(50.0, 24).ToArray(); + weatherDay[WeatherDataType.CloudCover] = Enumerable.Repeat(0.7, 24).ToArray(); + + EpwRow e = SingleHour(weatherDay, 51.5, -0.1, 0, 6, 21, 12); + + // 20 C, 50 %: dew point 9.26 C (psychrometric tables 9.3 C). + Assert.InRange(e.DewPoint, 9.2, 9.35); + Assert.Equal(7, e.TotalSkyCover); + Assert.Equal(7, e.OpaqueSkyCover); + } + + [Fact] + public void LeapYear_WritesOneRecordPerHourOf365Days_WithGeometryOfTheWrittenDate() + { + WeatherYear leap = new WeatherYear(2024); + WeatherYear common = new WeatherYear(2023); + for (int i = 0; i < 365; i++) + { + leap[i] = Day(global: 500, diffuse: 200); + common[i] = Day(global: 500, diffuse: 200); + } + + List leapRows = Rows(Located(leap).DataString().Split('\n')); + List commonRows = Rows(Located(common).DataString().Split('\n')); + + Assert.Equal(8760, leapRows.Count); + Assert.Equal(8760, commonRows.Count); + + // SAM WeatherYear has 365 days: in a leap year day index 59 is written as 29 February + // (existing labelling, unchanged here). The geometry follows the written date. + Assert.Contains(leapRows, x => x.Month == 2 && x.Day == 29); + EpwRow leapJune = leapRows.Single(x => x.Month == 6 && x.Day == 21 && x.Hour == 13); + EpwRow commonJune = commonRows.Single(x => x.Month == 6 && x.Day == 21 && x.Hour == 13); + Assert.InRange(leapJune.DirectNormal / commonJune.DirectNormal, 0.995, 1.005); + } + + #endregion + + #region Helpers + + private static WeatherDay Day(double global, double diffuse) + { + WeatherDay weatherDay = new WeatherDay(); + weatherDay[WeatherDataType.GlobalSolarRadiation] = Enumerable.Repeat(global, 24).ToArray(); + weatherDay[WeatherDataType.DiffuseSolarRadiation] = Enumerable.Repeat(diffuse, 24).ToArray(); + weatherDay[WeatherDataType.DryBulbTemperature] = Enumerable.Repeat(15.0, 24).ToArray(); + return weatherDay; + } + + private static WeatherData Located(WeatherYear weatherYear) + { + WeatherData weatherData = new WeatherData("Station", "Description", 51.5, -0.1, 10, weatherYear); + weatherData.SetValue(WeatherDataParameter.TimeZone, Core.Query.Description(Core.Query.UTC(0.0))); + return weatherData; + } + + private static EpwRow SingleHour(double latitude, double longitude, double timeZone, int month, int day, int hour, double global, double diffuse) + { + return SingleHour(Day(global, diffuse), latitude, longitude, timeZone, month, day, hour); + } + + /// + /// Exports one hour through the public WeatherData route (the route Convert.ToEPW uses). + /// + private static EpwRow SingleHour(WeatherDay weatherDay, double latitude, double longitude, double timeZone, int month, int day, int hour) + { + WeatherYear weatherYear = new WeatherYear(2025); + int index = new DateTime(2025, month, day).DayOfYear - 1; + weatherYear[index] = weatherDay; + + WeatherData weatherData = new WeatherData("Station", "Description", latitude, longitude, 10, weatherYear); + weatherData.SetValue(WeatherDataParameter.TimeZone, Core.Query.Description(Core.Query.UTC(timeZone))); + + return Rows(weatherData.DataString().Split('\n')).Single(x => x.Month == month && x.Day == day && x.Hour == hour + 1); + } + + private static string WeatherDirectory() + { + return Path.Combine(AppContext.BaseDirectory, "Resources", "Weather"); + } + + private static byte[] Decompress(string path) + { + using FileStream fileStream = File.OpenRead(path); + using GZipStream gZipStream = new GZipStream(fileStream, CompressionMode.Decompress); + using MemoryStream memoryStream = new MemoryStream(); + gZipStream.CopyTo(memoryStream); + return memoryStream.ToArray(); + } + + private static string ReadReference(string name) + { + return System.Text.Encoding.UTF8.GetString(Decompress(Path.Combine(WeatherDirectory(), name))); + } + + private static List Rows(IEnumerable lines) + { + List result = new List(); + foreach (string line in lines) + { + string[] values = line.Trim().Split(','); + if (values.Length < 24 || !int.TryParse(values[0], out _) || !int.TryParse(values[1], out _)) + continue; + + double Value(int index) => double.Parse(values[index], CultureInfo.InvariantCulture); + + result.Add(new EpwRow + { + FieldCount = values.Length, + Month = int.Parse(values[1]), + Day = int.Parse(values[2]), + Hour = int.Parse(values[3]), + DryBulb = Value(6), + DewPoint = Value(7), + ExtraterrestrialHorizontal = Value(10), + ExtraterrestrialDirectNormal = Value(11), + HorizontalInfrared = Value(12), + GlobalHorizontal = Value(13), + DirectNormal = Value(14), + DiffuseHorizontal = Value(15), + TotalSkyCover = Value(22), + OpaqueSkyCover = Value(23), + }); + } + + return result; + } + + private sealed class Comparison + { + public double ReferenceDirectNormal; + public double ExportedDirectNormal; + public double HorizontalBeam; + public double ExportedBeamOnHorizontal; + public double ReferenceClosureError; + public double ReferenceClosureErrorShifted; + public double MedianRelativeError; + public double Percentile90AbsoluteRelativeError; + public double[] ReferenceFacadeBeam = new double[4]; + public double[] ExportedFacadeBeam = new double[4]; + } + + /// + /// Imports a reference EPW into SAM, exports it again and compares the exported direct normal + /// with the reference's own, using the independent geometry. Closure error: sum over hours of + /// |DNI_ref * meanCosZenith - (GHI - DHI)| / sum (GHI - DHI), on the EPW time base and with the + /// geometry 30 minutes earlier. + /// + private Comparison Compare(string name) + { + List lines = ReadReference(name).Split('\n').Select(x => x.TrimEnd('\r')).ToList(); + List reference = Rows(lines); + WeatherData weatherData = Weather.Convert.ToSAM(lines); + List exported = Rows(Weather.Convert.ToEPW(weatherData).Split('\n')); + Assert.Equal(reference.Count, exported.Count); + + double timeZone = Core.Query.Double(Core.Query.UTC(weatherData.GetValue(WeatherDataParameter.TimeZone))); + + Comparison result = new Comparison(); + double closure = 0; + double closureShifted = 0; + List relativeErrors = new List(); + for (int i = 0; i < reference.Count; i++) + { + EpwRow r = reference[i]; + EpwRow e = exported[i]; + + AlmanacSun.HourGeometry geometry = AlmanacSun.Hour(2018, r.Month, r.Day, r.Hour, weatherData.Latitude, weatherData.Longitude, timeZone, 0, 120); + AlmanacSun.HourGeometry shifted = AlmanacSun.Hour(2018, r.Month, r.Day, r.Hour, weatherData.Latitude, weatherData.Longitude, timeZone, -30, 60); + double beam = System.Math.Max(0, r.GlobalHorizontal - r.DiffuseHorizontal); + + result.ReferenceDirectNormal += r.DirectNormal; + result.ExportedDirectNormal += e.DirectNormal; + result.HorizontalBeam += beam; + result.ExportedBeamOnHorizontal += e.DirectNormal * geometry.MeanCosZenith; + closure += System.Math.Abs(r.DirectNormal * geometry.MeanCosZenith - beam); + closureShifted += System.Math.Abs(r.DirectNormal * shifted.MeanCosZenith - beam); + + for (int k = 0; k < 4; k++) + { + result.ReferenceFacadeBeam[k] += r.DirectNormal * geometry.MeanCosIncidence[k]; + result.ExportedFacadeBeam[k] += e.DirectNormal * geometry.MeanCosIncidence[k]; + } + + // Sun above ~15 degrees for the whole hour and a material beam: the conversion is well conditioned. + if (geometry.MeanCosZenith > 0.26 && r.DirectNormal > 100) + relativeErrors.Add((e.DirectNormal - r.DirectNormal) / r.DirectNormal); + } + + relativeErrors.Sort(); + result.ReferenceClosureError = closure / result.HorizontalBeam; + result.ReferenceClosureErrorShifted = closureShifted / result.HorizontalBeam; + result.MedianRelativeError = relativeErrors[relativeErrors.Count / 2]; + result.Percentile90AbsoluteRelativeError = relativeErrors.Select(System.Math.Abs).OrderBy(x => x).ElementAt((int)(relativeErrors.Count * 0.9)); + + output.WriteLine($"{name}: lat {weatherData.Latitude}, lon {weatherData.Longitude}, tz {timeZone}"); + output.WriteLine($" reference closure |DNI*cosZ - (GHI-DHI)|: {result.ReferenceClosureError * 100:F2} % (geometry 30 min earlier: {result.ReferenceClosureErrorShifted * 100:F2} %)"); + output.WriteLine($" annual DNI: reference {result.ReferenceDirectNormal / 1000:F1}, exported {result.ExportedDirectNormal / 1000:F1} kWh/m2 ({(result.ExportedDirectNormal / result.ReferenceDirectNormal - 1) * 100:+0.00;-0.00;0.00} %)"); + output.WriteLine($" annual horizontal beam: GHI-DHI {result.HorizontalBeam / 1000:F1}, exported DNI*cosZ {result.ExportedBeamOnHorizontal / 1000:F1} kWh/m2 ({(result.ExportedBeamOnHorizontal / result.HorizontalBeam - 1) * 100:+0.00;-0.00;0.00} %)"); + output.WriteLine($" sun-up hours n={relativeErrors.Count}: median DNI error {result.MedianRelativeError * 100:+0.00;-0.00;0.00} %, 90th percentile |error| {result.Percentile90AbsoluteRelativeError * 100:F2} %"); + for (int k = 0; k < 4; k++) + output.WriteLine($" facade {AlmanacSun.FacadeNames[k]} annual beam: reference {result.ReferenceFacadeBeam[k] / 1000:F1}, exported {result.ExportedFacadeBeam[k] / 1000:F1} kWh/m2"); + + return result; + } + + private sealed class EpwRow + { + public int FieldCount; + public int Month; + public int Day; + public int Hour; + public double DryBulb; + public double DewPoint; + public double ExtraterrestrialHorizontal; + public double ExtraterrestrialDirectNormal; + public double HorizontalInfrared; + public double GlobalHorizontal; + public double DirectNormal; + public double DiffuseHorizontal; + public double TotalSkyCover; + public double OpaqueSkyCover; + } + + /// + /// Independent solar position: Michalsky (1988), "The Astronomical Almanac's algorithm for + /// approximate solar position (1950-2050)", Solar Energy 40(3) - mean longitude, anomaly, + /// ecliptic longitude, right ascension and Greenwich mean sidereal time. No equation of time. + /// + private static class AlmanacSun + { + public static readonly string[] FacadeNames = { "N", "E", "S", "W" }; + + public sealed class HourGeometry + { + /// Hour mean of max(0, cos zenith). + public double MeanCosZenith; + + /// Hour mean of max(0, cos incidence) on vertical facades facing N, E, S, W, with the sun up. + public double[] MeanCosIncidence = new double[4]; + } + + public static double MeanCosZenith(int year, int month, int day, int epwHour, double latitude, double longitude, double timeZone) + { + return Hour(year, month, day, epwHour, latitude, longitude, timeZone).MeanCosZenith; + } + + /// + /// Geometry averaged over the EPW hour ending at :00 local + /// standard time (shifted by ), from sub-steps. + /// + public static HourGeometry Hour(int year, int month, int day, int epwHour, double latitude, double longitude, double timeZone, double shiftMinutes = 0, int steps = 600) + { + DateTime start = new DateTime(year, month, day, 0, 0, 0, DateTimeKind.Utc).AddHours(epwHour - 1 - timeZone).AddMinutes(shiftMinutes); + HourGeometry result = new HourGeometry(); + for (int i = 0; i < steps; i++) + { + (double cosZenith, double east, double north) = Sun(start.AddHours((i + 0.5) / steps), latitude, longitude); + if (cosZenith <= 0) + continue; + + result.MeanCosZenith += cosZenith / steps; + + // Outward normals of vertical facades: N (0,1), E (1,0), S (0,-1), W (-1,0). + result.MeanCosIncidence[0] += System.Math.Max(0, north) / steps; + result.MeanCosIncidence[1] += System.Math.Max(0, east) / steps; + result.MeanCosIncidence[2] += System.Math.Max(0, -north) / steps; + result.MeanCosIncidence[3] += System.Math.Max(0, -east) / steps; + } + + return result; + } + + /// + /// Unit sun vector (up, east, north) at the given instant. + /// + private static (double, double, double) Sun(DateTime utc, double latitude, double longitude) + { + double n = (utc - new DateTime(2000, 1, 1, 12, 0, 0, DateTimeKind.Utc)).TotalDays; + double hours = utc.TimeOfDay.TotalHours; + + double meanLongitude = Normalise(280.460 + 0.9856474 * n); + double meanAnomaly = Radians(Normalise(357.528 + 0.9856003 * n)); + double eclipticLongitude = Radians(Normalise(meanLongitude + 1.915 * System.Math.Sin(meanAnomaly) + 0.020 * System.Math.Sin(2 * meanAnomaly))); + double obliquity = Radians(23.439 - 0.0000004 * n); + + double rightAscension = System.Math.Atan2(System.Math.Cos(obliquity) * System.Math.Sin(eclipticLongitude), System.Math.Cos(eclipticLongitude)); + double declination = System.Math.Asin(System.Math.Sin(obliquity) * System.Math.Sin(eclipticLongitude)); + + // n includes the fraction of the day, which supplies the sidereal excess of the UT hours. + double greenwichMeanSiderealTime = (6.697375 + 0.0657098242 * n + hours) % 24.0; + double localMeanSiderealTime = greenwichMeanSiderealTime + longitude / 15.0; + double hourAngle = Radians(localMeanSiderealTime * 15.0) - rightAscension; + + double latitudeRadians = Radians(latitude); + double up = System.Math.Sin(latitudeRadians) * System.Math.Sin(declination) + System.Math.Cos(latitudeRadians) * System.Math.Cos(declination) * System.Math.Cos(hourAngle); + double east = -System.Math.Cos(declination) * System.Math.Sin(hourAngle); + double north = System.Math.Cos(latitudeRadians) * System.Math.Sin(declination) - System.Math.Sin(latitudeRadians) * System.Math.Cos(declination) * System.Math.Cos(hourAngle); + return (up, east, north); + } + + private static double Normalise(double degrees) + { + double result = degrees % 360.0; + return result < 0 ? result + 360.0 : result; + } + + private static double Radians(double degrees) + { + return degrees * System.Math.PI / 180.0; + } + } + + #endregion + } +} diff --git a/SAM/SAM.Weather/Query/DataString.cs b/SAM/SAM.Weather/Query/DataString.cs index 526b36a32..603e4aeeb 100644 --- a/SAM/SAM.Weather/Query/DataString.cs +++ b/SAM/SAM.Weather/Query/DataString.cs @@ -3,6 +3,7 @@ using System; using System.Collections.Generic; +using System.Globalization; using System.Linq; namespace SAM.Weather @@ -23,9 +24,11 @@ public static string DataString(this WeatherData weatherData) if (years == null) return null; + double timeZone = weatherData.EPWTimeZone(); + List values = new List(); foreach (int year in years) - values.Add(DataString(year, weatherData[year], weatherData.Elevtion)); + values.Add(DataString(year, weatherData[year], weatherData.Elevtion, weatherData.Latitude, weatherData.Longitude, timeZone)); return string.Join("\n", values); } @@ -47,11 +50,28 @@ private static double StandardAtmosphericPressure(double elevation) /// /// Generates a string of data for a given year and WeatherYear object. + /// Without a location the direct-normal radiation cannot be derived from the horizontal + /// components and is written as missing; use the overload with location to derive it. /// /// The year to generate data for. /// The WeatherYear object to generate data from. /// A string of data for the given year and WeatherYear object. public static string DataString(int year, WeatherYear weatherYear, double elevation = double.NaN) + { + return DataString(year, weatherYear, elevation, double.NaN, double.NaN, double.NaN); + } + + /// + /// Generates a string of data for a given year and WeatherYear object. + /// + /// The year to generate data for. + /// The WeatherYear object to generate data from. + /// Site elevation [m] used for the standard-atmosphere pressure fallback; NaN falls back to sea level. + /// Site latitude [degrees], north positive. + /// Site longitude [degrees], east positive. + /// Standard time zone [hours from UTC], east positive. + /// A string of data for the given year and WeatherYear object. + public static string DataString(int year, WeatherYear weatherYear, double elevation, double latitude, double longitude, double timeZone) { if (weatherYear == null) return null; @@ -66,7 +86,7 @@ public static string DataString(int year, WeatherYear weatherYear, double elevat DateTime dateTime = new DateTime(year, 1, 1); dateTime = dateTime.AddDays(i); - values[i] = DataString(year, dateTime.Month, dateTime.Day, weatherDay, elevation); + values[i] = DataString(year, dateTime.Month, dateTime.Day, weatherDay, elevation, latitude, longitude, timeZone); } return string.Join("\n", values); @@ -74,6 +94,7 @@ public static string DataString(int year, WeatherYear weatherYear, double elevat /// /// Generates a string of data for a given day and weatherDay. + /// Without a location the direct-normal radiation is written as missing. /// /// The year of the data. /// The month of the data. @@ -81,6 +102,23 @@ public static string DataString(int year, WeatherYear weatherYear, double elevat /// The weatherDay of the data. /// A string of data for the given day and weatherDay. public static string DataString(int year, int month, int day, WeatherDay weatherDay, double elevation = double.NaN) + { + return DataString(year, month, day, weatherDay, elevation, double.NaN, double.NaN, double.NaN); + } + + /// + /// Generates a string of data for a given day and weatherDay. + /// + /// The year of the data. + /// The month of the data. + /// The day of the data. + /// The weatherDay of the data. + /// Site elevation [m] used for the standard-atmosphere pressure fallback; NaN falls back to sea level. + /// Site latitude [degrees], north positive. + /// Site longitude [degrees], east positive. + /// Standard time zone [hours from UTC], east positive. + /// A string of data for the given day and weatherDay. + public static string DataString(int year, int month, int day, WeatherDay weatherDay, double elevation, double latitude, double longitude, double timeZone) { if (weatherDay == null) return null; @@ -88,7 +126,7 @@ public static string DataString(int year, int month, int day, WeatherDay weather List values = new List(); for (int i = 0; i < 24; i++) { - string value = DataString(year, month, day, i, weatherDay, elevation); + string value = DataString(year, month, day, i, weatherDay, elevation, latitude, longitude, timeZone); values.Add(value == null ? string.Empty : value); } @@ -97,6 +135,7 @@ public static string DataString(int year, int month, int day, WeatherDay weather /// /// Generates a data string from the given year, month, day, hour, and WeatherDay object. + /// Without a location the direct-normal radiation is written as missing. /// /// The year. /// The month. @@ -106,6 +145,42 @@ public static string DataString(int year, int month, int day, WeatherDay weather /// Site elevation [m] used for the standard-atmosphere pressure fallback; NaN falls back to sea level. /// A data string. public static string DataString(int year, int month, int day, int hour, WeatherDay weatherDay, double elevation = double.NaN) + { + return DataString(year, month, day, hour, weatherDay, elevation, double.NaN, double.NaN, double.NaN); + } + + /// + /// Generates the EPW data record for one hour. + /// + /// Radiation: SAM stores horizontal components (global = direct + diffuse on the horizontal, + /// see ). The record carries global + /// horizontal, direct normal and diffuse horizontal radiation; direct normal is the horizontal + /// beam divided by the hour-mean cosine of the solar zenith at the site (EPW hour = the hour + /// ending at hour + 1, local standard time), so it needs latitude, longitude and time + /// zone. It never exceeds what a pure Rayleigh atmosphere transmits, which only matters with the + /// sun near the horizon, where global - diffuse is ill-conditioned; the stated global and + /// diffuse are written unchanged. The extraterrestrial fields are astronomical. A component that + /// is neither stated nor derivable is written as the EPW missing value 9999, never as 0. + /// + /// + /// Sky: SAM does not hold horizontal infrared radiation, so it is written as missing (9999); + /// EnergyPlus then derives it from dry-bulb, dew-point and opaque sky cover. SAM cloud cover + /// [0-1] is written in tenths to total and opaque sky cover (SAM holds one cloud amount; the + /// opaque value is that amount, not a separate observation). A dew point that is not stated is + /// derived from dry-bulb and relative humidity. + /// + /// + /// The year. + /// The month. + /// The day. + /// The hour. + /// The WeatherDay object. + /// Site elevation [m] used for the standard-atmosphere pressure fallback; NaN falls back to sea level. + /// Site latitude [degrees], north positive. + /// Site longitude [degrees], east positive. + /// Standard time zone [hours from UTC], east positive. + /// A data string. + public static string DataString(int year, int month, int day, int hour, WeatherDay weatherDay, double elevation, double latitude, double longitude, double timeZone) { if (weatherDay == null) return null; @@ -117,27 +192,48 @@ public static string DataString(int year, int month, int day, int hour, WeatherD double relativeHumidity = weatherDay[WeatherDataType.RelativeHumidity, hour]; double atmosphericPressure = weatherDay[WeatherDataType.AtmosphericPressure, hour]; + // EPW field 8 is the dew point; SAM keeps the EPW dew point in WetBulbTemperature + // (see TryGetData). Sources without it (TAS weather) get it derived from DBT and RH, + // because EnergyPlus uses it for the sky emissivity when the infrared field is missing. + double dewPointTemperature = wetBulbTemperature; + if (double.IsNaN(dewPointTemperature)) + dewPointTemperature = EPWDewPointTemperature(dryBulbTemperature, relativeHumidity); + + double globalSolarRadiation = weatherDay[WeatherDataType.GlobalSolarRadiation, hour]; double directSolarRadiation = weatherDay[WeatherDataType.DirectSolarRadiation, hour]; double diffuseSolarRadiation = weatherDay[WeatherDataType.DiffuseSolarRadiation, hour]; - if (double.IsNaN(directSolarRadiation) && !double.IsNaN(diffuseSolarRadiation)) - { - double globalSolarRadiation = weatherDay[WeatherDataType.GlobalSolarRadiation, hour]; - if (!double.IsNaN(globalSolarRadiation)) - directSolarRadiation = globalSolarRadiation - diffuseSolarRadiation; - } + // Complete the horizontal triple from global = direct + diffuse; stated values win. + if (double.IsNaN(directSolarRadiation) && !double.IsNaN(globalSolarRadiation) && !double.IsNaN(diffuseSolarRadiation)) + directSolarRadiation = globalSolarRadiation - diffuseSolarRadiation; + + if (double.IsNaN(diffuseSolarRadiation) && !double.IsNaN(globalSolarRadiation) && !double.IsNaN(directSolarRadiation)) + diffuseSolarRadiation = globalSolarRadiation - directSolarRadiation; + + if (double.IsNaN(globalSolarRadiation) && !double.IsNaN(directSolarRadiation) && !double.IsNaN(diffuseSolarRadiation)) + globalSolarRadiation = directSolarRadiation + diffuseSolarRadiation; + + double extraterrestrialHorizontalRadiation = double.NaN; + double extraterrestrialDirectNormalRadiation = double.NaN; + double directNormalRadiation = double.NaN; + + double pressure = double.IsNaN(atmosphericPressure) ? StandardAtmosphericPressure(elevation) : atmosphericPressure; - if (double.IsNaN(diffuseSolarRadiation) && !double.IsNaN(directSolarRadiation)) + DateTime date = new DateTime(year, month, day); + if (TryGetEPWSolarGeometry(date.DayOfYear, DateTime.IsLeapYear(year) ? 366 : 365, hour, latitude, longitude, timeZone, pressure, out double meanCosZenith, out double extraterrestrialNormalIrradiance, out double limitDirectNormalIrradiance)) { - double globalSolarRadiation = weatherDay[WeatherDataType.GlobalSolarRadiation, hour]; - if (!double.IsNaN(globalSolarRadiation)) - diffuseSolarRadiation = globalSolarRadiation - directSolarRadiation; + extraterrestrialHorizontalRadiation = extraterrestrialNormalIrradiance * meanCosZenith; + extraterrestrialDirectNormalRadiation = meanCosZenith > 0 ? extraterrestrialNormalIrradiance : 0; + directNormalRadiation = EPWDirectNormalRadiation(directSolarRadiation, meanCosZenith, limitDirectNormalIrradiance); } double windDirection = weatherDay[WeatherDataType.WindDirection, hour]; double windSpeed = weatherDay[WeatherDataType.WindSpeed, hour]; double cloudCover = weatherDay[WeatherDataType.CloudCover, hour]; + // EPW sky cover is in tenths (0-10); SAM cloud cover is 0-1 (TryGetData divides by 10). + double skyCover = double.IsNaN(cloudCover) ? double.NaN : Math.Round(Math.Max(0, Math.Min(1, cloudCover)) * 10.0, 6); + values.Add(year.ToString()); //Year values.Add(month.ToString()); //Month values.Add(day.ToString()); //Day @@ -145,23 +241,50 @@ public static string DataString(int year, int month, int day, int hour, WeatherD values.Add(0.ToString()); //Minute values.Add("C9C9C9C9*0?9?9?9?9?9?9?9A7A7A7A7A7A7*0E8*0*0"); //Flags values.Add(double.IsNaN(dryBulbTemperature) ? 0.ToString() : dryBulbTemperature.ToString()); //Dry Bulb Temperature - values.Add(double.IsNaN(wetBulbTemperature) ? 0.ToString() : wetBulbTemperature.ToString()); //Wet Bulb Temperature + values.Add(EPWValue(dewPointTemperature, "99.9")); //Dew Point Temperature values.Add(double.IsNaN(relativeHumidity) ? 0.ToString() : relativeHumidity.ToString()); //Relative Humidity // EnergyPlus rejects station pressure outside [31000, 120000] Pa with a severe // ReadEPlusWeatherForDay error — a missing value must never be written as 0. // Fall back to the standard-atmosphere pressure at the site elevation (invariant // culture: a decimal-comma locale would corrupt the comma-separated row). - values.Add(double.IsNaN(atmosphericPressure) ? StandardAtmosphericPressure(elevation).ToString(System.Globalization.CultureInfo.InvariantCulture) : atmosphericPressure.ToString()); //Atmospheric Pressure - values.Add(string.Join(",", Enumerable.Repeat(0.ToString(), 4))); // Solar - values.Add(double.IsNaN(directSolarRadiation) ? 0.ToString() : directSolarRadiation.ToString()); //Direct Solar Radiation - values.Add(double.IsNaN(diffuseSolarRadiation) ? 0.ToString() : diffuseSolarRadiation.ToString()); //Diffuse Solar Radiation + values.Add(double.IsNaN(atmosphericPressure) ? StandardAtmosphericPressure(elevation).ToString(CultureInfo.InvariantCulture) : atmosphericPressure.ToString()); //Atmospheric Pressure + values.Add(EPWValue(extraterrestrialHorizontalRadiation, "9999")); //Extraterrestrial Horizontal Radiation + values.Add(EPWValue(extraterrestrialDirectNormalRadiation, "9999")); //Extraterrestrial Direct Normal Radiation + values.Add("9999"); //Horizontal Infrared Radiation Intensity - not held by SAM; EnergyPlus derives it + values.Add(EPWRadiationValue(globalSolarRadiation)); //Global Horizontal Radiation + values.Add(EPWRadiationValue(directNormalRadiation)); //Direct Normal Radiation + values.Add(EPWRadiationValue(diffuseSolarRadiation)); //Diffuse Horizontal Radiation values.Add(string.Join(",", Enumerable.Repeat(0.ToString(), 4))); // Illumination values.Add(double.IsNaN(windDirection) ? 0.ToString() : windDirection.ToString()); //Wind Direction values.Add(double.IsNaN(windSpeed) ? 0.ToString() : windSpeed.ToString()); //Wind Speed - values.Add(double.IsNaN(cloudCover) ? 0.ToString() : cloudCover.ToString()); //Cloud Cover - values.Add(string.Join(",", Enumerable.Repeat(0.ToString(), 9))); // Sky, Precipiatation, Snow etc. + values.Add(EPWValue(skyCover, "99")); //Total Sky Cover + values.Add(EPWValue(skyCover, "99")); //Opaque Sky Cover + values.Add(string.Join(",", Enumerable.Repeat(0.ToString(), 8))); // Visibility, Ceiling, Present Weather, Precipitable Water, Aerosol, Snow return string.Join(",", values); } + + /// + /// EPW text of a value (invariant culture, at most 2 decimals) or the missing code when NaN. + /// + private static string EPWValue(double value, string missing) + { + if (double.IsNaN(value) || double.IsInfinity(value)) + return missing; + + return Math.Round(value, 2).ToString(CultureInfo.InvariantCulture); + } + + /// + /// EPW text of a radiation value [Wh/m2]: negative values (global - diffuse rounding) are 0, + /// NaN is the missing code 9999. + /// + private static string EPWRadiationValue(double value) + { + if (double.IsNaN(value) || double.IsInfinity(value)) + return "9999"; + + return EPWValue(Math.Max(0, value), "9999"); + } } } diff --git a/SAM/SAM.Weather/Query/EPWSolarGeometry.cs b/SAM/SAM.Weather/Query/EPWSolarGeometry.cs new file mode 100644 index 000000000..0b250abdb --- /dev/null +++ b/SAM/SAM.Weather/Query/EPWSolarGeometry.cs @@ -0,0 +1,188 @@ +// SPDX-License-Identifier: LGPL-3.0-or-later +// Copyright (c) 2020–2026 Michal Dengusiak & Jakub Ziolkowski and contributors + +using System; + +namespace SAM.Weather +{ + public static partial class Query + { + /// + /// Solar constant [W/m2] used for the extraterrestrial irradiance of the EPW record + /// (the value the EnergyPlus weather converter and the TMY2/TMY3 EPW sources use). + /// + private const double EPWSolarConstant = 1367.0; + + /// + /// Number of equal sub-intervals used to average the solar geometry over one EPW hour. + /// + private const int EPWSolarSubIntervals = 60; + + /// + /// Hour-averaged solar geometry of one EPW record, used to convert SAM horizontal beam + /// radiation into EPW direct-normal radiation. + /// + /// An EPW hour h (1-24) holds the radiation received during the hour that ends at + /// h:00 local standard time (EnergyPlus Auxiliary Programs, EPW data dictionary: + /// "during the number of minutes preceding the time indicated"); SAM stores it at hour index + /// h - 1 (see ). The interval is therefore + /// [hour, hour + 1) local standard time for SAM . + /// + /// + /// is the mean over that interval of max(0, cos zenith), so a + /// constant direct-normal irradiance DNI over the hour puts DNI * meanCosZenith + /// on the horizontal - the energy relation the conversion inverts. Sun position: NOAA + /// general solar position equations (Spencer 1971 Fourier series for declination and + /// equation of time), no refraction; accuracy is a few arc-minutes, well below the + /// sub-hour averaging. + /// + /// + /// is the largest direct-normal irradiance + /// that can reach the ground in the hour: the beam of a pure Rayleigh (aerosol-, water- and + /// ozone-free) atmosphere, I0 exp(-dR(m') m'), with the Kasten and Young (1989) + /// relative air mass m, m' = m p / 101325 and the Kasten (1996) integral Rayleigh + /// optical thickness dR, averaged over the hour with cos-zenith weights (the constant + /// direct normal that would give that beam on the horizontal). A real atmosphere only + /// attenuates more, so the value is a physical bound, not a clear-sky model. + /// + /// + /// Day of year, 1-based. + /// Days in the calendar year of the record (365 or 366). + /// SAM hour index 0-23 (EPW hour minus one). + /// Latitude [degrees], north positive. + /// Longitude [degrees], east positive. + /// Standard time zone [hours from UTC], east positive (EPW LOCATION convention). + /// Station pressure [Pa] for the air-mass correction; NaN uses 101325 Pa. + /// Hour mean of max(0, cos zenith) [-]. + /// Extraterrestrial direct-normal irradiance for the day [W/m2]. + /// Rayleigh-atmosphere bound on the hour's direct normal [W/m2]; 0 when the sun is down all hour. + /// False when the location or time zone is not stated (NaN) or out of range. + private static bool TryGetEPWSolarGeometry(int dayOfYear, int daysInYear, int hour, double latitude, double longitude, double timeZone, double pressure, out double meanCosZenith, out double extraterrestrialNormalIrradiance, out double limitDirectNormalIrradiance) + { + meanCosZenith = double.NaN; + extraterrestrialNormalIrradiance = double.NaN; + limitDirectNormalIrradiance = double.NaN; + + if (double.IsNaN(latitude) || double.IsNaN(longitude) || double.IsNaN(timeZone)) + return false; + + if (latitude < -90 || latitude > 90 || longitude < -180 || longitude > 180 || timeZone < -14 || timeZone > 14) + return false; + + if (hour < 0 || hour > 23 || dayOfYear < 1 || dayOfYear > daysInYear) + return false; + + double pressureRatio = double.IsNaN(pressure) || pressure <= 0 ? 1.0 : pressure / 101325.0; + + double dayAngle = 2.0 * Math.PI * (dayOfYear - 1) / daysInYear; + extraterrestrialNormalIrradiance = EPWSolarConstant * (1.000110 + 0.034221 * Math.Cos(dayAngle) + 0.001280 * Math.Sin(dayAngle) + 0.000719 * Math.Cos(2 * dayAngle) + 0.000077 * Math.Sin(2 * dayAngle)); + + double latitudeRadians = latitude * Math.PI / 180.0; + double sumCosZenith = 0; + double sumLimitBeam = 0; + for (int i = 0; i < EPWSolarSubIntervals; i++) + { + // Mid-point of the sub-interval, in hours of local standard time. + double time = hour + (i + 0.5) / EPWSolarSubIntervals; + + double gamma = 2.0 * Math.PI / daysInYear * (dayOfYear - 1 + (time - 12.0) / 24.0); + + double equationOfTime = 229.18 * (0.000075 + 0.001868 * Math.Cos(gamma) - 0.032077 * Math.Sin(gamma) - 0.014615 * Math.Cos(2 * gamma) - 0.040849 * Math.Sin(2 * gamma)); + + double declination = 0.006918 - 0.399912 * Math.Cos(gamma) + 0.070257 * Math.Sin(gamma) - 0.006758 * Math.Cos(2 * gamma) + 0.000907 * Math.Sin(2 * gamma) - 0.002697 * Math.Cos(3 * gamma) + 0.00148 * Math.Sin(3 * gamma); + + // True solar time [min]: standard meridian is 15 degrees per hour of time zone. + double trueSolarTime = time * 60.0 + equationOfTime + 4.0 * longitude - 60.0 * timeZone; + double hourAngle = (trueSolarTime / 4.0 - 180.0) * Math.PI / 180.0; + + double cosZenith = Math.Sin(latitudeRadians) * Math.Sin(declination) + Math.Cos(latitudeRadians) * Math.Cos(declination) * Math.Cos(hourAngle); + if (cosZenith <= 0) + continue; + + sumCosZenith += cosZenith; + + double airMass = RelativeAirMass(cosZenith) * pressureRatio; + sumLimitBeam += extraterrestrialNormalIrradiance * Math.Exp(-RayleighOpticalThickness(airMass) * airMass) * cosZenith; + } + + meanCosZenith = sumCosZenith / EPWSolarSubIntervals; + limitDirectNormalIrradiance = sumCosZenith > 0 ? sumLimitBeam / sumCosZenith : 0; + return true; + } + + /// + /// Relative optical air mass, Kasten and Young (1989), from the cosine of the solar zenith + /// (sun above the horizon). + /// + private static double RelativeAirMass(double cosZenith) + { + double zenith = Math.Acos(Math.Min(1.0, cosZenith)) * 180.0 / Math.PI; + return 1.0 / (cosZenith + 0.50572 * Math.Pow(96.07995 - zenith, -1.6364)); + } + + /// + /// Integral Rayleigh optical thickness, Kasten (1996), for the pressure-corrected air mass. + /// + private static double RayleighOpticalThickness(double airMass) + { + if (airMass <= 20) + return 1.0 / (6.6296 + 1.7513 * airMass - 0.1202 * Math.Pow(airMass, 2) + 0.0065 * Math.Pow(airMass, 3) - 0.00013 * Math.Pow(airMass, 4)); + + return 1.0 / (10.4 + 0.718 * airMass); + } + + /// + /// Direct-normal radiation [Wh/m2] for one EPW hour from the horizontal beam (global minus + /// diffuse) and the hour-averaged solar geometry: DNI = beam / meanCosZenith, limited + /// to the Rayleigh-atmosphere bound of the hour. Horizontal beam in an hour where the sun is + /// below the horizon for the whole interval cannot be represented as a direct-normal value + /// and gives 0. NaN when the horizontal beam or the geometry is not known (written as the EPW + /// missing value). + /// + private static double EPWDirectNormalRadiation(double horizontalBeamRadiation, double meanCosZenith, double limitDirectNormalIrradiance) + { + if (double.IsNaN(horizontalBeamRadiation) || double.IsNaN(meanCosZenith)) + return double.NaN; + + if (horizontalBeamRadiation <= 0 || meanCosZenith <= 0) + return 0; + + double result = horizontalBeamRadiation / meanCosZenith; + if (!double.IsNaN(limitDirectNormalIrradiance) && result > limitDirectNormalIrradiance) + result = limitDirectNormalIrradiance; + + return result; + } + + /// + /// Dew-point temperature [C] from dry-bulb temperature [C] and relative humidity [%] + /// (Magnus form over water, Alduchov and Eskridge 1996 coefficients; within about 0.35 K + /// from -40 C to 50 C). NaN when either input is missing or the humidity is not positive. + /// + private static double EPWDewPointTemperature(double dryBulbTemperature, double relativeHumidity) + { + if (double.IsNaN(dryBulbTemperature) || double.IsNaN(relativeHumidity) || relativeHumidity <= 0) + return double.NaN; + + const double a = 17.625; + const double b = 243.04; + + double value = Math.Log(Math.Min(relativeHumidity, 100.0) / 100.0) + a * dryBulbTemperature / (b + dryBulbTemperature); + return b * value / (a - value); + } + + /// + /// Standard time zone [hours from UTC] stated on the weather data, or NaN when it is not stated. + /// + private static double EPWTimeZone(this WeatherData weatherData) + { + if (weatherData == null) + return double.NaN; + + if (!weatherData.TryGetValue(WeatherDataParameter.TimeZone, out string timeZone) || string.IsNullOrWhiteSpace(timeZone)) + return double.NaN; + + return Core.Query.Double(Core.Query.UTC(timeZone)); + } + } +} diff --git a/THIRD_PARTY.md b/THIRD_PARTY.md index 9232d2da7..975267b15 100644 --- a/THIRD_PARTY.md +++ b/THIRD_PARTY.md @@ -34,3 +34,10 @@ listing comments and termination text. It has no GenOpt banner, copyright header, console transcript or code. - Distribution: no GenOpt code, `genopt.jar` or Java runtime is included in this repository or its builds. The oracle tool `SAM/SAM.Tests.GenOptOracle` runs a locally supplied `genopt.jar`, which is never committed. + +## TMYx weather files (test fixtures only) +- Files: `SAM/SAM.Tests/Resources/Weather/*.epw.gz`, five TMYx EPW files from https://climate.onebuilding.org/, + gzip-compressed without modification. The SHA-256 of each original `.epw` is in `SHA256SUMS.txt` next to them, and + provenance per file is in that folder's `README.md`. +- License: free to use with attribution: "Climatic data provided by climate.onebuilding.org". +- Use: reference data for the EPW writer tests (`WeatherEPWExportTests`). They are not part of any build output. diff --git a/documentation/EPW-Weather-Radiation-F1-PR.md b/documentation/EPW-Weather-Radiation-F1-PR.md new file mode 100644 index 000000000..5e1878150 --- /dev/null +++ b/documentation/EPW-Weather-Radiation-F1-PR.md @@ -0,0 +1,243 @@ +# EPW writer: physically correct radiation and sky fields (F1) + +PR record. PR https://github.com/SAM-BIM/SAM/pull/212, branch `fix/epw-weather-radiation-q4`, base `sow/2026-Q4` +(`618db819`). Tested content SHA: `a2cb16cf`. Production code is unchanged since `b721c5f4`: `a2cb16cf` adds one +test and evidence, and the commits after it change only this record. Not merged. + +## Status +Final review complete (2026-10-19, see "Final review"). No blocking issues; recommended for owner approval and +merge. TAS was not re-run (see Risks). + +## Defect (SAM_UI PR #239, finding F1) +`SAM.Weather.Query.DataString`, which `Convert.ToEPW` and SAM_OpenStudio `EmbeddedAnnualWeatherPath` use, wrote: +- fields 11-14 (extraterrestrial horizontal, extraterrestrial normal, horizontal IR, **global horizontal**) as `0`; +- field 15 (**direct normal**) as `DirectSolarRadiation`, i.e. global minus diffuse, which is the horizontal beam; +- field 23 (total sky cover, in tenths) as the SAM 0-1 fraction, and field 24 (opaque sky cover) as `0`; +- field 8 (dew point) as `0` when the source has none (TAS weather), and missing radiation as `0`. + +Measured on the box model's embedded weather, which comes from the Budapest TMYx EPW: GHI 1276 -> 0 kWh/m2 and +DNI 1305 -> 735 kWh/m2 (the horizontal beam). + +## Source-field semantics (verified in code) +- SAM `WeatherDataType.GlobalSolarRadiation`, `DiffuseSolarRadiation` and `DirectSolarRadiation` are all + **horizontal**: `WeatherHour.CalculatedGlobalSolarRadiation` gives global = direct + diffuse. They are + hourly values in W/m2, which equal Wh/m2 for one hour. +- EPW import (`Query.TryGetData`) reads GHI (field 14) and DHI (field 16) and does **not** read DNI. It stores + the EPW dew point (field 8) under `WetBulbTemperature`, and total sky cover / 10 under `CloudCover`. SAM hour + index `i` is EPW hour `i + 1`. +- TAS weather (SAM_Tas `Convert.ToSAM(TWD.WeatherYear)`) carries global, diffuse, cloud cover (0-1), dry bulb, + RH and wind, with no dew point and no pressure. TBD hour `i` (1-based) maps to SAM index `i - 1`. +- `WeatherYear` has exactly 365 days. Time zone is the `WeatherDataParameter.TimeZone` string ("UTC+01:00"). +- SAM has no solar-position code of its own. The one in SAM_SolarCalculator cannot be referenced, because that + repository depends on SAM. The geometry is therefore implemented privately in SAM.Weather. + +## Transformations (EnergyPlus EPW data dictionary conventions) +| EPW field | Written value | +|---|---| +| 8 dew point | stated value; else Magnus (Alduchov-Eskridge) from DBT and RH; else 99.9 | +| 11 / 12 extraterrestrial horizontal / normal | `I0 * meanCosZ` / `I0` (Spencer eccentricity, 1367 W/m2); 9999 without a location | +| 13 horizontal IR | 9999 (missing). SAM holds no IR; EnergyPlus derives it (Clark-Allen) from DBT, dew point and opaque sky cover | +| 14 global horizontal | stated value; else direct + diffuse; else 9999 | +| 15 direct normal | `beam / meanCosZ`, bounded by the Rayleigh-atmosphere limit; 0 if the sun is down all hour; 9999 without beam or location | +| 16 diffuse horizontal | stated value; else global - direct; else 9999 | +| 23 / 24 total / opaque sky cover | cloud cover x 10 (tenths) to both; 99 if missing | + +- **Time base:** EPW hour `h` holds the radiation of the hour ending at `h:00` local standard time. + `meanCosZ` is the mean of max(0, cos zenith) over `[h-1, h)` from 60 sub-steps. This makes the conversion + energy-consistent, so `DNI * meanCosZ` returns the stated beam, and it handles the sunrise and sunset hours. + Daylight saving does not apply: the writer declares none, and EPW is standard time. +- **Sun position:** NOAA general solar position (Spencer 1971 declination and equation of time), with the + longitude and the time-zone meridian. No refraction. +- **Low-sun bound:** global - diffuse is ill-conditioned at grazing sun. A few W/m2 of beam divided by + `meanCosZ` of about 0.02 gives 1000-1400 W/m2 of DNI. EnergyPlus then puts that DNI on east or west facades + at near-normal incidence. The bound is the beam of a pure Rayleigh atmosphere: Kasten-Young (1989) air mass, + pressure-corrected, with the Kasten (1996) integral Rayleigh optical thickness, averaged with cos-zenith + weights. A real atmosphere only attenuates more, so this is a physical bound rather than a tuned clear-sky + model. It does not bind at high sun: about 1210 W/m2 at air mass 1. The stated global and diffuse are always + written unchanged. + +## Decisions and assumptions +- **Missing is not zero.** A missing value is written with its EPW missing code. EnergyPlus sets missing DNI + and DHI to 0 and reports the count, and it derives a missing IR. +- **IR is not computed in SAM.** Writing 9999 lets EnergyPlus apply its own documented algorithm. The previous + 0 was also treated by EnergyPlus as invalid and derived. Because opaque sky cover was 0 (clear sky), that + gave a sky about 3 K too cold (see the box evidence). It was not a 0 W/m2 sky, as the PR #239 record says. +- **Opaque sky cover = SAM cloud cover.** SAM holds one cloud amount. TAS uses that amount for long-wave sky + exchange, which is the role opaque cover has in EnergyPlus. For EPW-sourced data it is the total cover, and + the opaque value of the original file cannot be recovered. This is a stated proxy. +- **No time zone means no DNI.** The time zone is not guessed from longitude. EPW and TAS imports both set it. +- **Public API is unchanged.** All four existing `DataString` overloads keep their signatures. Without a + location they now write DNI as 9999 instead of the horizontal beam. New overloads take latitude, longitude + and time zone. `DataString(WeatherData)` passes them. +- **Versioning:** none in the repository. `SAMVersion` is stamped by CI (`Directory.Build.props`, `build.yml`), + and there is no Grasshopper component change. + +## Independent verification +Five unmodified TMYx EPWs are committed under `SAM/SAM.Tests/Resources/Weather` (gzip, SHA-256 of the +originals, provenance README, attribution in `THIRD_PARTY.md`): Budapest (incident, ERA5-based), Fairbanks +(64.8 N), Sydney (southern hemisphere), Singapore (tropics), Phoenix (clear). Each is imported into SAM and +re-exported. The exported DNI is compared with the dataset's own DNI, using an independent geometry in the test: +the Michalsky (1988) Astronomical Almanac algorithm, which goes through sidereal time and has no equation of +time. Full numbers: `documentation/evidence/epw-weather-f1/reference-epw-statistics.txt`. + +| Reference | Own closure (EPW time base / 30 min earlier) | Annual DNI exported vs reference | Sun-up median / P90 error | Facade beam N/E/S/W | +|---|---|---|---|---| +| Sydney | 0.40 % / 8.15 % | +0.22 % | +0.24 % / 0.59 % | within 0.3 % | +| Singapore | 0.54 % / 7.78 % | +0.28 % | +0.27 % / 0.57 % | within 0.4 % | +| Fairbanks | 1.32 % / 6.07 % | +0.23 % | +0.17 % / 0.61 % | within 0.4 % | +| Phoenix | **8.85 % / 2.34 %** | +6.7 % | -1.49 % / 20.0 % | E -10 %, W +27 % | +| Budapest | **9.31 % / 4.26 %** | +14.7 % | -1.36 % / 20.9 % | E -7 %, W +39 % | + +- The writer's time base is the one the self-consistent datasets use. Its extraterrestrial-horizontal values + match the Sydney and Phoenix files' own field hour by hour (for example 788.18 against 788). +- **Phoenix and Budapest are internally inconsistent.** Their DNI closes with their own GHI - DHI only when the + geometry is taken about 30 minutes earlier, and Budapest's own extraterrestrial column is a further 30 minutes + later. No conversion that follows the EPW convention can reproduce such a DNI from GHI - DHI. The tests + measure this rather than assume it. For these files the tests require energy closure instead: the exported + `DNI * cosZ` equals GHI - DHI within -0.4 %, which is the beam above the Rayleigh bound. This is a property of + the source data, not of the writer. +- Synthetic hours: a known DNI at solar noon, at low sun, in winter and summer, in both hemispheres, 13-14 + degrees off the zone meridian and with a half-hour zone (New Delhi) is recovered within 1 %, or within 3 % in + hours with the sun below about 5 degrees, where the two algorithms' few arc-minutes are a percent of the + altitude. Night, sun-below-horizon beam, the low-sun bound, high sun untouched, missing data, partial + components, TAS-style weather, leap year (8760 records, geometry of the written date) and export, import and + export stability are also covered. +- **Fail-first:** with the original `DataString.cs` and the same tests, 34 of 38 fail. The 4 that pass are + fixture hashes, the unchanged importer, night zeros and the leap-year row count. + +## Box-model check (machine-local, EnergyPlus 25.2.0 from Ladybug Tools) +The producer's `in.idf` from PR #239 experiment e0 is used, with only `Version` changed to 25.2 and seven site +and surface output variables added. Three weather files: the old export, the corrected export (the box model's +own embedded `WeatherData` through the new writer) and the explicit Budapest EPW. Committed summaries: +`documentation/evidence/epw-weather-f1/` (`eplus-box-*.txt`, `weather-field-summary.txt`, `run-inputs.sha256`). +The old-export run reproduces the PR #239 E+ 25.1 values exactly (2810.0183 / 3832.5030 kWh), and the explicit run +reproduces e2. + +| Annual (box) | Old export | Corrected export | Explicit EPW | +|---|---|---|---| +| Window transmitted solar [kWh] | 2592.6 (-28.5 %) | 3755.1 (+3.6 %) | 3623.7 | +| Site direct solar, mean [W/m2] | 83.6 | 162.1 (+8.9 %) | 148.8 | +| Site horizontal IR, mean [W/m2] | 296.4 | 319.6 (+2.9 %) | 310.5 | +| Site sky temperature, mean [C] | -4.83 | 0.28 | -1.87 | +| Ideal Loads heating [kWh] | 2810.0 (+29.3 %) | 2081.9 (-4.2 %) | 2173.9 | +| Ideal Loads sensible cooling [kWh] | 3432.6 (-17.1 %) | 4311.7 (+4.1 %) | 4140.8 | + +The remaining embedded-vs-explicit difference comes from the Budapest file's internal half-hour inconsistency +(direct) and from Clark-Allen IR with opaque = total cover against the file's own IR column. It is not an +acceptance criterion, and engine agreement between TAS and EnergyPlus is not claimed. + +Reproduce: run `energyplus.exe -w .epw -d . in.idf` in each `C:\TasOut\f1-epw-weather-20261019\run-*` +folder, then `eso-summary.ps1 -Eso eplusout.eso` (PR #239 tool). These are machine-local. Copy that folder and +`C:\TasOut\licensed-e2e-20261019` to continue elsewhere. + +## Final review (2026-10-19) + +### 1. Rayleigh bound: how often it binds and what it changes +Test `RayleighBound_BindsOnlyAtGrazingSun_AndRemovesLittleEnergy`, evidence +`documentation/evidence/epw-weather-f1/rayleigh-bound-activation.txt`. It reads the writer's own hour geometry back +from the export (ETR horizontal / ETR normal), so "binds" counts exactly the hours the bound changed (more than 1 %). +It compares against the same conversion without the bound. Facade beam uses the independent incidence geometry. + +| Reference | Hours bound / beam hours | Highest sun when bound | Annual DNI | Horizontal beam | Facade beam N / E / S / W | +|---|---|---|---|---|---| +| Sydney, Singapore, Fairbanks | **0** | - | 0.00 % | 0.00 % | 0.00 % on all four | +| Budapest | 170 / 4217 (4.0 %) | ~4.7 deg | -51.9 % | -0.18 % | -23.1 / 0.0 / -4.5 / -15.5 % | +| Phoenix | 251 / 4217 (6.0 %) | ~20.7 deg | -3.6 % | -0.38 % | -11.9 / 0.0 / -1.4 / -6.3 % | + +- For the three self-consistent datasets the bound never binds, so the exported DNI is the unbounded energy-consistent + value. The test asserts this. +- It binds only for the two files whose DNI and GHI/DHI are half an hour apart. There, the stated beam at low sun + implies impossible DNI: up to 2450 W/m2 at ~21 deg in Phoenix, where a pure Rayleigh atmosphere allows about + 1080 W/m2, and up to 4.7e5 W/m2 in Budapest in the minutes after sunrise. Without the bound, Budapest's annual DNI + would be 3109 instead of 1495 kWh/m2, and its west facade beam 567 instead of 479 kWh/m2 (the dataset's own DNI + gives 345). The bound removes at most 0.38 % of the stated horizontal beam. +- East facades are unaffected (0.00 %): the binding hours are afternoon and sunset hours. +- Separately, by rule, beam stated while the sun is below the horizon all hour becomes DNI 0. That is 163 hours and + 0.51 kWh/m2 in Budapest, and none elsewhere. +- **What the -51.9 % means:** it compares Budapest's exported DNI with the same conversion *without* the bound + (3109 kWh/m2). That unbounded total is inflated by grazing-sun hours in which a few W/m2 of beam is divided by a + near-zero cosine. It is **not** a change against the source data: the exported annual DNI (1495 kWh/m2) is + +14.7 % against the file's own DNI (1305), and the stated horizontal beam is kept within 0.18 %. +- **The bound is a physical safeguard, not a correction of the half-hour mismatch.** It never shifts, rescales or + redistributes radiation in time. It only limits an hour's DNI to what a pure Rayleigh atmosphere can transmit at + that sun position. The mismatch itself is not corrected, and it remains in the exported DNI as the east/west + asymmetry recorded in "Independent verification". The source `WeatherData` is never modified: the writer only + reads it through side-effect-free getters, and the stated global and diffuse are written unchanged. +- **Conclusion:** the bound is inert for physically consistent data and removes only physically impossible beam + from inconsistent data. Accepted. + +### 2. Opaque-sky-cover stand-in +The stand-in is documented in the code (`DataString` remarks), in "Decisions and assumptions" above and in the PR +description. It was isolated with two more EnergyPlus box runs (`eplus-box-corrected-refopaque.txt`, +`eplus-box-corrected-refIR.txt`): the corrected export with the reference file's real opaque cover, and with the +reference file's own IR column. + +| Box, annual | Stand-in (opaque = total) | Real opaque cover | Dataset IR | Explicit EPW | +|---|---|---|---|---| +| Horizontal IR, mean [W/m2] | 319.6 | 322.2 | 310.5 | 310.5 | +| Sky temperature, mean [C] | 0.28 | 0.88 | -1.87 | -1.87 | +| Heating [kWh] | 2081.9 | 2066.9 (-0.7 %) | 2101.1 (+0.9 %) | 2173.9 | +| Sensible cooling [kWh] | 4311.7 | 4327.9 (+0.4 %) | 4238.4 (-1.7 %) | 4140.8 | + +- Using the stand-in instead of the real opaque cover changes IR by 2.6 W/m2 and loads by less than 1 %. In this + ERA5-based file the opaque cover is on average larger than the total (6.56 against 5.80 tenths), so the "real" + value is itself unreliable. +- EnergyPlus's Clark-Allen IR against the dataset's IR column is a larger effect (9 W/m2, up to 1.7 % of loads), and + it is EnergyPlus's documented behaviour for a missing IR. +- With the dataset's IR the residual against the explicit EPW is -3.4 % heating and +2.4 % cooling. That is the + solar part, from the file's half-hour DNI inconsistency. +- **Conclusion:** documented, small, and better than the alternatives (opaque 0 means a clear sky; 99 means an + unknown EnergyPlus substitute). Accepted. + +### 3. CI, tests and mergeability +- PR #212 at `90822fac`: GitHub `build (Release)`, `test (Release)` and `spdx` pass. Mergeable, state CLEAN, base + `618db819` unchanged, no reviews or comments, ready for review. The final head's checks are in the PR description. +- Local at `a2cb16cf`: `SAM.Tests` 3726/3726 in Release and in Debug, which includes the 43 `WeatherEPWExportTests` + cases. Production code is identical to `b721c5f4`, whose sandboxed `SAM.sln` Release rebuild had 0 errors. + +**Blocking issues: none.** Recommended for owner approval and merge. + +## Files +- `SAM/SAM.Weather/Query/DataString.cs`: field mapping, missing codes and location overloads. +- `SAM/SAM.Weather/Query/EPWSolarGeometry.cs` (new): hour-mean geometry, Rayleigh bound, dew point and time zone. +- `SAM/SAM.Tests/WeatherEPWExportTests.cs` (new, 43 test cases) and `SAM/SAM.Tests/SAM.Tests.csproj` (copies fixtures). +- `SAM/SAM.Tests/Resources/Weather/` (new): 5 `.epw.gz`, `SHA256SUMS.txt` and `README.md`. +- `THIRD_PARTY.md`: TMYx attribution. +- `documentation/evidence/epw-weather-f1/` (new) and this record. + +## Validation (local, Windows, .NET 8 SDK) +- At `a2cb16cf`: `dotnet test SAM/SAM.Tests/SAM.Tests.csproj -c Release`: **3726/3726** passed. With `-c Debug`: + 3726/3726. +- At `b721c5f4` (same production code): `dotnet test` gave 3721/3721 in Release and Debug, before the bound test was + added. `MSBuild SAM.sln /t:Restore;Rebuild /p:Configuration=Release /m:1 /p:UseSharedCompilation=false`: 0 errors. + It ran with `APPDATA`/`USERPROFILE` redirected to a sandbox so the Grasshopper post-build copies could not touch + the installed SAM. +- `git diff --check` is clean, and the new `.cs` files carry the SPDX header. +- GitHub checks: green at `90822fac`; see the PR description for the final head. + +## Unresolved issues and risks +- **Consumers see changed weather.** Embedded-weather OpenStudio and Both runs get more solar and a warmer sky. + Results from embedded weather (including PR #239 B1a) change and are not comparable with earlier runs. The + explicit-EPW route is unchanged: the file is passed through, and the SAM importer is untouched (test + `ReferenceEpw_Import_IsUnchanged`). +- **TAS embedded path not re-run.** No TAS licence run was done here. Whether the TAS route reads this EPW + (PR #239: "Tas then simulates with no solar gain") is to be confirmed in the next licensed run. +- **Internally inconsistent sources** (above): DNI cannot be reconstructed better than the source's own + consistency allows, and east/west facade beam can differ by tens of percent from the dataset's own DNI. + Mitigation: use the explicit EPW when the original file is available, since SAM does not retain DNI. +- **Leap years:** `WeatherYear` holds 365 days. In a leap year, day index 59 is written as 29 February and + 31 December is not written (existing labelling, unchanged). The header declares no leap day. +- **Out of scope, existing behaviour:** DBT, RH, wind and pressure use the current culture, so a decimal-comma + locale breaks the row. Missing DBT, RH and wind are still written as 0. The importer does not recognise EPW + missing codes (a 9999 GHI is imported as 9999). Illuminance, visibility and present-weather fields are 0. +- The SAM_OpenStudio `R7/WeatherImportTests` and `C4/EmbeddedWeatherDesignDaysTests` export through this writer + but assert no radiation fields. They were not run here, because they need the OpenStudio SDK. + +## Recommended next step +1. Owner approval and merge of PR #212. Then the `PROJECT_PROGRESS.md` closeout with the merge SHA, and cleanup + of the task branch. +2. F2 (EPW culture and missing-code correctness for the other fields and the importer), in the same session, on a + new branch from the merged `sow/2026-Q4`. +3. Later: re-run the PR #239 box and HungaryHouse embedded-weather comparisons with a licensed TAS, and record + the new embedded-vs-explicit residuals. Suggested settings: Opus 5.5, High effort, a new session in Build mode, + because it is a licensed validation run. diff --git a/documentation/evidence/epw-weather-f1/eplus-box-corrected-refIR.txt b/documentation/evidence/epw-weather-f1/eplus-box-corrected-refIR.txt new file mode 100644 index 000000000..73a3be0d2 --- /dev/null +++ b/documentation/evidence/epw-weather-f1/eplus-box-corrected-refIR.txt @@ -0,0 +1,11 @@ +environments=1 hours=61320 +Environment,Site Outdoor Air Dewpoint Temperature [C] | n=8760 | mean=5.636 | max=28.583 @h5030 +Environment,Site Sky Temperature [C] | n=8760 | mean=-1.872 | max=27.231 @h5030 +Environment,Site Horizontal Infrared Radiation Rate per Area [W/m2] | n=8760 | mean=310.477 | max=461.583 @h5030 +Environment,Site Diffuse Solar Radiation Rate per Area [W/m2] | n=8760 | mean=61.176 | max=454.167 @h4594 +Environment,Site Direct Solar Radiation Rate per Area [W/m2] | n=8760 | mean=162.071 | max=866.884 @h4689 +SAM_SPACE_OFFICE_CB68BED1,Enclosure Windows Total Transmitted Solar Radiation Energy [J] | n=8760 | sum=3755.1314 kWh | max=2268.92 W @h901 +SAM_IDEALLOADS_OFFICE_CB68BED1,Zone Ideal Loads Supply Air Sensible Heating Energy [J] | n=8760 | sum=2101.1184 kWh | max=3209.97 W @h8455 +SAM_IDEALLOADS_OFFICE_CB68BED1,Zone Ideal Loads Supply Air Total Heating Energy [J] | n=8760 | sum=2101.1184 kWh | max=3209.97 W @h8455 +SAM_IDEALLOADS_OFFICE_CB68BED1,Zone Ideal Loads Supply Air Sensible Cooling Energy [J] | n=8760 | sum=4238.4259 kWh | max=3544.42 W @h4381 +SAM_IDEALLOADS_OFFICE_CB68BED1,Zone Ideal Loads Supply Air Total Cooling Energy [J] | n=8760 | sum=4656.8346 kWh | max=5388.94 W @h5029 diff --git a/documentation/evidence/epw-weather-f1/eplus-box-corrected-refopaque.txt b/documentation/evidence/epw-weather-f1/eplus-box-corrected-refopaque.txt new file mode 100644 index 000000000..acb62aaae --- /dev/null +++ b/documentation/evidence/epw-weather-f1/eplus-box-corrected-refopaque.txt @@ -0,0 +1,11 @@ +environments=1 hours=61320 +Environment,Site Outdoor Air Dewpoint Temperature [C] | n=8760 | mean=5.636 | max=28.583 @h5030 +Environment,Site Sky Temperature [C] | n=8760 | mean=0.876 | max=28.062 @h5031 +Environment,Site Horizontal Infrared Radiation Rate per Area [W/m2] | n=8760 | mean=322.199 | max=466.710 @h5031 +Environment,Site Diffuse Solar Radiation Rate per Area [W/m2] | n=8760 | mean=61.176 | max=454.167 @h4594 +Environment,Site Direct Solar Radiation Rate per Area [W/m2] | n=8760 | mean=162.071 | max=866.884 @h4689 +SAM_SPACE_OFFICE_CB68BED1,Enclosure Windows Total Transmitted Solar Radiation Energy [J] | n=8760 | sum=3755.1314 kWh | max=2268.92 W @h901 +SAM_IDEALLOADS_OFFICE_CB68BED1,Zone Ideal Loads Supply Air Sensible Heating Energy [J] | n=8760 | sum=2066.8763 kWh | max=3141.50 W @h8455 +SAM_IDEALLOADS_OFFICE_CB68BED1,Zone Ideal Loads Supply Air Total Heating Energy [J] | n=8760 | sum=2066.8763 kWh | max=3141.50 W @h8455 +SAM_IDEALLOADS_OFFICE_CB68BED1,Zone Ideal Loads Supply Air Sensible Cooling Energy [J] | n=8760 | sum=4327.9049 kWh | max=3599.83 W @h3902 +SAM_IDEALLOADS_OFFICE_CB68BED1,Zone Ideal Loads Supply Air Total Cooling Energy [J] | n=8760 | sum=4748.0915 kWh | max=5390.59 W @h5029 diff --git a/documentation/evidence/epw-weather-f1/eplus-box-embedded-before.txt b/documentation/evidence/epw-weather-f1/eplus-box-embedded-before.txt new file mode 100644 index 000000000..4906fe959 --- /dev/null +++ b/documentation/evidence/epw-weather-f1/eplus-box-embedded-before.txt @@ -0,0 +1,11 @@ +environments=1 hours=61320 +Environment,Site Outdoor Air Dewpoint Temperature [C] | n=8760 | mean=5.636 | max=28.583 @h5030 +Environment,Site Sky Temperature [C] | n=8760 | mean=-4.834 | max=19.811 @h4382 +Environment,Site Horizontal Infrared Radiation Rate per Area [W/m2] | n=8760 | mean=296.384 | max=417.637 @h4382 +Environment,Site Diffuse Solar Radiation Rate per Area [W/m2] | n=8760 | mean=61.176 | max=454.167 @h4594 +Environment,Site Direct Solar Radiation Rate per Area [W/m2] | n=8760 | mean=83.616 | max=735.250 @h4690 +SAM_SPACE_OFFICE_CB68BED1,Enclosure Windows Total Transmitted Solar Radiation Energy [J] | n=8760 | sum=2592.6255 kWh | max=1531.90 W @h2316 +SAM_IDEALLOADS_OFFICE_CB68BED1,Zone Ideal Loads Supply Air Sensible Heating Energy [J] | n=8760 | sum=2810.0183 kWh | max=3343.73 W @h919 +SAM_IDEALLOADS_OFFICE_CB68BED1,Zone Ideal Loads Supply Air Total Heating Energy [J] | n=8760 | sum=2810.0183 kWh | max=3343.73 W @h919 +SAM_IDEALLOADS_OFFICE_CB68BED1,Zone Ideal Loads Supply Air Sensible Cooling Energy [J] | n=8760 | sum=3432.5972 kWh | max=3335.32 W @h4381 +SAM_IDEALLOADS_OFFICE_CB68BED1,Zone Ideal Loads Supply Air Total Cooling Energy [J] | n=8760 | sum=3832.5030 kWh | max=5206.16 W @h5029 diff --git a/documentation/evidence/epw-weather-f1/eplus-box-embedded-corrected.txt b/documentation/evidence/epw-weather-f1/eplus-box-embedded-corrected.txt new file mode 100644 index 000000000..5702b5d0a --- /dev/null +++ b/documentation/evidence/epw-weather-f1/eplus-box-embedded-corrected.txt @@ -0,0 +1,11 @@ +environments=1 hours=61320 +Environment,Site Outdoor Air Dewpoint Temperature [C] | n=8760 | mean=5.636 | max=28.583 @h5030 +Environment,Site Sky Temperature [C] | n=8760 | mean=0.275 | max=28.062 @h5031 +Environment,Site Horizontal Infrared Radiation Rate per Area [W/m2] | n=8760 | mean=319.601 | max=466.710 @h5031 +Environment,Site Diffuse Solar Radiation Rate per Area [W/m2] | n=8760 | mean=61.176 | max=454.167 @h4594 +Environment,Site Direct Solar Radiation Rate per Area [W/m2] | n=8760 | mean=162.071 | max=866.884 @h4689 +SAM_SPACE_OFFICE_CB68BED1,Enclosure Windows Total Transmitted Solar Radiation Energy [J] | n=8760 | sum=3755.1314 kWh | max=2268.92 W @h901 +SAM_IDEALLOADS_OFFICE_CB68BED1,Zone Ideal Loads Supply Air Sensible Heating Energy [J] | n=8760 | sum=2081.8983 kWh | max=3189.84 W @h8455 +SAM_IDEALLOADS_OFFICE_CB68BED1,Zone Ideal Loads Supply Air Total Heating Energy [J] | n=8760 | sum=2081.8983 kWh | max=3189.84 W @h8455 +SAM_IDEALLOADS_OFFICE_CB68BED1,Zone Ideal Loads Supply Air Sensible Cooling Energy [J] | n=8760 | sum=4311.7222 kWh | max=3604.29 W @h3902 +SAM_IDEALLOADS_OFFICE_CB68BED1,Zone Ideal Loads Supply Air Total Cooling Energy [J] | n=8760 | sum=4731.6249 kWh | max=5390.59 W @h5029 diff --git a/documentation/evidence/epw-weather-f1/eplus-box-explicit-Budapest.txt b/documentation/evidence/epw-weather-f1/eplus-box-explicit-Budapest.txt new file mode 100644 index 000000000..4b107f602 --- /dev/null +++ b/documentation/evidence/epw-weather-f1/eplus-box-explicit-Budapest.txt @@ -0,0 +1,11 @@ +environments=1 hours=61320 +Environment,Site Outdoor Air Dewpoint Temperature [C] | n=8760 | mean=5.636 | max=28.583 @h5030 +Environment,Site Sky Temperature [C] | n=8760 | mean=-1.872 | max=27.231 @h5030 +Environment,Site Horizontal Infrared Radiation Rate per Area [W/m2] | n=8760 | mean=310.477 | max=461.583 @h5030 +Environment,Site Diffuse Solar Radiation Rate per Area [W/m2] | n=8760 | mean=61.176 | max=454.167 @h4594 +Environment,Site Direct Solar Radiation Rate per Area [W/m2] | n=8760 | mean=148.809 | max=864.000 @h2507 +SAM_SPACE_OFFICE_CB68BED1,Enclosure Windows Total Transmitted Solar Radiation Energy [J] | n=8760 | sum=3623.6820 kWh | max=2184.43 W @h924 +SAM_IDEALLOADS_OFFICE_CB68BED1,Zone Ideal Loads Supply Air Sensible Heating Energy [J] | n=8760 | sum=2173.9185 kWh | max=3259.18 W @h8455 +SAM_IDEALLOADS_OFFICE_CB68BED1,Zone Ideal Loads Supply Air Total Heating Energy [J] | n=8760 | sum=2173.9185 kWh | max=3259.18 W @h8455 +SAM_IDEALLOADS_OFFICE_CB68BED1,Zone Ideal Loads Supply Air Sensible Cooling Energy [J] | n=8760 | sum=4140.7887 kWh | max=3552.50 W @h4381 +SAM_IDEALLOADS_OFFICE_CB68BED1,Zone Ideal Loads Supply Air Total Cooling Energy [J] | n=8760 | sum=4556.2262 kWh | max=5332.62 W @h5029 diff --git a/documentation/evidence/epw-weather-f1/rayleigh-bound-activation.txt b/documentation/evidence/epw-weather-f1/rayleigh-bound-activation.txt new file mode 100644 index 000000000..3497db294 --- /dev/null +++ b/documentation/evidence/epw-weather-f1/rayleigh-bound-activation.txt @@ -0,0 +1,39 @@ +HUN_BU_Budapest.Met.Center.128400_TMYx.epw.gz: bound binds in 170 of 4217 beam hours (4.03 %), all with hour-mean cos zenith <= 0.082 (sun below ~4.7 deg); largest unbounded DNI 471630 W/m2 + beam with the sun below the horizon all hour (DNI 0 by rule): 163 hours, 0.51 kWh/m2 + of the bound hours, 2 have the sun up for only about a minute (geometry not read back; excluded from the energy sums) + annual DNI: unbounded 3108.8, bounded 1495.2 kWh/m2 (-51.91 %) + annual horizontal beam: stated 734.6, exported 733.3 kWh/m2 (-0.18 %) + facade N annual beam: unbounded 41.0, bounded 31.5 kWh/m2 (-23.14 %) + facade E annual beam: unbounded 302.2, bounded 302.2 kWh/m2 (0.00 %) + facade S annual beam: unbounded 667.6, bounded 637.5 kWh/m2 (-4.50 %) + facade W annual beam: unbounded 567.2, bounded 479.0 kWh/m2 (-15.54 %) +SGP_SG_Changi.Intl.AP.486980_TMYx.2004-2018.epw.gz: bound binds in 0 of 3685 beam hours (0.00 %), all with hour-mean cos zenith <= 0.000 (sun below ~0.0 deg); largest unbounded DNI 0 W/m2 + beam with the sun below the horizon all hour (DNI 0 by rule): 0 hours, 0.00 kWh/m2 + of the bound hours, 0 have the sun up for only about a minute (geometry not read back; excluded from the energy sums) + annual DNI: unbounded 468.1, bounded 468.1 kWh/m2 (0.00 %) + annual horizontal beam: stated 360.5, exported 360.5 kWh/m2 (0.00 %) + facade N annual beam: unbounded 48.2, bounded 48.2 kWh/m2 (0.00 %) + facade E annual beam: unbounded 90.8, bounded 90.8 kWh/m2 (0.00 %) + facade S annual beam: unbounded 68.7, bounded 68.7 kWh/m2 (0.00 %) + facade W annual beam: unbounded 127.4, bounded 127.4 kWh/m2 (0.00 %) +AUS_NSW_Sydney.Intl.AP.947670_TMYx.2004-2018.epw.gz: bound binds in 0 of 3819 beam hours (0.00 %), all with hour-mean cos zenith <= 0.000 (sun below ~0.0 deg); largest unbounded DNI 0 W/m2 + annual DNI: unbounded 1482.6, bounded 1482.6 kWh/m2 (0.00 %) + annual horizontal beam: stated 951.8, exported 951.8 kWh/m2 (0.00 %) + facade N annual beam: unbounded 594.1, bounded 594.1 kWh/m2 (0.00 %) + facade E annual beam: unbounded 345.8, bounded 345.8 kWh/m2 (0.00 %) + facade S annual beam: unbounded 15.1, bounded 15.1 kWh/m2 (0.00 %) + facade W annual beam: unbounded 374.4, bounded 374.4 kWh/m2 (0.00 %) +USA_AK_Fairbanks.Intl.AP.702610_TMYx.2004-2018.epw.gz: bound binds in 0 of 2428 beam hours (0.00 %), all with hour-mean cos zenith <= 0.000 (sun below ~0.0 deg); largest unbounded DNI 0 W/m2 + annual DNI: unbounded 385.2, bounded 385.2 kWh/m2 (0.00 %) + annual horizontal beam: stated 184.7, exported 184.7 kWh/m2 (0.00 %) + facade N annual beam: unbounded 4.9, bounded 4.9 kWh/m2 (0.00 %) + facade E annual beam: unbounded 108.9, bounded 108.9 kWh/m2 (0.00 %) + facade S annual beam: unbounded 218.4, bounded 218.4 kWh/m2 (0.00 %) + facade W annual beam: unbounded 87.2, bounded 87.2 kWh/m2 (0.00 %) +USA_AZ_Phoenix-Sky.Harbor.Intl.AP.722780_TMYx.2004-2018.epw.gz: bound binds in 251 of 4217 beam hours (5.95 %), all with hour-mean cos zenith <= 0.353 (sun below ~20.7 deg); largest unbounded DNI 2450 W/m2 + annual DNI: unbounded 2911.3, bounded 2807.3 kWh/m2 (-3.57 %) + annual horizontal beam: stated 1620.1, exported 1613.9 kWh/m2 (-0.38 %) + facade N annual beam: unbounded 77.3, bounded 68.1 kWh/m2 (-11.93 %) + facade E annual beam: unbounded 535.9, bounded 535.9 kWh/m2 (0.00 %) + facade S annual beam: unbounded 987.1, bounded 973.4 kWh/m2 (-1.38 %) + facade W annual beam: unbounded 1078.9, bounded 1011.4 kWh/m2 (-6.25 %) diff --git a/documentation/evidence/epw-weather-f1/reference-epw-statistics.txt b/documentation/evidence/epw-weather-f1/reference-epw-statistics.txt new file mode 100644 index 000000000..ea6a5113b --- /dev/null +++ b/documentation/evidence/epw-weather-f1/reference-epw-statistics.txt @@ -0,0 +1,45 @@ +SGP_SG_Changi.Intl.AP.486980_TMYx.2004-2018.epw.gz: lat 1.35, lon 103.994, tz 8 + reference closure |DNI*cosZ - (GHI-DHI)|: 0.54 % (geometry 30 min earlier: 7.78 %) + annual DNI: reference 466.8, exported 468.1 kWh/m2 (+0.28 %) + annual horizontal beam: GHI-DHI 360.5, exported DNI*cosZ 360.4 kWh/m2 (-0.01 %) + sun-up hours n=631: median DNI error +0.27 %, 90th percentile |error| 0.57 % + facade N annual beam: reference 48.1, exported 48.2 kWh/m2 + facade E annual beam: reference 90.5, exported 90.8 kWh/m2 + facade S annual beam: reference 68.5, exported 68.7 kWh/m2 + facade W annual beam: reference 127.0, exported 127.4 kWh/m2 +USA_AK_Fairbanks.Intl.AP.702610_TMYx.2004-2018.epw.gz: lat 64.804, lon -147.876, tz -9 + reference closure |DNI*cosZ - (GHI-DHI)|: 1.32 % (geometry 30 min earlier: 6.07 %) + annual DNI: reference 384.3, exported 385.2 kWh/m2 (+0.23 %) + annual horizontal beam: GHI-DHI 184.7, exported DNI*cosZ 185.3 kWh/m2 (+0.30 %) + sun-up hours n=661: median DNI error +0.17 %, 90th percentile |error| 0.61 % + facade N annual beam: reference 4.9, exported 4.9 kWh/m2 + facade E annual beam: reference 108.6, exported 108.9 kWh/m2 + facade S annual beam: reference 217.7, exported 218.4 kWh/m2 + facade W annual beam: reference 87.1, exported 87.2 kWh/m2 +AUS_NSW_Sydney.Intl.AP.947670_TMYx.2004-2018.epw.gz: lat -33.946, lon 151.177, tz 10 + reference closure |DNI*cosZ - (GHI-DHI)|: 0.40 % (geometry 30 min earlier: 8.15 %) + annual DNI: reference 1479.4, exported 1482.6 kWh/m2 (+0.22 %) + annual horizontal beam: GHI-DHI 951.8, exported DNI*cosZ 952.0 kWh/m2 (+0.02 %) + sun-up hours n=2186: median DNI error +0.24 %, 90th percentile |error| 0.59 % + facade N annual beam: reference 592.5, exported 594.1 kWh/m2 + facade E annual beam: reference 345.1, exported 345.8 kWh/m2 + facade S annual beam: reference 15.1, exported 15.1 kWh/m2 + facade W annual beam: reference 373.6, exported 374.4 kWh/m2 +USA_AZ_Phoenix-Sky.Harbor.Intl.AP.722780_TMYx.2004-2018.epw.gz: lat 33.43333, lon -112.0167, tz -7 + reference closure |DNI*cosZ - (GHI-DHI)|: 8.85 % (geometry 30 min earlier: 2.34 %) + annual DNI: reference 2632.0, exported 2807.3 kWh/m2 (+6.66 %) + annual horizontal beam: GHI-DHI 1620.1, exported DNI*cosZ 1613.9 kWh/m2 (-0.38 %) + sun-up hours n=3243: median DNI error -1.49 %, 90th percentile |error| 19.96 % + facade N annual beam: reference 49.6, exported 68.1 kWh/m2 + facade E annual beam: reference 598.2, exported 535.9 kWh/m2 + facade S annual beam: reference 952.6, exported 973.4 kWh/m2 + facade W annual beam: reference 796.6, exported 1011.4 kWh/m2 +HUN_BU_Budapest.Met.Center.128400_TMYx.epw.gz: lat 47.517, lon 19.033, tz 1 + reference closure |DNI*cosZ - (GHI-DHI)|: 9.31 % (geometry 30 min earlier: 4.26 %) + annual DNI: reference 1304.6, exported 1496.2 kWh/m2 (+14.69 %) + annual horizontal beam: GHI-DHI 735.2, exported DNI*cosZ 733.3 kWh/m2 (-0.25 %) + sun-up hours n=2505: median DNI error -1.36 %, 90th percentile |error| 20.91 % + facade N annual beam: reference 16.5, exported 31.5 kWh/m2 + facade E annual beam: reference 324.7, exported 302.2 kWh/m2 + facade S annual beam: reference 613.0, exported 637.5 kWh/m2 + facade W annual beam: reference 345.1, exported 479.1 kWh/m2 diff --git a/documentation/evidence/epw-weather-f1/run-inputs.sha256 b/documentation/evidence/epw-weather-f1/run-inputs.sha256 new file mode 100644 index 000000000..c6bddfb25 --- /dev/null +++ b/documentation/evidence/epw-weather-f1/run-inputs.sha256 @@ -0,0 +1,10 @@ +c3e1054fe8bcb4f44844d4560c811bc5c209952c073f6e107f46a0a5db2f6e8f *weather/corrected-refIR.epw +50fc4ffdf7b28b0c34c6bb1f29759776f97b30b884ee625bd7b3c7f456099d46 *weather/corrected-refopaque.epw +e1353bd2600a401f76262d52a8e16176a15dc45a4967e58d9107c374704abf7c *weather/embedded-before.epw +5afbeb329bea4ce807f402f095fbc695f7fce8a2a29aa6c7b0d93a31a0063082 *weather/embedded-corrected.epw +472ffb4e2f10867600949a392ce46891fc1bd05f2fc8d8917a99bca7d26209c3 *weather/explicit-Budapest.epw +0331f55b13dd76cba401543272e29f690b66c9862c8ff05d73d6a4918d90f0ec *run-corrected-refIR/in.idf +0331f55b13dd76cba401543272e29f690b66c9862c8ff05d73d6a4918d90f0ec *run-corrected-refopaque/in.idf +0331f55b13dd76cba401543272e29f690b66c9862c8ff05d73d6a4918d90f0ec *run-embedded-before/in.idf +0331f55b13dd76cba401543272e29f690b66c9862c8ff05d73d6a4918d90f0ec *run-embedded-corrected/in.idf +0331f55b13dd76cba401543272e29f690b66c9862c8ff05d73d6a4918d90f0ec *run-explicit-Budapest/in.idf diff --git a/documentation/evidence/epw-weather-f1/weather-field-summary.txt b/documentation/evidence/epw-weather-f1/weather-field-summary.txt new file mode 100644 index 000000000..a34a043eb --- /dev/null +++ b/documentation/evidence/epw-weather-f1/weather-field-summary.txt @@ -0,0 +1,3 @@ +embedded-before n=8760 DBT 10.321 DPT 5.638 ETRH 0.0 kWh IR 0.0 kWh (missing 0 h) GHI 0.0 DNI 735.2 DHI 541.2 kWh/m2 sky tot 0.58 opq 0.00 tenths +embedded-corrected n=8760 DBT 10.321 DPT 5.638 ETRH 2601.4 kWh IR 0.0 kWh (missing 8760 h) GHI 1276.4 DNI 1496.2 DHI 541.2 kWh/m2 sky tot 5.80 opq 5.80 tenths +explicit-Budapest n=8760 DBT 10.321 DPT 5.638 ETRH 2598.6 kWh IR 2719.8 kWh (missing 0 h) GHI 1276.4 DNI 1304.6 DHI 541.2 kWh/m2 sky tot 5.80 opq 6.56 tenths