如何在 Cartopy 投影中绘制圆的半径

我将使用以下函数来完成，该函数实现了一些基本的三角函数来获取圆：

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

def draw_circle(ax, lon, lat, radius, resolution=100):
    """
    Draws a circle on a Cartopy map

    Parameters:
    - ax: Cartopy axes
    - lon, lat: Center coordinates of the circle
    - radius: Radius of the circle in degrees
    - resolution: Number of points to use for drawing the circle (default is 100)
    """
    theta = np.linspace(0, 2*np.pi, resolution)
    circle_lon = lon + np.cos(theta) * radius
    circle_lat = lat + np.sin(theta) * radius

    ax.plot(circle_lon, circle_lat, transform=ccrs.PlateCarree(), color='red', label='Circle')

下面是一个用法示例：

center_lon, center_lat = 1, 0  # Center of the circle
radius_degrees = 10  # Radius of the circle in degrees

# Create a Cartopy map
fig, ax = plt.subplots(
    subplot_kw={'projection': ccrs.PlateCarree()}
) # See https://scitools.org.uk/cartopy/docs/v0.15/crs/projections.html for projections
ax.stock_img()

# Draw the circle
draw_circle(ax, center_lon, center_lat, radius_degrees)

# Some additional matplotlib commands
ax.set_title('Circle on Cartopy Projection')
ax.legend()

plt.show()

您可能需要访问文档以获取所需的可视化格式： https://scitools.org.uk/cartopy/docs/v0.15/matplotlib/intro.html

此代码应提供通过 X 标记的所需半径线。

代码：

from shapely.geometry import Polygon
from cartopy.geodesic import Geodesic  # from geographiclib
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
import cartopy
import cartopy.feature as cfeature

fig = plt.figure(figsize=(10, 6))
ax = fig.add_subplot(1, 1, 1, projection=ccrs.Mercator())

# Avoid error AttributeError: 'GeoAxesSubplot' object has no attribute '_autoscaleXon'
ax._autoscaleXon = False
ax._autoscaleYon = False

ax.set_extent([18, 28, 59.5, 64.1], crs=ccrs.PlateCarree())
ax.coastlines(linewidth=.5)

# Add the radar distance circle
lon_ika = 23.076
lat_ika = 61.76
ika = [lon_ika, lat_ika]
radius = 250     # km
n_samples = 80

circles = Polygon(Geodesic().circle(lon_ika, lat_ika, radius*1000., n_samples=n_samples))
feature = cfeature.ShapelyFeature([circles], ccrs.PlateCarree(), fc='None', ec="black", lw=1, linestyle="-")

circle = ax.add_feature(feature)

# Adding red dot and name of the radar station to the plot
ax.plot(lon_ika, lat_ika, "o", c='r', transform=ccrs.PlateCarree(), markersize=6, label="Ikaalinen", zorder=30)

# Adding red cross and name of IOP location to the plot
lon_hyy = 24.3
lat_hyy = 61.83
hyy = [lon_hyy, lat_hyy]

# Get (geod_distance, forward and backward azimuth) between 2 points
dist_m, fw_azim, bw_azim  = Geodesic().inverse(ika, hyy).T
# Get (long, lat, forward_azimuth) of target point using direct problem solver
px_lon, px_lat, fwx_azim = Geodesic().direct(ika, fw_azim, radius*1000).T

ax.plot(lon_hyy, lat_hyy, "x", c='r', transform=ccrs.PlateCarree(), markersize=6, label="Hyytiälä", zorder=10)

# Plot the target point on the circle's perimeter
ax.plot(px_lon, px_lat, "x", c='g', transform=ccrs.PlateCarree(), markersize=12, label="Target")

# Plot great-circle arc from circle center to the target point
ax.plot([ika[0], px_lon], [ika[1], px_lat], '.-', color='blue',  transform=ccrs.Geodetic(), zorder=5 )

gl = ax.gridlines(draw_labels=True)

#ax.set_aspect(1)

# Add labels
plt.legend(loc='upper left', fontsize=12, framealpha=1, edgecolor='black')

plt.show()

输出图：