Nightshade w/ refraction=-6.0 yields ValueError: GEOSGeom_createLinearRing_r returned a NULL pointer

[captcha1]

Description

I'm trying use Nightshade to draw twilight lines (and sub solar point).

If I do Nightshade w/ refraction=-6.0, I get ValueError: GEOSGeom_createLinearRing_r returned a NULL pointer

Nightshade w/ refraction=-6.0001 works OK

Also, refraction=89 causes extraneous lines.

Code to reproduce

import matplotlib.pyplot as plt
fig = plt.figure()

import cartopy
ax = fig.add_subplot(projection=cartopy.crs.PlateCarree())
ax.coastlines()

import datetime
date = datetime.datetime(1999, 12, 31, 12)

from cartopy.feature.nightshade import Nightshade

ax.add_feature(Nightshade(date, color='none'))
ax.add_feature(Nightshade(date, refraction=-6.0, color='none')) # Civil Twilight
#ax.add_feature(Nightshade(date, refraction=-6.0001, color='none')) # Civil Twilight
ax.add_feature(Nightshade(date, refraction=-12, color='none')) # Nautical Twilight

ax.add_feature(Nightshade(date, refraction=89, color='none')) # sub solar circle ...

plt.show()

Traceback

python /tmp/nightshade-bug.py 
/usr/lib64/python3.10/site-packages/cartopy/feature/nightshade.py:80: RuntimeWarning: invalid value encountered in arccos
  omega0 = np.rad2deg(np.arccos(np.sin(np.deg2rad(refraction)) /
IllegalArgumentException: Points of LinearRing do not form a closed linestring
Traceback (most recent call last):
  File "/tmp/nightshade-bug.py", line 15, in <module>
    ax.add_feature(Nightshade(date, refraction=-6, color='none')) # Civil Twilight
  File "/usr/lib64/python3.10/site-packages/cartopy/feature/nightshade.py", line 93, in __init__
    geom = sgeom.Polygon(np.column_stack((x, y)))
  File "/usr/lib64/python3.10/site-packages/shapely/geometry/polygon.py", line 243, in __init__
    ret = geos_polygon_from_py(shell, holes)
  File "/usr/lib64/python3.10/site-packages/shapely/geometry/polygon.py", line 509, in geos_polygon_from_py
    ret = geos_linearring_from_py(shell)
  File "shapely/speedups/_speedups.pyx", line 408, in shapely.speedups._speedups.geos_linearring_from_py
ValueError: GEOSGeom_createLinearRing_r returned a NULL pointer

Full environment definition

Operating system

Fedora 35

Cartopy version

python3-cartopy-0.20.1-1.fc35.x86_64

python3-shapely-1.7.1-12.fc35.x86_64
python3-matplotlib-3.5.1-3.fc35.x86_64

Thanks !

[greglucas]

This was an unfortunate floating point issue in the math underlying the Nightshade calculations. The -6 value should be fixed by #1989.

I think what is likely going on with the subsolar point extra lines is that there is some more floating point imprecision going on in the opposite direction of the other issue. I don't think "refraction" is a good use for that, as the light curve corrections probably aren't what you actually want? I assume you want an actual circle and not a curve, so I'd suggest using a patch or single scatter point instead, which would avoid all of the extra curved calculations.