diff --git a/docs/specification.md b/docs/specification.md index 95951ae..fa6846c 100644 --- a/docs/specification.md +++ b/docs/specification.md @@ -140,15 +140,15 @@ cell spacing `sqrt(4πR² / (12 · 4^order))`. The table below is regenerated from these formulas and pinned by `mortie/tests/test_spec_page.py` so it cannot drift. -**Note — code vs. page.** These are the **normative, sphere-derived** -values. `mortie.tools.order2res` in *code* retains the historical flat -constant `111 km/deg × 58.6323 × 0.5^order` (an implied sphere `R ≈ 6366 -km`) for **behavioral compatibility** — it is consumed by `res2display` and -by buffer-pad computations (e.g. `tests/test_coverage_boundary.py`, which -scales the coverage pad by `order2res`) whose outputs would shift if the -constant changed. The page's cell-scale column is therefore ~0.2% larger -than `order2res` at every order. Unifying `order2res` onto this sphere is a -behavioral change tracked separately in +**Note — code and page unified.** These are the **normative, sphere-derived** +values, and `mortie.tools.order2res` now derives from the same sphere: +`order2res(order) = sqrt(4πR² / (12 · 4^order))` with the single +`mortie.tools.EARTH_RADIUS_KM = 6371.0088` constant. Its consumers +(`res2display` and the buffer-pad computation in +`tests/test_coverage_boundary.py`) therefore read the cell-scale column +below directly. This replaced the historical flat constant +`111 km/deg × 58.6323 × 0.5^order` (an implied sphere `R ≈ 6366 km`), a +behavioral change of ~0.2% at every order, per [mortie #119](https://github.com/espg/mortie/issues/119). diff --git a/mortie/tests/test_spec_page.py b/mortie/tests/test_spec_page.py index f6bad51..affff5f 100644 --- a/mortie/tests/test_spec_page.py +++ b/mortie/tests/test_spec_page.py @@ -5,10 +5,9 @@ **Both** columns derive from the exact HEALPix sphere at mean radius ``EARTH_RADIUS_KM``: every order-k cell has identical area ``4*pi*R**2 / (12 * 4**k)``, and the cell scale is the square root of that -area (RMS cell spacing). This is the sphere-derived normative value, not the -historical ``order2res`` constant kept in ``mortie.tools`` for behavioral -compatibility (spec §3 code-vs-page note; unification tracked as a mortie -follow-up issue). This test compares every regenerated row literally against +area (RMS cell spacing). ``mortie.tools.order2res`` is derived from this same +``EARTH_RADIUS_KM`` sphere (issue #119), so the code and page now share one +Earth model. This test compares every regenerated row literally against the rows between the ``table:order2res`` markers; to refresh the doc after a deliberate formula change, paste the rows this module's ``table_rows()`` produces. @@ -17,13 +16,13 @@ import math from pathlib import Path -from mortie.tools import MAX_ORDER +from mortie.tools import EARTH_RADIUS_KM, MAX_ORDER SPEC_PAGE = Path(__file__).resolve().parents[2] / "docs" / "specification.md" BEGIN = "" END = "" # Exact HEALPix sphere, mean Earth radius (km); drives BOTH table columns. -EARTH_RADIUS_KM = 6371.0088 +# order2res is now unified onto this same sphere (issue #119). def format_row(order): diff --git a/mortie/tests/test_tools.py b/mortie/tests/test_tools.py index 2f1f0aa..991c1a7 100644 --- a/mortie/tests/test_tools.py +++ b/mortie/tests/test_tools.py @@ -11,6 +11,8 @@ - Tests focus on consistency, determinism, and structural validation """ +import math + import pytest import numpy as np from numpy.testing import assert_array_equal, assert_allclose @@ -18,25 +20,30 @@ from mortie import tools +def _sphere_res(order): + """Reference RMS cell spacing (km): sqrt of the equal-area cell area.""" + R = tools.EARTH_RADIUS_KM + return math.sqrt(4 * math.pi * R**2 / (12 * 4**order)) + + class TestOrder2Res: """Test order to resolution conversion""" def test_order2res_basic(self): """Test basic order to resolution calculations""" - # Order 0 should be largest resolution + # Order 0 is the RMS cell spacing on the unified HEALPix sphere. res0 = tools.order2res(0) - assert_allclose(res0, 111 * 58.6323, rtol=1e-10) + assert_allclose(res0, _sphere_res(0), rtol=1e-10) - # Order 1 should be half of order 0 + # Each order halves the cell scale (area drops by 4). res1 = tools.order2res(1) assert_allclose(res1, res0 / 2.0, rtol=1e-10) def test_order2res_range(self): - """Test full range of valid orders""" + """Test full range of valid orders against the sphere formula""" for order in range(20): res = tools.order2res(order) - expected = 111 * 58.6323 * (0.5 ** order) - assert_allclose(res, expected, rtol=1e-10) + assert_allclose(res, _sphere_res(order), rtol=1e-10) def test_order2res_decreasing(self): """Test that resolution decreases with order""" @@ -72,12 +79,12 @@ def test_unit_ladder_km_m_cm(self, capsys): tools.res2display() lines = capsys.readouterr().out.strip().split('\n') by_order = {int(line.rsplit(' ', 1)[1]): line for line in lines} - # order 12 = 1.589 km, order 13 = 794.456 m (the issue's examples) - assert by_order[12] == '1.589 km at tessellation order 12' - assert by_order[13] == '794.456 m at tessellation order 13' + # order 12 = 1.592 km, order 13 = 795.852 m (unified sphere, issue #119) + assert by_order[12] == '1.592 km at tessellation order 12' + assert by_order[13] == '795.852 m at tessellation order 13' # finest orders drop to cm rather than tiny km/m fractions - assert by_order[25] == '19.396 cm at tessellation order 25' - assert by_order[29] == '1.212 cm at tessellation order 29' + assert by_order[25] == '19.43 cm at tessellation order 25' + assert by_order[29] == '1.214 cm at tessellation order 29' def test_rounds_within_bracket(self, capsys): """Values are rounded to three decimals inside the chosen unit""" diff --git a/mortie/tools.py b/mortie/tools.py index 3c80971..b847caa 100644 --- a/mortie/tools.py +++ b/mortie/tools.py @@ -7,6 +7,11 @@ from . import _healpix as hp from . import _rustie +# Mean Earth radius (km), the exact HEALPix sphere the spec page's resolution +# table (docs/specification.md §3) derives from. order2res is the RMS cell +# spacing on this sphere; unified with the page in issue #119. +EARTH_RADIUS_KM = 6371.0088 + # Rust-native geo2mort (uses healpix crate, no Python HEALPix backend) _rust_geo2mort = _rustie.rust_geo2mort # Packed-word kernel bridge: morton <-> HEALPix NESTED (vectorized). @@ -20,8 +25,16 @@ def order2res(order): - res = 111 * 58.6323 * .5**order - return res + """Approximate cell scale (km) at a HEALPix tessellation ``order``. + + The exact RMS cell spacing on the mean-radius HEALPix sphere: every + order-k cell has identical area ``4*pi*R**2 / (12 * 4**order)`` (HEALPix is + equal-area), and the cell scale is the square root of that area. Derived + from :data:`EARTH_RADIUS_KM` so code and the spec page (§3) share one Earth + model (issue #119). + """ + area = 4 * np.pi * EARTH_RADIUS_KM**2 / (12 * 4**order) # km2 + return float(np.sqrt(area)) def res2display(max_order=MAX_ORDER): @@ -32,7 +45,7 @@ def res2display(max_order=MAX_ORDER): Each resolution is rendered in the largest sensible unit -- km at coarse orders, m once it drops below 1 km, cm once it drops below 1 m -- and rounded to three decimals within that bracket, so fine orders read - naturally (e.g. order 12 -> ``1.589 km``, order 13 -> ``794.456 m``) + naturally (e.g. order 12 -> ``1.592 km``, order 13 -> ``795.852 m``) rather than as tiny km fractions. ``max_order`` must lie in 0..MAX_ORDER, the order range the packed-u64