python利用IBTrACS的nc版本数据进行区域每年台风路径的绘制

在绘制图像之前我们进行一个对数据的分析，此处我们选取的是第四版的数据，然后我们会发现这一版的数据和之前的第三版有很多的不同，我们先看看数据里面有什么？

1、分析数据

from datetime import time
from os import times
import numpy as np
import pandas as pd
import xarray as xr
import matplotlib.pyplot as plt
import matplotlib as mpl
import matplotlib.patches as patches
import cartopy.crs as ccrs
import cartopy.feature as cfeature
from cartopy.mpl.ticker import LongitudeFormatter,LatitudeFormatter
tcFile =r'C:\\Users\\sha\\Desktop\\IBTrACS.ALL.v04r00.nc'
infile = xr.open_dataset(tcFile)

阅读infile的内部：

infile

我们可以看出这里的参数变量非常多，多达147个，而我们正需要从这里面抽取我们需要的参数，那么我们就打开看看里面的内容是什么：

 根据说明书的内容，我们大概可以知道这些参数对应的意义大概是这些：

SID-由IBTrACS算法分配的唯一风暴标识符（SID）

SEASON-风暴开始的季节

NUMBER-当年风暴的数量（每年重新开始为1次）

BASIN-当前风暴位置的流域

SUBBASIN-当前风暴位置的子流域

NAME-按来源提供的系统名称

 这里我们可以知道数据的来源有usa、jma、cma、hko、new、bom等等，同时还会发现每个地区都有对应的变量：

 我们会发现还是usa的数据会相对比较全，所以我们在这里用的一些参数都是取usa的数据，当然其实每种数据的lat和lon都是大差不差的。

其实还有其他很多的数据，这些数据大家可以选择用jupyter notebook 进行查看，因为尝试用spyder查看的时候就会发现数据太多只能看到部分内容，这样的话对数据没有办法看到全貌不方便进行处理。

2、读取需要的数据

# read variables
lat  = infile['usa_lat']          # storm center latitude
lon  = infile['usa_lon']         # storm center longitude          # 2 = WP - western north Pacific
stormYear = infile['season'].data              # year based on season
sid=infile['sid']
number=infile['number'].data
time2d = infile['iso_time'].data         # time step
name=infile['name'].data
basin=infile['basin'].data
basin1=basin[:,0]
tracktype=infile['track_type'].data

3、选取需要的年份，地区和轨道类型

比如此处我们选取的是2020年的，然后地区是WP，西太平洋地区的台风，然后轨道类型是main，为何需要选择main，因为我们看说明书可以发现，说明书中有说明，有时候会对台风进行过度计算，我们有时候会把spur-split也算在里面，但是这样的话个数就会偏多，所以我们需要去掉spur-split，而是只选取其中的main即可。

tcmask = np.where((stormYear==2019)&(basin1==b'WP')&(tracktype==b'main'))

4、读取相对应满足tcmask的数据

latselect  = lat[tcmask]
lonselect  = lon[tcmask]
nameselect = name[tcmask]
timeselect = time2d[tcmask]
numberselect = number[tcmask]
basinselect= basin[tcmask]
print(len(numberselect))
#此处我们可以输出一个数据，我们就可以准确知道满足年份地区和轨道类型的台风的具体个数，而不需要自己数
latselect1=np.array(latselect)
lonselect1=np.array(lonselect)

5、绘图

fig=plt.figure(figsize=(20,10)) #设置画布大小

ax = plt.axes(projection=ccrs.PlateCarree(central_longitude=180))   #设置投影方式
# Set figure extent & ticks
ax.coastlines()
ax.gridlines(draw_labels=True, dms=True, x_inline=False, y_inline=False)
ax.set_extent([95,180, 0, 60])  #设置纬度范围
ax.add_feature(cfeature.OCEAN)
ax.set_title('TCs in the Northwest Pacific in 2019'+'  ('+str(len(numberselect))+')',fontsize=20, loc='center')
ax.add_feature(cfeature.LAND, edgecolor='b')
for i in range(0,len(numberselect)):
    tclat=latselect1[i,:]
    tclon=lonselect1[i,:]
    ax.plot(tclon,tclat,color='k',linewidth=0.5,transform=ccrs.PlateCarree())
    cb = ax.scatter(tclon,tclat,s=5,transform=ccrs.PlateCarree())

 虽然不太好看，但是画出来啦！

如果要画其他地区的直接在

tcmask = np.where((stormYear==2019)&(basin1==b'WP')&(tracktype==b'main'))对时间和地区进行修改就好啦，如果不限制地区可以把后面的删掉，还有我们要注意的一点是可能我们的nc文件资料里的数据中，对2021年所有的台风的tracktype都是PROVISIONAL，临时性的一个定义，目前也不太清楚原因是为何，但是在绘制2021年的时候将(tracktype==b'main')改为(tracktype==b'PROVISIONAL')就好啦！