How to use the cartopy.feature.OCEAN function in Cartopy

# ---------------------------------Overview Location Map ------------------------
#
# set up index map 20% height, left 16% of figure
left = 0.03
bottom = 0
width = 0.16
height = 0.2
rect = [left, bottom, width, height]

ax2 = plt.axes(rect, projection=ccrs.PlateCarree(), )
ax2.set_extent((110, 160, -45, -10))
#  ax2.set_global()  will show the whole world as context

ax2.coastlines(resolution='110m', zorder=2)
ax2.add_feature(cfeature.LAND)
ax2.add_feature(cfeature.OCEAN)

ax2.gridlines()

lon0, lon1, lat0, lat1 = ax.get_extent()
box_x = [lon0, lon1, lon1, lon0, lon0]
box_y = [lat0, lat0, lat1, lat1, lat0]

plt.plot(box_x, box_y, color='red', transform=ccrs.Geodetic())

# -------------------------------- Title -----------------------------
# set up map title top 4% of figure, right 80% of figure

left = 0.2
bottom = 0.95
width = 0.8
height = 0.04

def main():
    fig = plt.figure()
    ax = fig.add_subplot(1, 1, 1, projection=ccrs.Orthographic(-10, 45))

    ax.add_feature(cfeature.OCEAN, zorder=0)
    ax.add_feature(cfeature.LAND, zorder=0, edgecolor='black')

    ax.set_global()
    ax.gridlines()

    x, y, u, v, vector_crs = sample_data()
    ax.quiver(x, y, u, v, transform=vector_crs)

    plt.show()

def main():
    fig = plt.figure()
    ax = fig.add_subplot(1, 1, 1, projection=ccrs.PlateCarree())
    ax.set_extent([-20, 60, -40, 45], crs=ccrs.PlateCarree())

    ax.add_feature(cfeature.LAND)
    ax.add_feature(cfeature.OCEAN)
    ax.add_feature(cfeature.COASTLINE)
    ax.add_feature(cfeature.BORDERS, linestyle=':')
    ax.add_feature(cfeature.LAKES, alpha=0.5)
    ax.add_feature(cfeature.RIVERS)

    plt.show()

def get_ax(self, data):
        """Create and populate new axes object."""
        timezone = list([x for x in data if 'UTC' in x])

        timezone_start = tuple((x/255 for x in (0, 255, 0, 100)))
        country_start = tuple((x/255 for x in (0, 100, 0)))
        # We ignore some countries, as they are too big and need a higher
        # resolution for precise timezone assignment.
        ignored_countries = ['United States', 'Australia', 'Brazil', 'Canada']

        ax = plt.axes(projection=ccrs.PlateCarree())

        # Print countries and state borders
        ax.add_feature(cartopy.feature.LAND)
        ax.add_feature(cartopy.feature.OCEAN)
        ax.add_feature(cartopy.feature.COASTLINE)
        ax.add_feature(cartopy.feature.BORDERS)
        for state in self.states:
            ax.add_geometries(
                state.geometry,
                ccrs.PlateCarree(),
                facecolor=np.array((240, 240, 220)) / 256,
                edgecolor='black',
                label=state.attributes['name'],
            )

        collected_countries = []
        collected_timezones = []
        collected_states = []

        timezones_to_draw = []

custom_layout.add_barb('eastward_wind', 'northward_wind', units='knots')
custom_layout.add_value('NW', 'air_temperature', fmt='.1f', units='degF', color='darkred')
custom_layout.add_value('SW', 'dew_point_temperature', fmt='.1f', units='degF',
                        color='darkgreen')

# Also, we'll add a field that we don't have in our dataset. This will be ignored
custom_layout.add_value('E', 'precipitation', fmt='0.2f', units='inch', color='blue')

# Create the figure and an axes set to the projection
fig = plt.figure(figsize=(20, 10))
add_metpy_logo(fig, 1080, 290, size='large')
ax = fig.add_subplot(1, 1, 1, projection=proj)

# Add some various map elements to the plot to make it recognizable
ax.add_feature(cfeature.LAND)
ax.add_feature(cfeature.OCEAN)
ax.add_feature(cfeature.LAKES)
ax.add_feature(cfeature.COASTLINE)
ax.add_feature(cfeature.STATES)
ax.add_feature(cfeature.BORDERS, linewidth=2)

# Set plot bounds
ax.set_extent((-118, -73, 23, 50))

#
# Here's the actual station plot
#

# Start the station plot by specifying the axes to draw on, as well as the
# lon/lat of the stations (with transform). We also the fontsize to 12 pt.
stationplot = StationPlot(ax, data['longitude'], data['latitude'],
                          transform=ccrs.PlateCarree(), fontsize=12)

def navtimeseries(nav: xarray.Dataset):
    if not isinstance(nav, xarray.Dataset):
        return

    svs = nav.sv.values

    if cartopy is not None:
        ax = figure().gca(projection=cartopy.crs.PlateCarree())

        ax.add_feature(cpf.LAND)
        ax.add_feature(cpf.OCEAN)
        ax.add_feature(cpf.COASTLINE)
        ax.add_feature(cpf.BORDERS, linestyle=':')
    else:
        ax = figure().gca()

    for sv in svs:
        if sv[0] == 'S':
            lat, lon, alt = pm.ecef2geodetic(nav.sel(sv=sv)['X'].dropna(dim='time', how='all'),
                                             nav.sel(sv=sv)['Y'].dropna(
                dim='time', how='all'),
                nav.sel(sv=sv)['Z'].dropna(dim='time', how='all'))

            if ((alt < 35.7e6) | (alt > 35.9e6)).any():
                logging.warning('unrealistic geostationary satellite altitudes')

            if ((lat < -1) | (lat > 1)).any():

###########################################
# The payoff
# ----------

# Change the DPI of the resulting figure. Higher DPI drastically improves the
# look of the text rendering.
plt.rcParams['savefig.dpi'] = 255

# Create the figure and an axes set to the projection.
fig = plt.figure(figsize=(20, 10))
add_metpy_logo(fig, 1080, 290, size='large')
ax = fig.add_subplot(1, 1, 1, projection=proj)

# Add some various map elements to the plot to make it recognizable.
ax.add_feature(cfeature.LAND)
ax.add_feature(cfeature.OCEAN)
ax.add_feature(cfeature.LAKES)
ax.add_feature(cfeature.COASTLINE)
ax.add_feature(cfeature.STATES)
ax.add_feature(cfeature.BORDERS)

# Set plot bounds
ax.set_extent((-118, -73, 23, 50))

#
# Here's the actual station plot
#

# Start the station plot by specifying the axes to draw on, as well as the
# lon/lat of the stations (with transform). We also the fontsize to 12 pt.
stationplot = StationPlot(ax, data['lon'].values, data['lat'].values, clip_on=True,
                          transform=ccrs.PlateCarree(), fontsize=12)

"""
    if ax is None and kwargs is not None and crs is None:
        if 'subplot_kw' not in kwargs:
            kwargs['subplot_kw'] = {'projection': ccrs.PlateCarree()}
        elif 'subplot_kw' in kwargs and 'projection' not in kwargs['subplot_kw']:
            kwargs['subplot_kw']['projection'] = ccrs.PlateCarree()
        f, ax = plt.subplots(figsize=(10, 5), **kwargs)
    elif ax is None and crs is not None:
        f, ax = plt.subplots(figsize=(10, 5), subplot_kw={'projection': crs})
    else:
        f, ax = plt.subplots(
            figsize=(10, 5), subplot_kw={'projection': ccrs.PlateCarree()})
    if natural_earth:
        #~ ax.stock_img()
        ax.add_feature(cfeature.OCEAN)
        ax.add_feature(cfeature.LAND)
        ax.add_feature(cfeature.LAKES)
        ax.add_feature(cfeature.RIVERS)

    if states:
        states_provinces = cfeature.NaturalEarthFeature(
            category='cultural',
            name='admin_1_states_provinces_lines',
            scale=resolution,
            facecolor='none')
        ax.add_feature(states_provinces, edgecolor='black')

    if coastlines:
        ax.coastlines(resolution)

    if countries:

# * Departure > 2 is Magenta
#
# This should help us visualize where the precipitation event was strongest!

proj = ccrs.LambertConformal(central_longitude=-105, central_latitude=0,
                             standard_parallels=[35])

fig = plt.figure(figsize=(20, 10))
ax = fig.add_subplot(1, 1, 1, projection=proj)

state_boundaries = feat.NaturalEarthFeature(category='cultural',
                                            name='admin_1_states_provinces_lines',
                                            scale='110m', facecolor='none')

ax.add_feature(feat.LAND, zorder=-1)
ax.add_feature(feat.OCEAN, zorder=-1)
ax.add_feature(feat.LAKES, zorder=-1)
ax.coastlines(resolution='110m', zorder=2, color='black')
ax.add_feature(state_boundaries, edgecolor='black')
ax.add_feature(feat.BORDERS, linewidth=2, edgecolor='black')

# Set plot bounds
ax.set_extent((-109.9, -101.8, 36.5, 41.3))

# Plot each station, labeling based on departure
for stn in range(len(pcpn)):
    if pcpn_dep[stn] >= 0 and pcpn_dep[stn] < 2:
        ax.plot(lon[stn], lat[stn], 'g+', markersize=7, transform=ccrs.Geodetic())
        ax.text(lon[stn], lat[stn], pcpn[stn], transform=ccrs.Geodetic())
    elif pcpn_dep[stn] >= 2:
        ax.plot(lon[stn], lat[stn], 'm+', markersize=7, transform=ccrs.Geodetic())
        ax.text(lon[stn], lat[stn], pcpn[stn], transform=ccrs.Geodetic())

ax2 = fig.add_subplot(1, 2, 2, projection=ccrs.SouthPolarStereo(),
                          sharex=ax1, sharey=ax1)
    fig.subplots_adjust(bottom=0.05, top=0.95,
                        left=0.04, right=0.95, wspace=0.02)

    # Limit the map to -60 degrees latitude and below.
    ax1.set_extent([-180, 180, -90, -60], ccrs.PlateCarree())

    ax1.add_feature(cfeature.LAND)
    ax1.add_feature(cfeature.OCEAN)

    ax1.gridlines()
    ax2.gridlines()

    ax2.add_feature(cfeature.LAND)
    ax2.add_feature(cfeature.OCEAN)

    # Compute a circle in axes coordinates, which we can use as a boundary
    # for the map. We can pan/zoom as much as we like - the boundary will be
    # permanently circular.
    theta = np.linspace(0, 2*np.pi, 100)
    center, radius = [0.5, 0.5], 0.5
    verts = np.vstack([np.sin(theta), np.cos(theta)]).T
    circle = mpath.Path(verts * radius + center)

    ax2.set_boundary(circle, transform=ax2.transAxes)

    plt.show()