From 5dbede1a4b284e150fb93a680db5608b4e1972d8 Mon Sep 17 00:00:00 2001 From: Don Cross Date: Tue, 28 May 2024 17:10:07 -0400 Subject: [PATCH] Looking at _CalcVsop as the possible origin of numeric errors. --- demo/python/astronomy.py | 5 +++++ generate/.gitignore | 1 + generate/commit_hook | 2 +- generate/commit_hook.bat | 3 ++- generate/template/astronomy.py | 5 +++++ source/python/README.md | 9 +++++++++ source/python/astronomy/astronomy.py | 5 +++++ 7 files changed, 28 insertions(+), 2 deletions(-) diff --git a/demo/python/astronomy.py b/demo/python/astronomy.py index afaca19a..ef5b859f 100644 --- a/demo/python/astronomy.py +++ b/demo/python/astronomy.py @@ -36,6 +36,7 @@ import enum import re import abc +from os import getenv from typing import Any, List, Tuple, Optional, Union, Callable, Dict def _cbrt(x: float) -> float: @@ -3180,12 +3181,16 @@ def _VsopSphereToRect(lon: float, lat: float, rad: float) -> _TerseVector: rad * math.sin(lat) ) +_ProblemDebug:bool = ('1' == getenv('ASTRONOMY_ENGINE_PYTHON_PROBLEM')) + def _CalcVsop(model: _vsop_model_t, time: Time) -> Vector: t = time.tt / _DAYS_PER_MILLENNIUM lon = _VsopFormula(model.lon, t, True) lat = _VsopFormula(model.lat, t, False) rad = _VsopFormula(model.rad, t, False) eclip = _VsopSphereToRect(lon, lat, rad) + if _ProblemDebug: + print('_CalcVsop: lon={:0.16g}, lat={:0.16g}, rad={:0.16g}, eclip=({:0.16g}, {:0.16g}, {:0.16g})'.format(lon, lat, rad, eclip.x, eclip.y, eclip.z)) return _VsopRotate(eclip).ToAstroVector(time) class _body_state_t: diff --git a/generate/.gitignore b/generate/.gitignore index 4077993a..503a4f94 100644 --- a/generate/.gitignore +++ b/generate/.gitignore @@ -13,3 +13,4 @@ obj/ profile/ private/ node_modules/ +precision_problem.txt diff --git a/generate/commit_hook b/generate/commit_hook index 666f21f8..5477cc52 100755 --- a/generate/commit_hook +++ b/generate/commit_hook @@ -1,6 +1,6 @@ #/bin/bash [[ -d generate ]] && cd generate -python3 test.py precision || exit $? +ASTRONOMY_ENGINE_PYTHON_PROBLEM=1 python3 test.py precision || exit $? rm -f output/vsop*.txt output/*.eph output/jupiter_moons.txt ./run || exit $? ./verify_clean || exit $? diff --git a/generate/commit_hook.bat b/generate/commit_hook.bat index db8d8eb0..1c8e4300 100644 --- a/generate/commit_hook.bat +++ b/generate/commit_hook.bat @@ -5,9 +5,10 @@ REM Astronomy Engine - GitHub Actions steps for Windows. REM This batch file is executed on every push to GitHub. REM -------------------------------------------------------------------------------- -REM **** HACK **** Just test the Python precision problem and quit. cd %~dp0 +set ASTRONOMY_ENGINE_PYTHON_PROBLEM=1 py test.py precision || exit /b 1 +set ASTRONOMY_ENGINE_PYTHON_PROBLEM= REM Change to project/repo root directory. cd %~dp0\.. diff --git a/generate/template/astronomy.py b/generate/template/astronomy.py index 5d995118..27d4a656 100644 --- a/generate/template/astronomy.py +++ b/generate/template/astronomy.py @@ -36,6 +36,7 @@ import enum import re import abc +from os import getenv from typing import Any, List, Tuple, Optional, Union, Callable, Dict def _cbrt(x: float) -> float: @@ -1865,12 +1866,16 @@ def _VsopSphereToRect(lon: float, lat: float, rad: float) -> _TerseVector: rad * math.sin(lat) ) +_ProblemDebug:bool = ('1' == getenv('ASTRONOMY_ENGINE_PYTHON_PROBLEM')) + def _CalcVsop(model: _vsop_model_t, time: Time) -> Vector: t = time.tt / _DAYS_PER_MILLENNIUM lon = _VsopFormula(model.lon, t, True) lat = _VsopFormula(model.lat, t, False) rad = _VsopFormula(model.rad, t, False) eclip = _VsopSphereToRect(lon, lat, rad) + if _ProblemDebug: + print('_CalcVsop: lon={:0.16g}, lat={:0.16g}, rad={:0.16g}, eclip=({:0.16g}, {:0.16g}, {:0.16g})'.format(lon, lat, rad, eclip.x, eclip.y, eclip.z)) return _VsopRotate(eclip).ToAstroVector(time) class _body_state_t: diff --git a/source/python/README.md b/source/python/README.md index 277cfdd4..5a7e928c 100644 --- a/source/python/README.md +++ b/source/python/README.md @@ -3567,3 +3567,12 @@ latitude, longitude, and elevation for that vector. The geographic latitude, longitude, and elevation above sea level that corresponds to the given equatorial vector. +--- + + +### getenv(key, default=None) + +Get an environment variable, return None if it doesn't exist. +The optional second argument can specify an alternate default. +key, default and the result are str. + diff --git a/source/python/astronomy/astronomy.py b/source/python/astronomy/astronomy.py index afaca19a..ef5b859f 100644 --- a/source/python/astronomy/astronomy.py +++ b/source/python/astronomy/astronomy.py @@ -36,6 +36,7 @@ import enum import re import abc +from os import getenv from typing import Any, List, Tuple, Optional, Union, Callable, Dict def _cbrt(x: float) -> float: @@ -3180,12 +3181,16 @@ def _VsopSphereToRect(lon: float, lat: float, rad: float) -> _TerseVector: rad * math.sin(lat) ) +_ProblemDebug:bool = ('1' == getenv('ASTRONOMY_ENGINE_PYTHON_PROBLEM')) + def _CalcVsop(model: _vsop_model_t, time: Time) -> Vector: t = time.tt / _DAYS_PER_MILLENNIUM lon = _VsopFormula(model.lon, t, True) lat = _VsopFormula(model.lat, t, False) rad = _VsopFormula(model.rad, t, False) eclip = _VsopSphereToRect(lon, lat, rad) + if _ProblemDebug: + print('_CalcVsop: lon={:0.16g}, lat={:0.16g}, rad={:0.16g}, eclip=({:0.16g}, {:0.16g}, {:0.16g})'.format(lon, lat, rad, eclip.x, eclip.y, eclip.z)) return _VsopRotate(eclip).ToAstroVector(time) class _body_state_t: