I am drawing a map using basemap from matplotlib. The data are spreaded all over the world, but I just want to retain all the data on the continent and drop those on the ocean. Is there a way that I can filter the data, or is there a way to draw the ocean again to cover the data?
The simplest way is to use basemap's maskoceans.
If for each lat, lon you have a data and you want to use contours: After meshgrid and interpolation:
from scipy.interpolate import griddata as gd
from mpl_toolkits.basemap import Basemap, cm, maskoceans
xi, yi = np.meshgrid(xi, yi)
zi = gd((mlon, mlat),
scores,
(xi, yi),
method=grid_interpolation_method)
#mask points on ocean
data = maskoceans(xi, yi, zi)
con = m.contourf(xi, yi, data, cmap=cm.GMT_red2green)
#note instead of zi we have data now.
Update (much faster than in_land or in_polygon solutions):
If for each lat, lon you don't have any data, and you just want to scatter the points only over land:
x, y = m(lons, lats)
samples = len(lons)
ocean = maskoceans(lons, lats, datain=np.arange(samples),
resolution='i')
ocean_samples = np.ma.count_masked(ocean)
print('{0} of {1} points in ocean'.format(ocean_samples, samples))
m.scatter(x[~ocean.mask], y[~ocean.mask], marker='.', color=colors[~ocean.mask], s=1)
m.drawcountries()
m.drawcoastlines(linewidth=0.7)
plt.savefig('a.png')
The HYRY's answer won't work on new versions of matplotlib (nxutils is deprecated). I've made a new version that works:
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
If you love us? You can donate to us via Paypal or buy me a coffee so we can maintain and grow! Thank you!
Donate Us With