python处理nc气象数据_气象数据处理——nc文件
數據說明NetCDF(network Common Data Form)網絡通用數據格式是一種面向數組型并適于網絡共享的數據的描述和編碼標準。目前,NetCDF廣泛應用于大氣科學、水文、海洋學、環境模擬、地球物理等諸多領域。用戶可以借助多種方式方便地管理和操作 NetCDF 數據集。
NetCDF全稱為network Common Data Format,中文譯法為“網絡通用數據格式”;netcdf文件開始的目的是用于存儲氣象科學中的數據,現在已經成為許多數據采集軟件的生成文件的格式。?從數學上來說,netcdf存儲的數據就是一個多自變量的單值函數。用公式來說就是f(x,y,z,…)=value;
?函數的自變量x,y,z等在netcdf中叫做維(dimension) 或坐標軸(axix),
?函數值value在netcdf中叫做變量(Variables).
一個Netcdf文件的結構包括以下對象:
?變量(Variables) :變量對應著真實的物理數據。
?維(dimension):一個維對應著函數中的某個自變量,或者說函數圖象中的一個坐標軸,在線性代數中就是一個N維向量的一個分量。
?屬性(Attribute) :屬性對變量值和維的具體物理含義的注釋或者說解釋。
(原文鏈接:CSDN-專業IT技術社區-登錄)Python讀取:使用netCDF4的Dataset方法讀入文件
# -*- coding: utf-8 -*-
from netCDF4 import Dataset
import numpy as np
import sys
import matplotlib.pyplot as plt
import numpy as np
from mpl_toolkits.basemap import Basemap
from pandas import DataFrame
#數據讀入
nc=Dataset('bj2016_18pm.nc')
print(nc.variables.keys())
#該文件是輻射資料,來自ECMWF網站
odict_keys(['longitude', 'latitude', 'time', 'ssrd', 'ssr', 'fdir', 'strd', 'str'])
#取出各variable的數據看看,數據格式為numpy數組
for var in nc.variables.keys():
data=nc.variables[var][:].data
print(var,data.shape)
# np.save(var+'.npy',data)
longitude (1,)
latitude (2,)
time (4392,)
ssrd (4392, 2, 1)
ssr (4392, 2, 1)
fdir (4392, 2, 1)
strd (4392, 2, 1)
str (4392, 2, 1)
#time variable查看,時間戳變換
#看出是逐時數據
import datetime
time=nc.variables['time'][:].data
print(time[:10])
for i in range(3):
tstamp=(time[i]-613608)*3600 #1900年1月1日零時距離1970年1月1日零時有613608個小時
date= datetime.datetime.utcfromtimestamp(tstamp)
print (date.strftime("%Y-%m-%d%H:%M:%S"))
[1016851 1016852 1016853 1016854 1016855 1016856 1016857 1016858 1016859
1016860]
2016-01-01 19:00:00
2016-01-01 20:00:00
2016-01-01 21:00:00
#查看longitude和 latitude這兩個variable
#這個數據比較小,只包含兩個格點,所以直接輸出了
print(nc.variables['longitude'][:].data,nc.variables['latitude'][:].data)
[116.5] [39.875 39.75 ]
#查看輻射數據&數據切片
ssrd=nc.variables['ssrd'][:].data
#numpy數組切片,取數
ssrd_timeseq=ssrd[:,1,0]#取出某個格點所有時間的ssrd值
ssrd_timei=ssrd[0,:,:]#取出某個時間點所有格點的ssrd值Python作圖:使用matplotlib和Basemap
nc=Dataset('sfc_201803_2mt.nc')#數據來自ECMWF,為2018年3月份的全球溫度的格點數據
print(nc.variables.keys())
data=nc.variables['t2m'][:]
print(data.shape)
#d1=data[0,:,:].data##取出某個時刻的溫度
#輸出結果表明有744個時間序列,經度、緯度的取值有720、361個
odict_keys(['longitude', 'latitude', 'time', 't2m'])
(744, 361, 720)
#查看數據經緯度范圍,經度0-360,其中0~180為西經0~180,180~360為東經0~180;緯度正為北緯,負為南緯
#格點分辨率為0.5度
long= nc.variables['longitude'][:]
lati= nc.variables['latitude'][:]
print(long[0],long[-1],lati[0],lati[-1])
print(long.shape,lati.shape)
0.0 359.5 90.0 -90.0
(720,) (361,)
#plt對某個時刻的全球溫度作圖。左半部分為西半球,右邊是東半球
#選了某個時間點116作圖
plt.contourf(long,lati,data[116,:,:]-273) #轉為攝氏度
plt.colorbar()
用Basemap畫地圖
#用Basemap畫地圖
def graph(lon,lat,target,levelT,colorT,title):
b_map=Basemap(resolution='l', area_thresh=10000, projection='cyl',
llcrnrlon=min(lon), urcrnrlon=max(lon), llcrnrlat=min(lat),urcrnrlat=max(lat))
#llcrnrlon=0, urcrnrlon=360, llcrnrlat=-90,urcrnrlat=90
print(type(target))
fig=plt.figure(figsize=(9, 6)) #plt.figure(figsize=(12, 8))
ax=fig.add_axes([0.1,0.1,0.8,0.8])
lon,lat=np.meshgrid(lon,lat)
x,y=b_map(lon,lat)
print(x.shape,y.shape,target.shape)
cs=b_map.contourf(x,y,target,levels=levelT,colors=colorT) #target[0,:,:]
b_map.colorbar(cs)
b_map.drawcoastlines(linewidth=1)
b_map.drawcountries(linewidth=1.5)
plt.title(title,size=20)
#plt.savefig('Rainf_0.png',dpi=300)
plt.show()
plt.close()截取中國范圍來作圖:
#首先查找中國四至范圍對應的索引
#4,54為北緯4,54;73+180,133+180為東經73,133
print(np.argwhere(lati==4),np.argwhere(lati==54),
np.argwhere(long==73+180),np.argwhere(long==133+180))
#可以看到,中國范圍 對應的緯度索引為72:172,,經度索引為506:626
[[172]] [[72]] [[506]] [[626]]
#對中國范圍的溫度作圖,設定graph的參數
title='2m_temperature'
level_Tair= [-20,-10,0,5,10,15,20,25,30,1000] #[0,2.6,5,8,16,50,100,120,1000]
#[0,2.6,5,8,10,20,25,300,1000] [0,210,225,240,255,260,300,305,310,1000]
colors = ['#FFFFFF', '#AAF0FF', '#C8DC32', '#FFBE14', '#FF780A',
'#FF5A0A', '#F02800', '#780A00', '#140A00']
#注意這里要對經度做變換,原來東半球經度在180~360區間,現在減去180,轉為0~180
long=nc.variables['longitude'][506:630]-180
#lati= np.flip(nc.variables['latitude'][72:172])
lati= nc.variables['latitude'][62:152] #[60:164]
#datai=np.flipud(data[30,72:172,506:630])-274 #轉換為攝氏度
datai=data[30,62:152,506:630]-274 #溫度數據切片,選擇第30個時間點的溫度;將原溫度轉換為攝氏度
graph(long,lati,datai,level_Tair,colors,title)
(90, 124) (90, 124) (90, 124)
全球溫度作圖
#全球溫度作圖,設定graph的參數
title='2m_temperature'
level_Tair= [-20,-10,0,5,10,15,20,25,30,1000] #[0,2.6,5,8,16,50,100,120,1000]
#[0,2.6,5,8,10,20,25,300,1000] [0,210,225,240,255,260,300,305,310,1000]
colors = ['#FFFFFF', '#AAF0FF', '#C8DC32', '#FFBE14', '#FF780A',
'#FF5A0A', '#F02800', '#780A00', '#140A00']
#注意這里要對經度做變換,原來東半球經度在180~360區間,現在減去180,轉為0~180
long=nc.variables['longitude'][:]-180
#lati= np.flip(nc.variables['latitude'][72:172])
lati= nc.variables['latitude'][:] #[60:164]
#datai=np.flipud(data[30,72:172,506:630])-274 #轉換為攝氏度
datai=data[30,:,:]-274 #將原溫度轉換為攝氏度
graph(long,lati,datai,level_Tair,colors,title)
(361, 720) (361, 720) (361, 720)
總結
以上是生活随笔為你收集整理的python处理nc气象数据_气象数据处理——nc文件的全部內容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: 多个DataSet数据合并
- 下一篇: python基础版课件_Python入门