【问题标题】:How to plot gridded data in basemap and find grid points along a specific azimuth如何在底图中绘制网格数据并沿特定方位查找网格点
【发布时间】:2013-02-05 11:12:45
【问题描述】:

更新:我正在尝试映射一些数据。我有一组从网格中的参考点测量的后向方位角(baz)。我想在网格上找到沿着 baz 的一个大圆圈会穿过的所有点。为此,我遍历网格中的每个点,计算该点与参考点之间的预期后向方位角,并与每个测量的 baz 进行比较。如果两者之间的差异很小(小于 2 度),我会权衡这一点。然后我把它全部放在地图上。我使用的代码如下,但结果看起来有点奇怪,有谁知道我哪里出错了,或者如果有更好的方法(更快)然后我做了什么??

from matplotlib.colorbar import ColorbarBase
import matplotlib.cm as cm
import matplotlib.pyplot as plt
import numpy as np
from mpl_toolkits.basemap import Basemap
import mpl_toolkits.basemap.pyproj as pyproj


llcrnrlon = -30.0
llcrnrlat = 45.0
urcrnrlon = 0.0
urcrnrlat = 65.0
lon_0 = (urcrnrlon + llcrnrlon) / 2.
lat_0 = (urcrnrlat + llcrnrlat) / 2.

lat = 51.58661577  # reference point
lon = -9.18822525


# Generate random back-azimuths.
baz = zeros((20))
for i in xrange(len(baz)):
  baz[i] = random.randint(200,230)


####################################################################
## Set up the map background.
m = Basemap(llcrnrlon=llcrnrlon,llcrnrlat=llcrnrlat,urcrnrlon=urcrnrlon,urcrnrlat=urcrnrlat,
    resolution='i',projection='lcc',lon_0=lon_0,lat_0=lat_0)
m.drawcoastlines()
m.fillcontinents() 

# draw parallels
m.drawparallels(np.arange(10,70,10),labels=[1,0,0,0])
# draw meridians
m.drawmeridians(np.arange(-80, 25, 10),labels=[0,0,0,1])

# Plot station locations.
x, y = m(lon, lat)            # array ref points
m.plot(x,y,'ro', ms=5)


####################################################################
## Set up the grids etc.
glons = np.linspace(llcrnrlon, urcrnrlon, 100)
glats = np.linspace(llcrnrlat, urcrnrlat, 100)

# Convert to map coords.
xlons, ylats = m(glons, glats)

# create grid for pcolormesh.
grid_lon, grid_lat = np.meshgrid(xlons, ylats)

# create weights for pcolormesh.
weights = np.zeros(np.shape(grid_lon))

# create grid of lat-lon coords for baz calculation.
gln, glt = np.meshgrid(glons, glats)


####################################################################
## calculate baz from grid_lon, grid_lat to lon, lat. If less 
## than error weight grid point.

# method for BAZ calculation via pyproj.
def get_baz(lon1, lat1, lon2, lat2):
  g = pyproj.Geod(ellps='WGS84')
  az, baz, dist = g.inv(lon1, lat1, lon2, lat2)
  return baz

# BAZ calcultion for each point in grid.
ll=0
for mBAZ in baz:
  for i in xrange(len(gln)):
    for k in xrange(len(gln[i])):
      nbaz = get_baz(lon, lat, gln[i][k], glt[i][k])
      nbaz += 180
      if abs(nbaz - mBAZ) < 2:
    weights[i][k] = 1
  ll+=1


# plot grid.
m.pcolormesh(grid_lon, grid_lat, weights, cmap=plt.cm.YlOrBr)
plt.colorbar()
plt.show()

以下原始问题,现已过时。

我正在尝试映射一些数据。我有一个数据集,它为每个方向提供一系列值(频率)。我想将它们绘制在网格上,以便沿特定方位角的每个网格点都由特定频率的功率加权。 我创建了一个带有底图的地图,并在其上绘制了一个网格,如下所示,

from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
import numpy as np
from shoot import *

llcrnrlon = -20.0
llcrnrlat = 45.0
urcrnrlon = 10.0
urcrnrlat = 65.0
lon_0 = (urcrnrlon + llcrnrlon) / 2.
lat_0 = (urcrnrlat + llcrnrlat) / 2.

m = Basemap(llcrnrlon=llcrnrlon,llcrnrlat=llcrnrlat,urcrnrlon=urcrnrlon,urcrnrlat=urcrnrlat,
        resolution='i',projection='lcc',lon_0=lon_0,lat_0=lat_0)

## Set up the grid.
glons = np.linspace(-20,10,50)
glats = np.linspace(45, 65, 50)
xlons, ylats = m(glons, glats)
grid_lon, grid_lat = np.meshgrid(xlons, ylats) 
pwr = np.zeros((50,50))

m.drawcoastlines()
m.fillcontinents() 

# draw parallels
m.drawparallels(np.arange(10,70,10),labels=[1,0,0,0])
# draw meridians
m.drawmeridians(np.arange(-80, 25, 10),labels=[0,0,0,1])

lats = [54.8639587, 51.5641564]
lons = [-8.1778180, -9.2754284]

x, y = m(lons, lats)            # array ref points

# Plot station locations.
m.plot(x,y,'ro', ms=5)
m.pcolormesh(grid_lon, grid_lat, pwr)

然后我使用我在这个不错的site 找到的一些功能射出我想要的大圆圈

glon1 = lons[0]
glat1 = lats[0]
azimuth = 280.
maxdist = 200.
great(m, glon1, glat1, azimuth, color='orange', lw=2.0)
plt.show()

但是,仅绘制线是不够的,我希望能够找到大圆穿过的网格点,以便为它们分配值。有谁知道这个怎么办??

【问题讨论】:

    标签: python map matplotlib matplotlib-basemap


    【解决方案1】:

    你能具体说明你指的是哪个交叉点吗?运行您的代码只返回一行...

    【讨论】:

    • 我需要在网格(grid_lon,grid_lat)中找到线交叉的所有网格点。
    猜你喜欢
    • 1970-01-01
    • 2019-03-25
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 2022-10-04
    • 1970-01-01
    相关资源
    最近更新 更多