pandas 如何将等高线图叠加在底图上

声明:本页面是StackOverFlow热门问题的中英对照翻译,遵循CC BY-SA 4.0协议,如果您需要使用它,必须同样遵循CC BY-SA许可,注明原文地址和作者信息,同时你必须将它归于原作者(不是我):StackOverFlow 原文地址: http://stackoverflow.com/questions/26872337/
Warning: these are provided under cc-by-sa 4.0 license. You are free to use/share it, But you must attribute it to the original authors (not me): StackOverFlow

提示:将鼠标放在中文语句上可以显示对应的英文。显示中英文
时间:2020-09-13 22:39:21  来源:igfitidea点击:

How can I get my contour plot superimposed on a basemap

pythonmatplotlibpandasgismatplotlib-basemap

提问by Zilore Mumba

This is a question I asked several months ago and am still struggling to come to a solution. My code gives me a basemap and a contour plot side by side (but printing to file only gives the contour plot), but I want them superimposed. The best solution would the one here https://gist.github.com/oblakeobjet/7546272but this does not show how to introduce the data, and it is difficult when you are learning from scratch online. Without tiring very kind people, I hope the solution is easy as changing a line of code and that someone can help. My code

这是我几个月前问过的一个问题,目前仍在努力寻找解决方案。我的代码并排提供了底图和等高线图(但打印到文件只提供等高线图),但我希望它们叠加。最好的解决方案是这里的https://gist.github.com/oblakeobjet/7546272但这并没有展示如何引入数据,并且当您在线从头开始学习时很难。在不让非常善良的人感到疲倦的情况下,我希望解决方案就像更改一行代码一样简单,并且有人可以提供帮助。我的代码

#!/usr/bin/python
# vim: set fileencoding=UTF8

import numpy as np
import pandas as pd
from matplotlib.mlab import griddata
from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt

#fig = plt.figure(figsize=(10,8))  #when uncommented draws map with colorbar but no contours

#prepare a basemap

m = Basemap(projection = 'merc',llcrnrlon = 21, llcrnrlat = -18, urcrnrlon = 34, urcrnrlat = -8, resolution='h')

# draw country outlines.

m.drawcountries(linewidth=0.5, linestyle='solid', color='k', antialiased=1, ax=None, zorder=None)
m.drawmapboundary(fill_color = 'white')
m.fillcontinents(color='coral',lake_color='blue')
parallels = np.arange(-18, -8, 2.)
m.drawparallels(np.arange(-18, -8, 2.), color = 'black', linewidth = 0.5)
m.drawparallels(parallels,labels=[True,False,False,False])
meridians = np.arange(22,34, 2.)
m.drawmeridians(np.arange(21,36, 2.), color = '0.25', linewidth = 0.5)
m.drawmeridians(meridians,labels=[False,False,False,True])

fig = plt.figure(figsize=(10,8))       # At this position or commented draws teo figures side by side

#-- Read the data.

data = pd.read_csv('../../data/meansr.txt', delim_whitespace=True)

#-- Now gridding data.  First making a regular grid to interpolate onto

numcols, numrows = 300, 300
xi = np.linspace(data.Lon.min(), data.Lon.max(), numcols)
yi = np.linspace(data.Lat.min(), data.Lat.max(), numrows)
xi, yi = np.meshgrid(xi, yi)

#-- Interpolating at the points in xi, yi

x, y, z = data.Lon.values, data.Lat.values, data.Z.values
zi = griddata(x, y, z, xi, yi)

#-- Display and write the results

m = plt.contourf(xi, yi, zi)
plt.scatter(data.Lon, data.Lat, c=data.Z, s=100,
       vmin=zi.min(), vmax=zi.max())
fig.colorbar(m)
plt.savefig("rainfall.jpg", format="jpg")

The plots I get look like this contour plotand the basemap

我得到的情节是这样的等高线图,并底图

and my data

和我的数据

32.6  -13.6   41
27.1  -16.9   43
32.7  -10.2   46
24.2  -13.6   33
28.5  -14.4   43
28.1  -12.6   33
27.9  -15.8   46
24.8  -14.8   44
31.1  -10.2   35
25.9  -13.5   24
29.1   -9.8   10
25.8  -17.8   39
33.2  -12.3   44
28.3  -15.4   46
27.6  -16.1   47
28.9  -11.1   31
31.3   -8.9   39
31.9  -13.3   45
23.1  -15.3   31
31.4  -11.9   39
27.1  -15.0   42
24.4  -11.8   15
28.6  -13.0   39
31.3  -14.3   44
23.3  -16.1   39
30.2  -13.2   38
24.3  -17.5   32
26.4  -12.2   23
23.1  -13.5   27

回答by urschrei

You're almost there, but Basemap can be temperamental, and you have to manage the z-order of plots / map details. Also, you have to transform your lon / lat coordinates to map projection coordinatesbefore you plot them using basemap.

您快到了,但 Basemap 可能会喜怒无常,您必须管理绘图/地图详细信息的 z 顺序。此外,在使用底图绘制它们之前,您必须将 lon / lat 坐标转换为映射投影坐标

Here's a complete solution, which gives the following output. I've changed some colours and linewidths in order to make the whole thing more legible, YMMV. I've also scaled the size of the scatter points by the normalised 'mean' value (data['Z']) – you can simply remove it and substitute e.g. 50if you'd prefer a constant size (it'll look like the largest marker).

这是一个完整的解决方案,它提供了以下输出。我已经更改了一些颜色和线宽,以使整个内容更清晰,YMMV。我还通过归一化的“平均值”值 ( data['Z'])缩放了散点的大小——您可以简单地将其删除并替换,例如,50如果您更喜欢恒定大小(它看起来像最大的标记)。

Please also specify the units of rainfall, and the duration of the measurement which generated the means, if possible:

如果可能,还请指定降雨量单位以及产生平均值的测量持续时间:

Interpolated rainfall data, scatter points scaled by value

插值降雨数据,按值缩放的散点

import numpy as np
import pandas as pd
from matplotlib.mlab import griddata
from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
from matplotlib.colors import Normalize
%matplotlib inline

# set up plot
plt.clf()
fig = plt.figure(figsize=(10, 8))
ax = fig.add_subplot(111, axisbg='w', frame_on=False)

# grab data
data = pd.read_csv('../../data/meansr.txt', delim_whitespace=True)
norm = Normalize()

# define map extent
lllon = 21
lllat = -18
urlon = 34
urlat = -8

# Set up Basemap instance
m = Basemap(
    projection = 'merc',
    llcrnrlon = lllon, llcrnrlat = lllat, urcrnrlon = urlon, urcrnrlat = urlat,
    resolution='h')

# transform lon / lat coordinates to map projection
data['projected_lon'], data['projected_lat'] = m(*(data.Lon.values, data.Lat.values))

# grid data
numcols, numrows = 1000, 1000
xi = np.linspace(data['projected_lon'].min(), data['projected_lon'].max(), numcols)
yi = np.linspace(data['projected_lat'].min(), data['projected_lat'].max(), numrows)
xi, yi = np.meshgrid(xi, yi)

# interpolate
x, y, z = data['projected_lon'].values, data['projected_lat'].values, data.Z.values
zi = griddata(x, y, z, xi, yi)

# draw map details
m.drawmapboundary(fill_color = 'white')
m.fillcontinents(color='#C0C0C0', lake_color='#7093DB')
m.drawcountries(
    linewidth=.75, linestyle='solid', color='#000073',
    antialiased=True,
    ax=ax, zorder=3)

m.drawparallels(
    np.arange(lllat, urlat, 2.),
    color = 'black', linewidth = 0.5,
    labels=[True, False, False, False])
m.drawmeridians(
    np.arange(lllon, urlon, 2.),
    color = '0.25', linewidth = 0.5,
    labels=[False, False, False, True])

# contour plot
con = m.contourf(xi, yi, zi, zorder=4, alpha=0.6, cmap='RdPu')
# scatter plot
m.scatter(
    data['projected_lon'],
    data['projected_lat'],
    color='#545454',
    edgecolor='#ffffff',
    alpha=.75,
    s=50 * norm(data['Z']),
    cmap='RdPu',
    ax=ax,
    vmin=zi.min(), vmax=zi.max(), zorder=4)

# add colour bar and title
# add colour bar, title, and scale
cbar = plt.colorbar(conf, orientation='horizontal', fraction=.057, pad=0.05)
cbar.set_label("Mean Rainfall - mm")

m.drawmapscale(
    24., -9., 28., -13,
    100,
    units='km', fontsize=10,
    yoffset=None,
    barstyle='fancy', labelstyle='simple',
    fillcolor1='w', fillcolor2='#000000',
    fontcolor='#000000',
    zorder=5)

plt.title("Mean Rainfall")
plt.savefig("rainfall.png", format="png", dpi=300, transparent=True)
plt.show()

Using matplotlib's griddatamethod is convenient, but it can also be slow. As an alternative, you can use scipy's griddatamethods, which are both faster, and more flexible:

使用 matplotlib 的griddata方法很方便,但也可能很慢。作为替代方案,您可以使用 scipy 的griddata方法,它们更快、更灵活:

from scipy.interpolate import griddata as gd

zi = gd(
    (data[['projected_lon', 'projected_lat']]),
    data.Z.values,
    (xi, yi),
    method='linear')

If you use scipy's griddatamethod, you'll also have to determine which of the methods (nearest, linear, cubic) give the best resulting plot.

如果您使用 scipy 的griddata方法,您还必须确定哪种方法 ( nearest, linear, cubic) 给出了最好的结果图。

I should add that the interpolation methods demonstrated and discussed above are the simplest possible, and aren't necessarily valid for the interpolation of rainfall data. This articlegives a good overview of valid approaches and considerations for use in hydrology and hydrological modelling. The implementation of these (probably using Scipy) is left as an exercise &c.

我应该补充一点,上面演示和讨论的插值方法是最简单的方法,不一定对降雨数据的插值有效。本文很好地概述了用于水文和水文建模的有效方法和注意事项。这些(可能使用 Scipy)的实现留作练习 &c。

回答by heltonbiker

I have not everything installed here to run your code, but you should try plotting to the basemap myou created, like this:

我没有在这里安装所有东西来运行您的代码,但您应该尝试绘制到m您创建的底图,如下所示:

# fig = plt.figure(figsize=(10,8)) # omit this at line 28

(...)

m.contourf(xi, yi, zi)
m.scatter(data.Lon, data.Lat, c=data.Z, s=100,
   vmin=zi.min(), vmax=zi.max())

(please tell if this doesn't work)

(请告诉这是否不起作用)