isSunlit Function Potential Update

[pjimeco]

When comparing isSunlit Transition times between shadow/sunlight I found that when using the Skyfield library it assumes a sphere of earth_m radius instead of the showing the earth as an elipsoid. (The Intersect line and sphere function is used).

I created a function to collect transition times between sunlight/shadow and compared this to orekit transition times with and without flattening baked in.

Ii found that without flattening baked in (spherical) the Skyfield to Orekit transition time comparisons are pretty good <1 second difference.

However once you account for earth flattening in Orekit, the transition times in Skyfield seem to be accurate to only within 5-8 seconds.

Has there been any discussions on implementing a non-spherical assumption for the earth within the IsSunlit Function? This would help get a more accurate value of sunlit/not-sunlit when there is a 1 second resolution required.