Clipping issue with Geostationary/NearsidePerspective

[dopplershift]

I received this report from a user who was concerned by the points outside the circle of the globe:

import numpy as np
import cartopy.crs as ccrs
import cartopy.feature as cfeature
import matplotlib.pyplot as plt

lons, lats = np.meshgrid(np.linspace(0, 360, 100), np.linspace(-90, 90, 100))

data_projection = ccrs.Geodetic()
map_projection1 = ccrs.NearsidePerspective(central_longitude=0.0, central_latitude=0.0,     
                                           satellite_height=700000)

map_projection2 = ccrs.Geostationary(central_longitude=0.0, satellite_height=700000)

map_projections = [map_projection1, map_projection2]

for map_projection in map_projections:
    fig = plt.figure(figsize=(12, 8))
    ax = fig.add_subplot(1, 1, 1, projection=map_projection)

    # add some common geographic features
    ax.coastlines(resolution='10m', color='black')
    ax.add_feature(cfeature.STATES, edgecolor='black')
    ax.add_feature(cfeature.BORDERS)

    # add some lat/lon gridlines
    ax.gridlines()
    ax.set_global()
    ax.scatter(lons, lats, transform=data_projection)

produces these images for me.
NearsidePerspective:

Geostationary:

I'm a little confused at the moment as to why these are so different as well.

[QuLogic]

Probably related to #1081.

[dopplershift]

Nope, the problem in #1081 was fixed by matplotlib >=2.2.3. The problem here persists for me with current master and latest matplotlib release (currently 3.0.2).

[dopplershift]

So there are a few things going on here. First, Geostationary != NearSidePerspective@equator. Geostationary uses h*scan_angle for the projection coordinates, while Near Side Perspective (NSP) uses x/y distances on a plane tangent to the Earth right below the perspective point.

This implies that the shape of the boundary for NSP should be a circle, but we're using a slightly incorrect formula (PR incoming). I should also note that Proj's nsper only works on spheres, so I'll also warn about that.

For Geostationary, the shape of the boundary is more complicated than a circle. This is why we see points outside our boundary; our boundary should extend out further--more of a heavily rounded square I think. For the life of me, though, I can't manage to generate a sane formula, probably because it involves some horrific math on a sphere. In practice, it's not really a big deal, since the only data that actually comes in that projection is at a height of 35e6 m, at which point a circle is a fine approximation.