三探ECMWF（+NOAA）——ECMWF3.2

一、ECMWF的数据

from netCDF4 import Dataset
from netCDF4 import num2date

nc_obj = Dataset ( "12.5.nc" )

#取出时间
time=nc_obj['time'][:]

#将时间转化为人类可读时间 datetime.datetime
time1=num2date(time,units=nc_obj['time'].units)

#经度  128.5   12
longitude=nc_obj['longitude'][:]
print(longitude)

[127.    127.125 127.25  127.375 127.5   127.625 127.75  127.875 128.
 128.125 128.25  128.375 128.5   128.625 128.75  128.875 129.   ]

#纬度 45.5    4
latitude=nc_obj['latitude'][:]
print(latitude)

[46.    45.875 45.75  45.625 45.5   45.375 45.25  45.125 45.    44.875
 44.75  44.625 44.5   44.375 44.25  44.125 44.   ]

#取出变量blh  10 metre V wind component(10米V风分量)
v10=nc_obj['v10'][:]

#取出北纬45.5东经128.5的v10数值
a = [0]*124
for i in range(124):
    a[i]=nc_obj.variables['v10'][i][4][12]

#取出变量10 metre U wind component,(10米U风分量)
u10=nc_obj['u10'][:]

#取出北纬45.5东经128.5的u10数值
b = [0]*124
for i in range(124):
    b[i]=nc_obj.variables['u10'][i][4][12]

#对v10和u10平方和开平方
c = [0]*124
for i in range(124):
    c[i]=(a[i]**2+b[i]**2)**0.5

import matplotlib.pyplot as plt
import  pandas as pd

plt.rc('font', family='SimHei', size=15) #绘图中的中文显示问题，图表字体为SimHei，字号为15
#plt.rcParams['font.sans-serif']=['SimHei'] #用来正常显示中文标签
plt.rcParams['axes.unicode_minus']=False #用来正常显示负号

plt.figure(figsize=(30,8))
plt.title('ECMWF风速与时间的函数')  #有标题(风速与风向的函数)
plt.xlabel('时间') #横坐标的标题
plt.ylabel('风速（m/s）') #纵坐标的标题
plt.grid(color='#95a5a6',linestyle='--',linewidth=3,axis='both',alpha=0.4) #设置网格
plt.plot(time1,c, 'ro',time1,c, 'r')
plt.xticks(pd.date_range('2019-01-01','2019-01-31'))
# plt.xticks(rotation='vertical')
# plt.xticks(rotation=30)
plt.gcf().autofmt_xdate()
plt.savefig('ECMWF.png') #保存图片文件命名为3实测与设计功率曲线对比.png
plt.show() 

二、NOAA的数据（梁晨）

from datetime import datetime
import  pandas as pd

noaa=pd.read_csv("C:\\Users\\zzh\\Desktop\\ECMWF\\201901.CSV")

201901.CSV该数据存在问题，在用pandas读进来后，USAF变成行号，其他依次往前挪了一列

#取出时间
time2=noaa['WBAN']

# 201901041500，
# 每一个时间如上所示，为2019-01-04，15：00
# 我把分钟数截断，因为所有数据分钟数都是00
new_time=[]
for i in time2:
    # 201901041500
    new_time.append((str.split(str(i)[:10]))[0])

new_time2 =[]
for i in new_time:
    # 201901041500
    new_time2.append(datetime.strptime(i,'%Y%m%d%H'))

# 隔行选取数据
new_time3 =[]
for i in range(0,247,2):
    new_time3 .append(new_time2[i])

#取出数据，因为数据串列，风速命名为DIR
DIR=noaa['DIR']
# print(DIR)

#将英里每小时转化为米每秒
# 1英里每小时（mile/h）=0.44704米每秒(m/s)
spd=[]
for i in DIR:
    spd.append(i*0.44704)
# print(spd)

# 隔行选取数据
spd2 =[]
for i in range(0,247,2):
    spd2 .append(spd[i])

plt.figure(figsize=(30,8))
plt.title('NOAA风速与时间的函数')  #有标题(风速与风向的函数)
plt.xlabel('时间') #横坐标的标题
plt.ylabel('风速（m/s）') #纵坐标的标题
plt.grid(color='#95a5a6',linestyle='--',linewidth=3,axis='both',alpha=0.4) #设置网格
plt.plot(new_time3,spd2, 'go',new_time3,spd2, 'g')
plt.xticks(pd.date_range('2019-01-01','2019-01-31'))
# plt.xticks(rotation='vertical')
# plt.xticks(rotation=30)
plt.gcf().autofmt_xdate()
plt.savefig('NOAA.png') #保存图片文件命名为3实测与设计功率曲线对比.png
plt.show()

三、ECMWF和NOAA对比图

plt.figure(figsize=(30,8))
plt.title('风速与时间的函数')  #有标题(风速与风向的函数)
plt.xlabel('时间') #横坐标的标题
plt.ylabel('g风速（m/s）') #纵坐标的标题
plt.grid(color='#95a5a6',linestyle='--',linewidth=3,axis='both',alpha=0.4) #设置网格
plt.plot(time1,c, 'ro',time1,c, 'r',label='ECMWF')#ECMWF
plt.plot(new_time3,spd2, 'go',new_time3,spd2, 'g',label='NOAA')#NOAA
plt.xticks(pd.date_range('2019-01-01','2019-01-31'))
# plt.xticks(rotation='vertical')
# plt.xticks(rotation=30)
plt.gcf().autofmt_xdate()
plt.legend()
plt.savefig('ECMWF和NOAA.png') #保存图片文件命名为3实测与设计功率曲线对比.png
plt.show()