无法使用 cartopy 显示经纬度和海岸线

Question

我正在使用 netcdf 文件中的 GOES-16 卫星图像。我想裁剪文件，使其仅显示波多黎各和轮廓。我可以收割它，但不会显示经纬度，波多黎各的海岸线也不会显示。下面是我的示例代码。

import netCDF4 as nc
from netCDF4 import Dataset
import cartopy as ccrs
import cartopy.mpl.ticker as cticker
import matplotlib
import matplotlib.pyplot as plt
import matplotlib.ticker as mticker
#from mpl_toolkits.basemap import Basemap 
import metpy
import numpy as np
import cv2
from cv2 import remap
import xarray

fn = 'OR_ABI-L2-ACMC-M6_G16_s20211730001172_e20211730003545_c20211730004314.nc'

ds = xarray.open_dataset(fn)
nc = nc.Dataset(fn)

variable = 'BCM'
bcm = ds['BCM'][:].data

bcm = np.clip(bcm, 0, 1)

geo_extent = nc.variables['geospatial_lat_lon_extent']
print(geo_extent)
#set latitude and longitude values from metadata

min_lon = float(geo_extent.geospatial_westbound_longitude)
max_lon = float(geo_extent.geospatial_eastbound_longitude)
min_lat = float(geo_extent.geospatial_southbound_latitude)
max_lat = float(geo_extent.geospatial_northbound_latitude)

dat = ds.metpy.parse_cf('BCM')

geos = dat.metpy.cartopy_crs
x = dat.x
y = dat.y
print(x.min())

fig = plt.figure(figsize=(12, 15))


ax = fig.add_subplot(1, 1, 1, projection=geos)
ax.set_extent([-68.07, -64.65, 16.7, 19.79], crs=geos)

ax.imshow(bcm, origin='upper', extent=(min_lon, max_lon, min_lat, max_lat), transform=geos)

ax.coastlines(resolution='10m', color='black', linewidth=1)
ax.add_feature(ccrs.cartopy.feature.COASTLINE)

plt.title('Puerto Rico', loc='left', fontweight='bold', fontsize=15)

plt.show() 

这是地图的图片。 [1]：https://i.stack.imgur.com/Y0qbY.png

Answer 1

我想通了。我使用了不正确的语法。

dat = ds.metpy.parse_cf('BCM')

geos = dat.metpy.cartopy_crs
x = dat.x
y = dat.y


pc = ccrs.crs.PlateCarree()


fig = plt.figure(figsize=(12, 15))


ax = fig.add_subplot(1, 1, 1, projection = pc)

ax.set_extent([-68.07, -64.65, 16.7, 19.79], crs = pc) 


ax.imshow(bcm, origin='upper', extent=(x.min(), x.max(), y.min(), y.max()), transform=geos, interpolation = 'none')

ax.coastlines(resolution='10m', color='black', linewidth=1)
ax.add_feature(ccrs.cartopy.feature.COASTLINE)


ax.set_xticks((-68.07, -64.65), minor = 'True', crs= pc)

ax.set_yticks((16.7, 19.79), minor= 'True', crs= pc)



plt.title('Puerto Rico Clear-Sky Mask', loc='center', fontweight='bold', fontsize=15)


plt.show()