How to use the cartopy.crs.PlateCarree function in Cartopy

# Crater catalogue.
        self.craters = igen.ReadLROCHeadCombinedCraterCSV(
            filelroc="../catalogues/LROCCraters.csv",
            filehead="../catalogues/HeadCraters.csv",
            sortlat=True)

        # Long/lat limits
        self.cdim = [-180., 180., -60., 60.]

        # Coordinate systems.
        self.iglobe = ccrs.Globe(semimajor_axis=1737400,
                                 semiminor_axis=1737400,
                                 ellipse=None)
        self.geoproj = ccrs.Geodetic(globe=self.iglobe)
        self.iproj = ccrs.PlateCarree(globe=self.iglobe)

region = vd.get_region((data.longitude, data.latitude))
# The desired grid spacing in degrees (converted to meters using 1 degree approx. 111km)
spacing = 10 / 60
# Use Mercator projection because Spline is a Cartesian gridder
projection = pyproj.Proj(proj="merc", lat_ts=data.latitude.mean())
proj_coords = projection(data.longitude.values, data.latitude.values)

plt.figure(figsize=(7, 6))
ax = plt.axes(projection=ccrs.Mercator())
ax.set_title("Bathymetry from Baja California")
plt.scatter(
    data.longitude,
    data.latitude,
    c=data.bathymetry_m,
    s=0.1,
    transform=ccrs.PlateCarree(),
)
plt.colorbar().set_label("meters")
vd.datasets.setup_baja_bathymetry_map(ax)
plt.tight_layout()
plt.show()

########################################################################################
# We'll create a chain that applies a blocked median to the data, fits a polynomial
# trend, and then fits a standard gridder to the trend residuals.

chain = vd.Chain(
    [
        ("reduce", vd.BlockReduce(np.median, spacing * 111e3)),
        ("trend", vd.Trend(degree=1)),
        ("spline", vd.Spline()),
    ]

def main():
    rp = ccrs.RotatedPole(pole_latitude=45, pole_longitude=180)
    pc = ccrs.PlateCarree()
    
#    ax = plt.subplot(211, projection=rp)
#    ax.bluemarble()
#    ax.coastlines()
    # XXX The change in projection is a result of plt.plot making a poly or a line depending on some things...
    x, y = [-44, -44, 45, 45], [-45, 45, 45, -45]
#    x, y = [-44, -44, 45, 45, -44], [-45, 45, 45, -45, -45]
#    # XXX why is the horizontal line on the native projection not straight? BUG!
#    ax.plot(x, y, marker='o', transform=rp)
#    ax.gridlines()   
    
    
#    ax = plt.subplot(212, projection=pc)
    ax = plt.axes(projection=pc)
    ax = plt.axes(projection=ccrs.InterruptedGoodeHomolosine())
#    ax = plt.axes(projection=rp)

vmin = np.min(plot_vals)
            vmax = np.max(plot_vals)
        norm = normalize_colors(
            vmin, vmax, is_difference=use_cmap_RdBu, log_color_scale=log_color_scale
        )
    if xticklabels == []:
        xticklabels = ["{}$\degree$".format(x) for x in xtick_positions]

    if unit == "" and data_is_xr:
        unit = plot_vals.units.strip()

    if ax == None:
        if plot_type == "zonal_mean":
            ax = plt.axes()
        if plot_type == "single_level":
            ax = plt.axes(projection = ccrs.PlateCarree())

    if title == "fill" and data_is_xr:
        title = plot_vals.name

    # Create plot
    ax.set_title(title)
    if plot_type == "zonal_mean":
        if pedge.all() == -1:
            pedge = GEOS_72L_grid.p_edge()
        if pedge_ind.all() == -1:
            pedge_ind = np.where((pedge <= np.max(pres_range)) & (pedge >= np.min(pres_range)))
            pedge_ind = pedge_ind[0]
            # Pad edges if subset does not include surface or TOA so data spans entire subrange
            if min(pedge_ind) != 0:
                pedge_ind = np.append(min(pedge_ind) - 1, pedge_ind)
            if max(pedge_ind) != 72:

if shift is None:
            # flake B006
            shift = dict(units="dots", x=15)

        if text_kw is None:
            text_kw = {}
        else:
            # since we may modify it, let's make a copy
            text_kw = {**text_kw}

        if "projection" in ax.__dict__ and "transform" not in kwargs:
            from cartopy.crs import PlateCarree
            from matplotlib.transforms import offset_copy

            kwargs["transform"] = PlateCarree()
            geodetic_transform = PlateCarree()._as_mpl_transform(ax)
            text_kw["transform"] = offset_copy(geodetic_transform, **shift)

        if "color" not in kwargs:
            kwargs["color"] = "black"

        if "s" not in text_kw:
            if hasattr(self, "callsign"):
                text_kw["s"] = getattr(self, "callsign")  # noqa: B009
            if hasattr(self, "name"):
                text_kw["s"] = getattr(self, "name")  # noqa: B009

        cumul: List[Artist] = []
        cumul.append(ax.scatter(self.longitude, self.latitude, **kwargs))

        west, east, south, north = ax.get_extent(PlateCarree())
        if west <= self.longitude <= east and south <= self.latitude <= north:

#######################################
# **Plotting the Isentropic Analysis**

# Set up our projection
crs = ccrs.LambertConformal(central_longitude=-100.0, central_latitude=45.0)

# Coordinates to limit map area
bounds = [(-122., -75., 25., 50.)]
# Choose a level to plot, in this case 296 K
level = 0

fig = plt.figure(figsize=(17., 12.))
add_metpy_logo(fig, 120, 245, size='large')
ax = fig.add_subplot(1, 1, 1, projection=crs)
ax.set_extent(*bounds, crs=ccrs.PlateCarree())
ax.add_feature(cfeature.COASTLINE.with_scale('50m'), linewidth=0.75)
ax.add_feature(cfeature.STATES, linewidth=0.5)

# Plot the surface
clevisent = np.arange(0, 1000, 25)
cs = ax.contour(lon, lat, isentprs[level, :, :], clevisent,
                colors='k', linewidths=1.0, linestyles='solid', transform=ccrs.PlateCarree())
ax.clabel(cs, fontsize=10, inline=1, inline_spacing=7,
          fmt='%i', rightside_up=True, use_clabeltext=True)

# Plot RH
cf = ax.contourf(lon, lat, isentrh[level, :, :], range(10, 106, 5),
                 cmap=plt.cm.gist_earth_r, transform=ccrs.PlateCarree())
cb = fig.colorbar(cf, orientation='horizontal', extend='max', aspect=65, shrink=0.5, pad=0.05,
                  extendrect='True')
cb.set_label('Relative Humidity', size='x-large')

# markers = ['o', 's', 'v', '^', '>', '<', '*', 'p', 'D', 'h']

    feature_num = 0

    bounds = None

    for feature in featuredat['features']:
        geomtype = feature['geometry']['type']
        shape = shapely.geometry.shape(feature['geometry'])
        if(max_length > 0.0):
            shape = subdivide_geom(shape, geomtype, max_length)

        featurename = feature['properties']['name']
        print '  feature: {}'.format(featurename)

        refProjection = cartopy.crs.PlateCarree()

        color = colors[feature_num % len(colors)]
        # marker = markers[feature_num % len(markers)]

        if geomtype in ['Polygon', 'MultiPolygon']:
            props = {'linewidth': 2.0, 'edgecolor': color, 'alpha': 0.4,
                     'facecolor': color}
        elif geomtype in ['LineString', 'MultiLineString']:
            props = {'linewidth': 4.0, 'edgecolor': color, 'alpha': 1.,
                     'facecolor': 'none'}

        if bounds is None:
            bounds = list(shape.bounds)
        else:
            # expand the bounding box
            bounds[:2] = np.minimum(bounds[:2], shape.bounds[:2])

time = num2date(time_var[:].squeeze(), time_var.units)

# Combine 1D latitude and longitudes into a 2D grid of locations
lon_2d, lat_2d = np.meshgrid(lon, lat)

# Smooth mslp data
mslp = ndimage.gaussian_filter(mslp, sigma=3, order=0)

#####################################
# Begin making figure

# Create new figure
fig = plt.figure(figsize=(15, 12), facecolor='black')

# Add the map and set the extent
ax = plt.axes(projection=ccrs.PlateCarree())
ax.set_extent([-108., -91., 33., 45.])
ax.background_patch.set_fill(False)

# Add state boundaries to plot
ax.add_feature(cfeature.STATES, edgecolor='white', linewidth=2)

# Contour the MSLP
c = ax.contour(lon_2d, lat_2d, mslp, colors='lime', linewidths=6)
ax.clabel(c, fontsize=12, inline=1, inline_spacing=4, fmt='%i')

wslice = slice(1, None, 4)
# Plot 850-hPa wind vectors
vectors850 = ax.quiver(lon_2d[wslice, wslice], lat_2d[wslice, wslice],
                       u_wind850[wslice, wslice], v_wind850[wslice, wslice],
                       headlength=4, headwidth=3, angles='xy', scale_units='xy',
                       scale=12, color='gold', label='850mb wind')

# Plot #3: Contourf overlayed by coloured point data
    plt.figure()
    qplt.contourf(air_pressure)
    iplt.points(air_pressure, c=air_pressure.data)
    plt.gca().coastlines()
    iplt.show()

    # For the purposes of this example, add some bounds to the latitude
    # and longitude
    air_pressure.coord("grid_latitude").guess_bounds()
    air_pressure.coord("grid_longitude").guess_bounds()

    # Plot #4: Block plot
    plt.figure()
    plt.axes(projection=ccrs.PlateCarree())
    iplt.pcolormesh(air_pressure)
    plt.gca().stock_img()
    plt.gca().coastlines()
    iplt.show()