有没有办法使用matplotlib/cartopy绘制长度一致的直线?

mutmk8jj  于 2022-12-04  发布在  其他
关注(0)|答案(2)|浏览(182)

我正在使用matplotlib和cartopy在python中绘制覆盖在Map上的线。目前,我只是确定两个点的纬度,并在它们之间绘制一条线。由于我是在这些线上绘制横截面,我想找到一种方法使线的长度相同(比如说300公里长)无论我把它放在Map上的什么地方。这是不是可能的,而不只是使用反复试验和设置点,直到它们是所需的长度?
第一次,第二次,第一次,第二次= [34.5,36,-100,-97] x,y = [第一次,第二次],[第一次,第二次]
ax 1.绘图(x,y,颜色=“黑色”,标记=“o”,z阶数=3,转换= ccrs.PlateCarree(),线宽= 2.5)
下面是我现在使用的代码的相关部分。这是可行的,但我正在寻找一种方法来保持线长度不变,而不是改变端点“lat 1,lat 2,lon 1,lon 2”的值。我设想设置一个线长度,一个中点(lat/lon),以及一个围绕该点旋转的Angular 。我不知道这是否可能,但这就是我所想象的工作方式!
Example of a line that the cross section would be through

p1iqtdky

p1iqtdky1#

pyproj中的Geod类对于这些类型的操作非常方便。下面的示例基于方位角和距离(以米为单位)从给定纬度/经度沿直线移动。
既然你提到从中点开始,你可以为每个方向做两次,以得到完整的线。
在下面的例子中,我只是检索端点,并使用Cartopy(ccrs.Geodetic())对大圆进行插值。但是,您也可以使用相同的Geod对象(参见npts方法)自己进行插值,并沿直线采样给定数量的点。例如,如果您需要这些坐标来提取数据,后一种方法可能会很方便。

import matplotlib.pyplot as plt
import cartopy.crs as ccrs
from pyproj import Geod
import numpy as np

# start with a random point
np.random.seed(0)
lon = np.random.randint(-180,180)
lat  = np.random.randint(-90,90)

# and a random direction and distance
half_width = (1 + np.random.rand()) * 5000000 # meters
azimuth = np.random.randint(0,360)

# calculate the end points
geod = Geod(ellps="WGS84")
lon_end1, lat_end1, azi_rev1 = geod.fwd(lon, lat, azimuth, half_width)
lon_end2, lat_end2, azi_rev2 = geod.fwd(lon, lat, azimuth-180, half_width)

# visualize the result
fig, ax = plt.subplots(
    figsize=(8,4), dpi=86, layout="compressed", facecolor="w", 
    subplot_kw=dict(projection=ccrs.PlateCarree()),
)

ax.set_title(f"{lon=}, {lat=}, {azimuth=}, {half_width=:1.1f}m")
ax.plot(
    [lon_end1, lon, lon_end2], 
    [lat_end1, lat, lat_end2], 
    "ro-",
    transform=ccrs.Geodetic(),
)
ax.coastlines()
ax.set_global()

pbpqsu0x

pbpqsu0x2#

the easiest way to do this would probably be to re-project your data into an equidistant projection, such as azimuthal equidistant , then buffer the point by 300km. You could do this

import cartopy.crs as ccrs
import geopandas as gpd

point = gpd.GeoDataFrame(
    geometry=gpd.points_from_xy([-100], [34.5], crs="epsg:4326")
)

crs = ccrs.AzimuthalEquidistant(-100, 34.5)

circle = point.to_crs(crs).buffer(300000).boundary.to_crs("epsg:4326")

This creates an ellipsoid of points in lat/lon space (a circle in actual distance):

In [17]: circle.iloc[0]
Out[17]: <shapely.geometry.linestring.LineString at 0x18c3db7c0>

In [18]: circle.iloc[0].xy
Out[18]:
(array('d', [-96.73458302693649, -96.76051210175493, -96.81721890848735, -96.90389145413285, -97.0194601113924, -97.16261553916054, -97.33182721184983, -97.52536229026433, -97.74130463130324, -97.97757379088794, -98.23194392302969, -98.5020625175967, -98.78546895046141, -99.07961284313612, -99.38187224585855, -99.68957166961148, -100.0, -100.31042833038852, -100.61812775414145, -100.92038715686388, -101.21453104953859, -101.4979374824033, -101.76805607697031, -102.02242620911204, -102.25869536869676, -102.47463770973567, -102.66817278815017, -102.83738446083946, -102.98053988860761, -103.09610854586715, -103.18278109151265, -103.23948789824507, -103.26541697306351, -103.26003093177158, -103.22308261845816, -103.15462889147989, -103.05504203609833, -102.92501821713451, -102.7655823598689, -102.57808885099398, -102.3642174900262, -102.12596419991539, -101.86562612586265, -101.58578091250277, -101.28926014666801, -100.97911717692438, -100.65858975917035, -100.33105821410338, -100.0, -99.66894178589662, -99.34141024082963, -99.02088282307564, -98.710739853332, -98.41421908749726, -98.13437387413735, -97.87403580008461, -97.63578250997378, -97.42191114900605, -97.23441764013111, -97.07498178286552, -96.94495796390169, -96.84537110852011, -96.77691738154184, -96.73996906822842, -96.73458302693649]),
 array('d', [34.45635448578617, 34.191923718769814, 33.930839345631036, 33.67556892797355, 33.42850443127362, 33.191942126948454, 32.96806390193457, 32.75892001047773, 32.566413289704464, 32.39228485200862, 32.23810126231687, 32.105243205881486, 31.99489565146612, 31.90803951485474, 31.845444827911535, 31.80766541851577, 31.795035106337213, 31.80766541851577, 31.845444827911535, 31.90803951485474, 31.99489565146612, 32.105243205881486, 32.23810126231687, 32.392284852008615, 32.566413289704464, 32.75892001047773, 32.96806390193457, 33.191942126948454, 33.42850443127362, 33.67556892797355, 33.93083934563104, 34.191923718769814, 34.45635448578617, 34.72160994100843, 34.985136962504626, 35.244374905489536, 35.49678051263328, 35.739853647811024, 35.97116361010222, 36.18837573214065, 36.38927791396905, 36.57180669376817, 36.73407241410168, 36.874383010755274, 36.99126593481258, 37.08348772070849, 37.15007073604399, 37.19030669398129, 37.20376657546522, 37.19030669398129, 37.15007073604399, 37.08348772070849, 36.99126593481258, 36.87438301075528, 36.73407241410169, 36.57180669376818, 36.38927791396906, 36.18837573214065, 35.97116361010222, 35.739853647811024, 35.496780512633286, 35.244374905489536, 34.98513696250464, 34.72160994100843, 34.45635448578617]))

相关问题