12

matplotlib のベースマップを使用してマップを描画しています。データは世界中に散らばっていますが、すべてのデータを大陸に残し、それらを海に落としたいだけです。データをフィルタリングする方法はありますか、またはデータをカバーするために海を再度描画する方法はありますか?

4

5 に答える 5

6

matplotlib.basemapには次のメソッドがあります。 is_land(xpt, ypt)

True指定されたx、yポイント(投影座標内)が陸地上にある場合は戻ります。それ以外の場合は戻りますFalse。土地の定義は、クラスインスタンスに関連付けられたGSHHS海岸線ポリゴンに基づいています。陸域内の湖の上のポイントは、陸地ポイントとしてカウントされません。

詳細については、こちらを参照してください。

于 2012-12-10T07:43:14.447 に答える
6

is_land()すべてのポリゴンをループして、土地かどうかを確認します。データサイズが大きい場合、非常に遅くなります。points_inside_poly()from matplotlib を使用して、ポイントの配列をすばやく確認できます。これがコードです。lakepolygons湖のポイントを削除したい場合は、自分で追加できます。

私のPCで100000ポイントを確認するのに2.7秒かかりました。速度を上げたい場合は、ポリゴンをビットマップに変換できますが、これを行うのは少し困難です。次のコードがデータセットに対して十分に高速でない場合は教えてください。

from mpl_toolkits.basemap import Basemap
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.nxutils as nx

def points_in_polys(points, polys):
    result = []
    for poly in polys:
        mask = nx.points_inside_poly(points, poly)
        result.extend(points[mask])
        points = points[~mask]
    return np.array(result)

points = np.random.randint(0, 90, size=(100000, 2))
m = Basemap(projection='moll',lon_0=0,resolution='c')
m.drawcoastlines()
m.fillcontinents(color='coral',lake_color='aqua')
x, y = m(points[:,0], points[:,1])
loc = np.c_[x, y]
polys = [p.boundary for p in m.landpolygons]
land_loc = points_in_polys(loc, polys)
m.plot(land_loc[:, 0], land_loc[:, 1],'ro')
plt.show()
于 2012-12-11T00:37:41.073 に答える
5

HYRY の回答は、matplotlib の新しいバージョンでは機能しません (nxutils は非推奨です)。動作する新しいバージョンを作成しました:

from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
from matplotlib.path import Path
import numpy as np

map = Basemap(projection='cyl', resolution='c')

lons = [0., 0., 16., 76.]
lats = [0., 41., 19., 51.]

x, y = map(lons, lats)

locations = np.c_[x, y]

polygons = [Path(p.boundary) for p in map.landpolygons]

result = np.zeros(len(locations), dtype=bool) 

for polygon in polygons:

    result += np.array(polygon.contains_points(locations))

print result
于 2014-11-06T12:44:26.620 に答える