python | 降水数据分析（Ⅰ）绘制全国降水四季分布图_python画降水分布图

相关文章推荐

高兴的四季豆 · 【wrf-python】wrf模拟区域放大_ ...· 1 周前 ·

没有腹肌的排球 · 初次调用c++编写的dll文件_c调用c++ ...· 1 年前 ·

直爽的毛衣 · sql - ...· 1 年前 ·

直爽的马克杯 · Windows server 2016, ...· 1 年前 ·

温暖的拖把 · JavaScript 简介 - 知乎· 1 年前 ·

有腹肌的夕阳 · Springboot集成Camunda使用M ...· 1 年前 ·

上次交了作业（全国气温分布图）之后，老师开始让自己摸索处理数据的方法。为了这周先交一个成果上去 ~~（已经拖了两周了）~~ ，决定画一个四季分布图和一个时序图，用两篇文章来记录。

下面的0.可以忽略

0.批量修改文件后缀 （如txt文件 —>xlsx）（ 此步骤可忽略 ）

（此方法修改后的文件能用wps打开，但其他地方不一定能用）

用记事本创建一个新的txt文件，写入ren *.txt *.xlsx（注意是英文空格）

保存为.bat文件，保存在与需要修改后缀的文件同一目录下

保存后运行，即可看到后缀已修改

SOS：修改后的excel表格无法用office打开，显示已损坏，目前没有解决方案......

用python读取修改后的表格也会报错 Excel file format cannot be determined, you must specify an engine manually.

下面开始正文 —— 画2002.3-2003.2的全国站点降水四季分布图

1.数据处理

1-1 将三个月每天的站点数据合并到一个文本文件中

用记事本创建一个新的txt文件，写入

for %%i in (*.txt) do type %%i>>1.txt（注意是英文空格）（1.txt为合并后文件名称）

即for %%i in (所有.这种格式的文件) do type %%i>>导出后的文件名称.文件格式

与要合并的txt文件保存在同级目录下，运行就可以了

1-2 用excel读取、计算每个站点三个月降水的平均值（excel处理方法见上一篇文章）

2.绘制春夏秋冬分布图 直接上代码

from pickletools import long4
from pkgutil import extend_path
import matplotlib
import matplotlib.pyplot as plt
import numpy as np
from cartopy.mpl.ticker import LongitudeFormatter, LatitudeFormatter
import matplotlib.ticker as mticker
import cartopy.crs as ccrs
import cartopy.feature as cfeature
import pandas as pd
from cartopy.io import shapereader as shpreader
df1 = pd.read_excel(r"D:\春 终版.xlsx")
df2 = pd.read_excel(r"D:\夏 终.xlsx")
df3 = pd.read_excel(r"D:\秋 终.xlsx")
df4 = pd.read_excel(r"D:\冬 终.xlsx")
lat1 = df1['Lat']
lon1 = df1['Lon']
pre1 = df1['平均值项:PRE_Time_2020']
lat2 = df2['Lat']
lon2 = df2['Lon']
pre2 = df2['平均值项:PRE_Time_2020']
lat3 = df3['Lat']
lon3 = df3['Lon']
pre3 = df3['平均值项:PRE_Time_2020']
lat4 = df4['Lat']
lon4 = df4['Lon']
pre4 = df4['平均值项:PRE_Time_2020']
# 建立画布
fig2 = plt.figure(figsize=(20, 12))
norm = matplotlib.colors.Normalize(vmin=0, vmax=40)
proj = ccrs.PlateCarree(central_longitude=105)
leftlon, rightlon, lowerlat, upperlat = (70, 140, 15, 55)  # 根据上下限确定范围，至少为10°
# 春 ------------------------------------------------------------------------------------
#左 底 宽 高
f2_ax1 = fig2.add_axes([0.2, 0.35, 0.45, 0.225], projection=proj)
# 在画布的绝对坐标建立子图
f2_ax1.set_extent([leftlon, rightlon, lowerlat, upperlat],
                  crs=ccrs.PlateCarree())
# 海岸线，50m精度
f2_ax1.add_feature(cfeature.COASTLINE.with_scale('50m'))
# 湖泊数据(但是这个貌似只画了比较大的湖泊，比如贝湖巴湖)
f2_ax1.add_feature(cfeature.LAKES, alpha=0.5)
# 以下6条语句是定义地理坐标标签格式
f2_ax1.set_xticks(np.arange(leftlon, rightlon+10, 10),
                  crs=ccrs.PlateCarree())
f2_ax1.set_yticks(np.arange(lowerlat, upperlat+10, 10),
                  crs=ccrs.PlateCarree())
lon_formatter = LongitudeFormatter()
lat_formatter = LatitudeFormatter()
f2_ax1.xaxis.set_major_formatter(lon_formatter)
f2_ax1.yaxis.set_major_formatter(lat_formatter)
f2_ax1.set_title('Spr', loc='left', fontsize=10)
# 读取shp文件
china = shpreader.Reader(r"D:\bou2_4l.dbf").geometries()
# 绘制中国国界省界九段线等等
f2_ax1.add_geometries(china, ccrs.PlateCarree(),
                      facecolor='none', edgecolor='black', zorder=1)
# 添加南海，实际上就是新建一个子图覆盖在之前子图的右下角
f2_ax2 = fig2.add_axes([0.505, 0.353, 0.04, 0.07], projection=proj)
f2_ax2.set_extent([105, 125, 0, 25], crs=ccrs.PlateCarree())
f2_ax2.add_feature(cfeature.COASTLINE.with_scale('50m'))
china = shpreader.Reader(r"D:\bou2_4l.dbf").geometries()
f2_ax2.add_geometries(china, ccrs.PlateCarree(),
                      facecolor='none', edgecolor='black', zorder=1)
f1 = f2_ax1.scatter(lon1, lat1, s=2, c=pre1, cmap='seismic',
                    transform=ccrs.PlateCarree(), vmin=0, vmax=16)
f2 = f2_ax2.scatter(lon1, lat1, s=2, c=pre1, cmap='seismic',
                    transform=ccrs.PlateCarree(), vmin=0, vmax=16)
# 夏-------------------------------------------------------------------------------------
f2_ax3 = fig2.add_axes([0.5, 0.35, 0.45, 0.225], projection=proj)
# 在画布的绝对坐标建立子图
f2_ax3.set_extent([leftlon, rightlon, lowerlat, upperlat],
                  crs=ccrs.PlateCarree())
# 海岸线，50m精度
f2_ax3.add_feature(cfeature.COASTLINE.with_scale('50m'))
# 湖泊数据(但是这个貌似只画了比较大的湖泊，比如贝湖巴湖)
f2_ax3.add_feature(cfeature.LAKES, alpha=0.5)
# 以下6条语句是定义地理坐标标签格式
f2_ax3.set_xticks(np.arange(leftlon, rightlon+10, 10),
                  crs=ccrs.PlateCarree())
f2_ax3.set_yticks(np.arange(lowerlat, upperlat+10, 10),
                  crs=ccrs.PlateCarree())
lon_formatter = LongitudeFormatter()
lat_formatter = LatitudeFormatter()
f2_ax3.xaxis.set_major_formatter(lon_formatter)
f2_ax3.yaxis.set_major_formatter(lat_formatter)
f2_ax3.set_title('Sum', loc='left', fontsize=10)
# 读取shp文件
china = shpreader.Reader(r"D:\bou2_4l.dbf").geometries()
# 绘制中国国界省界九段线等等
f2_ax3.add_geometries(china, ccrs.PlateCarree(),
                      facecolor='none', edgecolor='black', zorder=1)
# 添加南海，实际上就是新建一个子图覆盖在之前子图的右下角
f2_ax4 = fig2.add_axes([0.805, 0.353, 0.04, 0.07], projection=proj)
f2_ax4.set_extent([105, 125, 0, 25], crs=ccrs.PlateCarree())
f2_ax4.add_feature(cfeature.COASTLINE.with_scale('50m'))
china = shpreader.Reader(r"D:\bou2_4l.dbf").geometries()
f2_ax4.add_geometries(china, ccrs.PlateCarree(),
                      facecolor='none', edgecolor='black', zorder=1)
f3 = f2_ax3.scatter(lon2, lat2, s=2, c=pre2, cmap='seismic',
                    transform=ccrs.PlateCarree(), vmin=0, vmax=16)
f4 = f2_ax4.scatter(lon2, lat2, s=2, c=pre2, cmap='seismic',
                    transform=ccrs.PlateCarree(), vmin=0, vmax=16)
# 秋-------------------------------------------------------------------------------------
f2_ax5 = fig2.add_axes([0.2, 0.05, 0.45, 0.225], projection=proj)
# 在画布的绝对坐标建立子图
f2_ax5.set_extent([leftlon, rightlon, lowerlat, upperlat],
                  crs=ccrs.PlateCarree())
# 海岸线，50m精度
f2_ax5.add_feature(cfeature.COASTLINE.with_scale('50m'))
# 湖泊数据(但是这个貌似只画了比较大的湖泊，比如贝湖巴湖)
f2_ax5.add_feature(cfeature.LAKES, alpha=0.5)
# 以下6条语句是定义地理坐标标签格式
f2_ax5.set_xticks(np.arange(leftlon, rightlon+10, 10),
                  crs=ccrs.PlateCarree())
f2_ax5.set_yticks(np.arange(lowerlat, upperlat+10, 10),
                  crs=ccrs.PlateCarree())
lon_formatter = LongitudeFormatter()
lat_formatter = LatitudeFormatter()
f2_ax5.xaxis.set_major_formatter(lon_formatter)
f2_ax5.yaxis.set_major_formatter(lat_formatter)
f2_ax5.set_title('Fal', loc='left', fontsize=10)
# 读取shp文件
china = shpreader.Reader(r"D:\bou2_4l.dbf").geometries()
# 绘制中国国界省界九段线等等
f2_ax5.add_geometries(china, ccrs.PlateCarree(),
                      facecolor='none', edgecolor='black', zorder=1)
# 添加南海，实际上就是新建一个子图覆盖在之前子图的右下角
f2_ax6 = fig2.add_axes([0.505, 0.053, 0.04, 0.07], projection=proj)
f2_ax6.set_extent([105, 125, 0, 25], crs=ccrs.PlateCarree())
f2_ax6.add_feature(cfeature.COASTLINE.with_scale('50m'))
china = shpreader.Reader(r"D:\bou2_4l.dbf").geometries()
f2_ax6.add_geometries(china, ccrs.PlateCarree(),
                      facecolor='none', edgecolor='black', zorder=1)
f5 = f2_ax5.scatter(lon3, lat3, s=2, c=pre3, cmap='seismic',
                    transform=ccrs.PlateCarree(), vmin=0, vmax=16)
f6 = f2_ax6.scatter(lon3, lat3, s=2, c=pre3, cmap='seismic',
                    transform=ccrs.PlateCarree(), vmin=0, vmax=16)
# 冬-------------------------------------------------------------------------------------
f2_ax7 = fig2.add_axes([0.5, 0.05, 0.45, 0.225], projection=proj)
# 在画布的绝对坐标建立子图
f2_ax7.set_extent([leftlon, rightlon, lowerlat, upperlat],
                  crs=ccrs.PlateCarree())
# 海岸线，50m精度
f2_ax7.add_feature(cfeature.COASTLINE.with_scale('50m'))
# 湖泊数据(但是这个貌似只画了比较大的湖泊，比如贝湖巴湖)
f2_ax7.add_feature(cfeature.LAKES, alpha=0.5)
# 以下6条语句是定义地理坐标标签格式
f2_ax7.set_xticks(np.arange(leftlon, rightlon+10, 10),
                  crs=ccrs.PlateCarree())
f2_ax7.set_yticks(np.arange(lowerlat, upperlat+10, 10),
                  crs=ccrs.PlateCarree())
lon_formatter = LongitudeFormatter()
lat_formatter = LatitudeFormatter()
f2_ax7.xaxis.set_major_formatter(lon_formatter)
f2_ax7.yaxis.set_major_formatter(lat_formatter)
f2_ax7.set_title('Win', loc='left', fontsize=10)
# 读取shp文件
china = shpreader.Reader(r"D:\bou2_4l.dbf").geometries()
# 绘制中国国界省界九段线等等
f2_ax7.add_geometries(china, ccrs.PlateCarree(),
                      facecolor='none', edgecolor='black', zorder=1)
# 添加南海，实际上就是新建一个子图覆盖在之前子图的右下角
f2_ax8 = fig2.add_axes([0.805, 0.053, 0.04, 0.07], projection=proj)
f2_ax8.set_extent([105, 125, 0, 25], crs=ccrs.PlateCarree())
f2_ax8.add_feature(cfeature.COASTLINE.with_scale('50m'))
china = shpreader.Reader(r"D:\bou2_4l.dbf").geometries()
f2_ax8.add_geometries(china, ccrs.PlateCarree(),
                      facecolor='none', edgecolor='black', zorder=1)
f7 = f2_ax7.scatter(lon4, lat4, s=2, c=pre4, cmap='seismic',
                    transform=ccrs.PlateCarree(), vmin=0, vmax=16)
f8 = f2_ax8.scatter(lon4, lat4, s=2, c=pre4, cmap='seismic',
                    transform=ccrs.PlateCarree(), vmin=0, vmax=16)
# colorbar 左 下 宽 高
l = 0.88
b = 0.06
w = 0.015
h = 0.5
# 对应 l,b,w,h；设置colorbar位置；
rect = [l, b, w, h]
cbar_ax = fig2.add_axes(rect)
cb = plt.colorbar(f4, cax=cbar_ax, extend='both')
plt.savefig(r'D:\桌面\Seasonal Precipitation.png', dpi=300)
plt.show()
3.一些参考和解释 
1.为了方便，直接用了第一次画图的代码，然后因为太懒了，就没有用循环来写代码。（实际上是学艺不精） 
2.与第一次画图的区别 
   为了给四幅图设置相同的colorbar，我在网上找了好多，但是都不是用plt.figure.add_axes方法画的图，就不适用于我这个图，只好在scatter中设置上下限。 
   关于为什么不用subplots等快捷的方法画这四幅图：因为图中还有南海子图，我没想明白怎么画子图套小图。。。 
   新增：对于colorbar添加了上下限（在scatter函数中，vmin和vmax）和尖角（在plt.colorbar中，extend=‘both’） 
3.下面是在画图过程中搜索到觉得值得收藏的几篇文章 
(13条消息) matplotlib contourf colorbar放到下方_Python空间绘图Colorbar详解_weixin_39790102的博客-CSDN博客 
(14条消息) [Python3] Matplotlib —— (九) 多子图_TreasureAI的博客-CSDN博客 
最全Python绘制多子图总结 - 知乎 (zhihu.com) 
(14条消息) Matplotlib学习手册A005_Figure的add_axes()方法_Python草堂的博客-CSDN博客_add_axes 
scatter函数参数(14条消息) Python中scatter函数参数详解_幸运六叶草的博客-CSDN博客_python scatter参数详解 
4.不得不感慨官网yyds matplotlib.axes.Axes.scatter — Matplotlib 3.5.2 documentation 
5.希望时序图顺利（合十）
本来想绘制如下的色斑图，但一开始不知到它叫这个名字，从等值线开始查起，发现等值线绘制是个比较大且难的问题，会出现等值点计算、等值点追踪、等值线裁剪等一些列的部分所组成，预想到最后还可能会出现效率问题，然，现今已有很多成熟和软件已集成（实现）了该功能，故在查找方法（方式）过程中小伙伴发现了它原来的真实名字，好了，废话到此结束。
一、色斑图绘制
1. 加载数据
做加载的数据包括，经度、维度、降水值三类数据，我把它们存在了一个csv文件中了，由于数据提前已进行了插值处理，因此这里不再进行插值计算。
				所用数据：2014年365天SURF_CHN_MUL_DAY_2014****格式数据，每个文件是每天所有站点的气象要素记录。
解决办法：遍历每天的站点，找到他属于哪个方格，先计算每天每方格内所有站点的平均，再把每天的降水二维数组加起来（365天）
难点：依次读取文件夹中数据。