Discontinuous RotatedPole parameters produce inverted coordinate system

reinderien commented 3 months ago

Description

From the beginning, the implementation of Nightshade has included

        if lat > 0:
            pole_lat = -90 + lat
            central_lon = 180
        else:
            pole_lat = 90 + lat
            central_lon = 0

I have code that depends on that routine and the resulting RotatedPole to do a series of ecliptic calculations. Unfortunately, I've had to copy-and-modify, because the original code - due to the if above - produces discontinuous behaviour when lat (i.e. solar declination) crosses zero. In one half of the annual domain, the behaviour is correct; and in the other half the following happen:

I believe the boundary of behaviour is at the September equinox
delta_sun changes sign, which in itself is not a problem
lat becomes positive; also not a problem
the if expression above abruptly starts evaluating to True
the entire coordinate system abruptly becomes inverted.

I have found that only the second case is correct (i.e. negative latitude). Copying the implementation and removing the conditional, preserving only the second case, produces correct, continuous behaviour.

I suspect that the reason this hasn't been noticed is that the nightshade curve is symmetric and so inversion is "invisible"; but anyone using the CRS owned by Nightshade for asymmetric curves is going to get a nasty surprise.

Code to reproduce

Compare behaviour before and after 3b77fe0, before and after an equinox.

Operating system

All

Cartopy version

0.22, but theoretically all versions

pip list

pip freeze
Cartopy==0.22.0
certifi==2024.2.2
charset-normalizer==3.3.2
contourpy==1.2.0
cycler==0.12.1
fonttools==4.49.0
idna==3.6
kiwisolver==1.4.5
matplotlib==3.8.3
mypy==1.8.0
mypy-extensions==1.0.0
numpy==1.26.4
packaging==23.2
pillow==10.2.0
pyparsing==3.1.1
pyproj==3.6.1
pyshp==2.3.1
python-dateutil==2.8.2
pytz==2024.1
requests==2.31.0
scipy==1.12.0
shapely==2.0.3
six==1.16.0
tomli==2.0.1
types-requests==2.31.0.20240311
typing_extensions==4.9.0
tzdata==2024.1
urllib3==2.2.1

lgolston commented 3 months ago

Can you post a screenshot of a specific day/time with incorrect behavior? As far as I understand from the description, this inversion twice a year during the equinox is working as it should, but could not follow the argument about asymmetric curves.

reinderien commented 3 months ago

@lgolston Thank you for your question. Here are the screenshots - Before the coordinate system fixup, before the equinox:

python -m isochrones -s -t 2024-09-15T12:00 -x -80 -y 40

After the equinox:

python -m isochrones -s -t 2024-09-30T12:00 -x -80 -y 40

In this 'isochrones' project it only makes sense for the isochrones from west to east to appear as orange, yellow, pink, purple and blue in that order. After the equinox they appear correct. Before the equinox, the coordinate system transformation essentially shows solar time in reverse.

greglucas commented 3 months ago

I am also not following exactly what is wrong/what is needed to fix here. Could you make a reproducible example with just Cartopy and the nightshade feature geometry/points explaining what you mean?

Is it that the geometry goes from something like [p0, p1, p2, ...] ascending points before September, and then goes to [... p2, p1, p0] after September and thus the coordinate order is reversed and that causes issues downstream? So it would help to always keep things in the same ascending/descending order within the geometry that is created.

SciTools / cartopy