爱气象,爱气象家园! 

气象家园

 找回密码
 立即注册

QQ登录

只需一步,快速开始

新浪微博登陆

只需一步, 快速开始

搜索
查看: 73357|回复: 56

[经验总结] 介绍一种绘图白化的方法

[复制链接]

新浪微博达人勋

发表于 2014-6-5 11:14:19 | 显示全部楼层 |阅读模式

登录后查看更多精彩内容~

您需要 登录 才可以下载或查看,没有帐号?立即注册 新浪微博登陆

x
本帖最后由 po_po1 于 2014-6-10 08:42 编辑

首先说明一下,这里的“白化”是指绘图区域限制在一个多边形区域,多边形外围没有图像。
一般来说,白化的程序是这样的,比如说绘制等值线填充图:
#####################################  
cs = plt.contourf(x,y,pr)#绘制等值线填充
for collection in cs.collections:
        collection.set_clip_path(my_path)#设置显示区域
#####################################

通常而言,只要是与matplotlib耦合的绘图程序都可以使用这个程序,包括Basemap、Cartopy等。
上面的程序有个问题,就是my_path是什么?其实就是Path元组,楼主通过两种方式获得:
方法一:(这个方法仅限于Cartopy,但其它方法可以借鉴这种思想)
import cartopy.io.shapereader as shpreader
from cartopy.mpl.patch import geos_to_path
file=r'D:/baiduyundownload/map/China boundary/province boundary/xibeidiqu/xibei.shp'#指定一条你自己的shapefile路径,我这里指定了一个世界地图
reader = shpreader.Reader(file)
countries = reader.records()
us_multipoly, = [country.geometry for country in countries
                          if country.attributes['NAME'] == 'China']#选择地图属性下'NAME'属性里名字是'China'的一条多边形
main_us_geom = sorted(us_multipoly.geoms, key=lambda geom: geom.area)[-1]
us_path, = geos_to_path(main_us_geom)

##########################
上面这些步骤是固定的,所以复制来直接使用就行
方法二
#读取一个描述多边形区域经纬度的坐标文件,比如楼主读取了一个txt文件
file=r'D:/Desktop/map/my location1.txt'

#然后
a=np.loadtxt(file)#具体参数自己填,获得了经纬度坐标
x=a[:,0]
y=a[:,1]

#比如x是经度坐标点集,y是纬度坐标点集,然后from matplotlib.path import Path
us_path=Path(zip(x,y))#这里的地图数据文件大家可以从软件MeteoInfo上弄,很方便的这样就生成了前面提到的us_path(名字不是重点,只是一个代称)





密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2019-1-9 11:20:48 | 显示全部楼层
这一步的时候显示错误,这是为啥呢
main_us_geom = sorted(us_multipoly.geoms, key=lambda geom: geom.area)[-1]

AttributeError: 'Polygon' object has no attribute 'geoms'
密码修改失败请联系微信:mofangbao
回复 支持 1 反对 0

使用道具 举报

新浪微博达人勋

 楼主| 发表于 2014-6-5 11:18:23 | 显示全部楼层
附带楼主上次提到的那个图和程序: figure_1.png

下面是程序:
import cartopy.io.shapereader as shpreader
from cartopy.mpl.patch import geos_to_path
file=r'D:/baiduyundownload/map/map/world map china line.shp'
reader = shpreader.Reader(file)
countries = reader.records()
us_multipoly, = [country.geometry for country in countries
                 if country.attributes['NAME'] == 'China']
main_us_geom = sorted(us_multipoly.geoms, key=lambda geom: geom.area)[-1]
us_path, = geos_to_path(main_us_geom)
############################
#以上地图设置阶段
############################
import matplotlib.pyplot as plt
import numpy as np
from matplotlib.patches import PathPatch
import cartopy.crs as ccrs
import matplotlib.ticker as mticker
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER , LATITUDE_FORMATTER
plt.figure(figsize=(10, 10))
ax = plt.axes(projection=ccrs.PlateCarree())
ax.set_extent([70,140,10,60])
ax.coastlines()
plate_carre_data_transform = ccrs.PlateCarree()._as_mpl_transform(ax)
ax.set_title('China Plot')
line=ax.gridlines(draw_labels=True,alpha=0.5,linestyle='--')
line.xlabels_top=False
line.ylabels_left=False
#line.xlines=False#只剩下x轴平行线
line.xlocator=mticker.FixedLocator(np.arange(70,140,20))#手动设置x轴刻度
line.xformatter= LONGITUDE_FORMATTER
line.yformatter= LATITUDE_FORMATTER
line.xlabel_style={'size':10,'color':'black'}
line.yalbel_style={'color':'red','weight':'bold'}
############################
#以上地图投影设置
############################
from netCDF4 import Dataset
q=Dataset(r'D:\Desktop\model\CUR\2006pr.nc')
lon=q.variables['lon'][:]
lat=q.variables['lat'][:]
pre=q.variables['pre'][4][:]
x,y=np.meshgrid(lon,lat)
############################
#以上读取nc数据
############################
#绘制等值线
quad_set = plt.contourf(x,y,pre,cmap='spectral',transform=ccrs.PlateCarree())#绘制等值线填充
for collection in quad_set.collections:
    collection.set_clip_path(us_path, plate_carre_data_transform)#设置显示区域
plt.colorbar(quad_set,orientation='horizontal',shrink=0.8,extend='both',pad=0.05)
plt.show()


评分

参与人数 1金钱 +20 收起 理由
何永利 + 20 很给力!

查看全部评分

密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2014-6-5 11:47:16 | 显示全部楼层
必须给32个赞,希望你能坚持下去,加油
密码修改失败请联系微信:mofangbao

新浪微博达人勋

 楼主| 发表于 2014-6-5 15:24:15 | 显示全部楼层
何永利 发表于 2014-6-5 11:47
必须给32个赞,希望你能坚持下去,加油

多谢何师兄夸奖
密码修改失败请联系微信:mofangbao

新浪微博达人勋

 楼主| 发表于 2014-6-5 15:31:14 | 显示全部楼层
本帖最后由 po_po1 于 2014-6-5 15:35 编辑

这是我用上述方法弄的地区图像,当然地图文件是我自己做的: 2000-08preciptation.png


密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2014-6-6 11:56:25 | 显示全部楼层
谢谢分享           
密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2014-6-20 19:18:55 | 显示全部楼层
好东西 牛逼
密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2014-7-14 17:08:56 | 显示全部楼层
楼主牛逼 32个赞
密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2015-3-10 12:25:35 | 显示全部楼层
po_po1 发表于 2014-6-5 15:31
这是我用上述方法弄的地区图像,当然地图文件是我自己做的:

楼主 咨询一下cartopy怎么引入的  你是在Windows平台下吗?
密码修改失败请联系微信:mofangbao

新浪微博达人勋

 楼主| 发表于 2015-3-12 14:52:55 | 显示全部楼层
nuist2006 发表于 2015-3-10 12:25
楼主 咨询一下cartopy怎么引入的  你是在Windows平台下吗?

是在windows平台
密码修改失败请联系微信:mofangbao
您需要登录后才可以回帖 登录 | 立即注册 新浪微博登陆

本版积分规则

Copyright ©2011-2014 bbs.06climate.com All Rights Reserved.  Powered by Discuz! (京ICP-10201084)

本站信息均由会员发表,不代表气象家园立场,禁止在本站发表与国家法律相抵触言论

快速回复 返回顶部 返回列表