|
最近在网上下载的NCEP的fnl数据，发现只有grib2文件格式。因为现在处理数据使用的更多的是python，matlab一类的。网上一搜python读取grib文件需要依赖pygrib库，pygrib是欧洲中期天气预报中心(ECMWF)的GRIG API C库的Python接口，通过这个库就可以将Grib数据读取出来。但这个依赖库只有Linux、Mac OS系统下可以使用，Windows系统没有这个依赖库，因此Windows下python不能直接读取grib文件。但python可以处理nc文件，想办法把grib文件转化为nc文件，就可以方便处理了。先说方法，主要有两种方案，一种是安装Cygwin，可以直接处理。另一种是安装wgrib2，转化数据。
1.读取grib2的两种方案
1.1下载安装Cygwin
下载安装Cygwin，Cygwin是一个在windows平台上运行的类UNIX模拟环境，就等于是一个Linux的模拟器，这样就可以下载pygrib依赖库处理grib2文件，本人不喜欢用这种模拟器，众所周知，模拟器容易卡卡卡。因此，主要讲述第二种方法。
1.2下载安装wgrib2
因为最新的数据都是grib2文件了，所以直接下载wgrib2。（如果需要处理的是grib文件，则下载wgrib） 下载wgrib2链接：ftp://ftp.cpc.ncep.noaa.gov/wd51we/wgrib2/Windows_64/这里直接是Windows 64的链接，如果需要别的版本可以下载对应的版本。(最近有人反映网站改成ftp链接打开了，所以我把我下载的上传成资源，wgrib2相关安装文件.rar)。 主要下载的东西：  你可以先下载wgrib2.exe，运行的时候会提示缺少啥，对应下载啥。 下载安装成功后，在cmd下对应的wgrib2文件夹下运行wgrib2.exe，就会看到对应得数据参数选项。  对应得所有功能如下，只是一个展示，不需要细看：
wgrib2 v0.2.0.9b1 3/2019 Wesley Ebisuzaki, Reinoud Bokhorst, John Howard, Jaakko Hyvätti, Dusan Jovic, Daniel Lee, Kristian Nilssen, Karl Pfeiffer, Pablo Romero, Manfred Schwarb, Gregor Schee, Arlindo da Silva, Niklas Sondell, Sam Trahan, Sergey Varlamov
stock build
-else else else
-elseif_fs elif X elseif X (fixed string) conditional execution
-elseif_n elif X elseif (inv numbers in range), X=(start:end:step)
-elseif_rec elif X elseif (record numbers in range), X=(start:end:step)
-elseif_reg elif X elseif rpn registers defined, X = A, A:B, A:B:C, etc A = register number
-elseif elif X elseif X (POSIX regular expression) conditional on match
-endif endif endif
-if if X if X (POSIX regular expression), conditional execution on match
-if_fs if X if X (fixed string), conditional execution on match
-if_n if X if (inv numbers in range), X=(start:end:step)
-if_rec if X if (record numbers in range), X=(start:end:step)
-if_reg if X if rpn registers defined, X = A, A:B, A:B:C, etc A = register number
-not_if if X not_if X (regular expression), conditional execution on not match
-not_if_fs if X if X (fixed string) does not match, conditional execution up to next output/fi
-0xSec inv X Hex dump of section X (0..8)
-aerosol_size inv optical properties of an aerosol
-aerosol_wavelength inv optical properties of an aerosol
-bitmap inv bitmap mode
-center inv center
-checksum inv X CRC checksum of section X (0..8), whole message (X = -1/message) or (X=data)
-disc inv discipline (code table 0.0)
-domain inv find rectangular domain for g2ctl/GrADS plots
-end_ft inv verf time = reference_time + forecast_time + stat. proc time (YYYYMMDDHH) (same as -vt)
-end_FT inv verf time = reference_time + forecast_time + stat. proc time (YYYYMMDDHHMMSS) (same as -VT)
-ens inv ensemble information
-ext_name inv extended name, var+qualifiers
-ftime inv either ftime1 or ftime2 dep on version_ftime
-ftime1 inv forecast time
-ftime2 inv timestamp -- will replace -ftime in the future TESTING
-ftn_api_fn0 inv n npnts nx ny msg_no submsg i11,5(1x,i11)
-full_name inv extended name, var+qualifiers
-gdt inv contents of Grid Definition Template (g2c)
-get_byte inv X Y Z get bytes in Section X, Octet Y, number of bytes Z (decimal format)
-get_hex inv X Y Z get bytes in Section X, Octet Y, number of bytes Z (bytes in hexadecimal format)
-get_ieee inv X Y Z get ieee float in Section X, Octet Y, number of floats Z
-get_int inv X Y Z get 4-byte ints in Section X, Octet Y, number of ints Z
-get_int2 inv X Y Z get 2-byte ints in Section X, Octet Y, number of ints Z
-grib_max_bits inv maximum bits used in grib encoding
-grid inv grid definition
-grid_id inv show values from grid_id
-hybrid inv shows vertical coordinate parameters from Sec4
-ij inv X Y value of field at grid(X,Y) X=1,..,nx Y=1,..,ny (WxText enabled)
-ijlat inv X Y lat,lon and grid value at grid(X,Y) X=1,..,nx Y=1,..,ny (WxText enabled)
-ilat inv X lat,lon and grid value at Xth grid point, X=1,..,npnts (WxText enabled)
-JMA inv inventory for JMA locally defined PDT
-lev inv level (code table 4.5)
-ll2i inv X Y x=lon y=lat, converts to (i), 1..ndata
-ll2ij inv X Y x=lon y=lat, converts lon-lat to (i,j) using gctpc
-lon inv X Y value at grid point nearest lon=X lat=Y (WxText enabled)
-match_inv inv inventory used by -match, -not, -if and -not_if
-Match_inv inv same as -match_inv except d=YYYYMMDDHH D=YYYYMMDDHHmmss
-max inv print maximum value
-min inv print minimum value
-misc inv variable name qualifiers like chemical, ensemble, probability, etc
-MM inv reference time MM
-model_version_date inv prints model date code
-n inv prints out inventory number
-N_ens inv number of ensemble members
-nl inv inserts new line into inventory
-nlons inv number of longitudes for each latitude
-npts inv number of grid points
-nxny inv nx and ny of grid
-packing inv shows the packing mode (use -v for more details)
-pdt inv Product Definition Table (Code Table 4.0)
-precision inv precision of packing
-print inv X inserts string (X) into inventory
-prob inv probability information
-process inv Process (code table 4.3)
-processid inv process id (locally defined)
-proj4_ij2ll inv X Y X=x Y=y, converts to (i,j) to lon-lat using proj.4 (experimental) we:sn
-proj4_ll2i inv X Y x=lon y=lat, converts to (i) using proj.4 (experimental) 1..ndata
-proj4_ll2ij inv X Y x=lon y=lat, converts lon-lat (i,j) using proj.4 (experimental)
-radius inv radius of Earth
-range inv print out location of record in bytes, 0 = first byte
-RT inv type of reference Time
-s inv simple inventory
-S inv simple inventory with minutes and seconds (subject to change)
-s2 inv simple inventory .. for testing ftime2
-scale inv scale for packing
-scaling inv scaling for packing (old format)
-scan inv scan order of grid
-Sec_len inv length of various grib sections
-Sec0 inv contents of section0
-Sec3 inv contents of section 3 (Grid Definition Section)
-Sec4 inv Sec 4 values (Product definition section)
-Sec5 inv Sec 5 values (Data representation section)
-Sec6 inv show bit-map section
-spatial_proc inv show spacial processing, pdt=4.15
-spectral_bands inv spectral bands for satellite, pdt=4.31 or 4.32
-start_ft inv verf time = reference_time + forecast_time (YYYYMMDDHH) : no stat. proc time
-start_FT inv verf time = reference_time + forecast_time (YYYYMMDDHHMMSS) - no stat. proc time
-stats inv statistical summary of data values
-subcenter inv subcenter
-t inv reference time YYYYMMDDHH, -v2 for alt format
-T inv reference time YYYYMMDDHHMMSS
-table inv parameter table
-timer inv reads OpenMP timer
-unix_time inv print unix timestamp for rt & vt
-V inv diagnostic output
-var inv short variable name
-varX inv raw variable name - discipline mastertab localtab center parmcat parmnum
-vector_dir inv grid or earth relative winds
-verf inv simple inventory using verification time
-vt inv verf time = reference_time + forecast_time, -v2 for alt format
-VT inv verf time = reference_time + forecast_time (YYYYMMDDHHMMSS)
-warn_old_g2 inv warn if old g2lib would have problem
-wave_partition inv ocean surface wave partition (pdt=4.52)
-YY inv reference time YYYY
-inv_f77 inv> X Y Z match inventory written to Z with character*(Y) and X=(bin,ieee)
-last inv> X write last inv item to file X
-last0 inv> X write last inv item to beginning of file X
-nl_out inv> X write new line in file X
-print_out inv> X Y prints string (X) in file (Y)
-s_out inv> X simple inventory written to X
-big_endian misc sets ieee output to big endian (default is big endian)
-box_ave misc X Y Z box average X=odd integer (lon) Y=odd integer (lat) critical_weight
-colon misc X replace item deliminator (:) with X
-config misc shows the configuration
-count misc prints count, number times this -count was processed
-end misc stop after first (sub)message (save time)
-error_final misc X Y Z error if at end X=count Y=ne,eq,le,lt,gt,ge Z=integer
-fix_CFSv2_fcst misc X Y Z fixes CFSv2 monthly fcst X=daily or 00/06/12/18 Y=pert no. Z=number ens fcsts v1.0
-fix_ncep misc fix ncep PDT=8 headers produced by cnvgrib
-gctpc misc X X=0,1 use gctpc library (default=1)
-grid_changes misc prints number of grid changes
-grid_def misc read lon and lat data from grib file -- experimental
-h misc help, shows common options
-header misc f77 header or nx-ny header in text output (default)
-help misc X help [search string|all], -help all, shows all options
-ijundefine misc X Y Z sets grid point values to undefined X=(in-box|out-box) Y=ix0:ix1 Z=iy0:iy1 ix=(1..nx) iy=(1..ny)
-import_bin misc X read binary file (X) for data
-import_grib misc X read grib2 file (X) for data
-import_grib_fs misc X Y read grib2 file (Y) sequentially for record that matches X (fixed string)
-import_ieee misc X read ieee file (X) for data
-import_netcdf misc X Y Z TESTING X=file Y=var Z=hyper-cube specification
-import_text misc X read text file (X) for data
-limit misc X stops after X fields decoded
-little_endian misc sets ieee output to little endian (default is big endian)
-mem_del misc X delete mem file X
-mem_final misc X Y write mem file X to file Y at cleanup step
-mem_init misc X Y read mem file X from file Y (on initialization)
-ndate misc X Y X=date Y=dt print date + dt
-ndates misc X Y Z X=date0 Y=(date1|dt1) Z=dt2 for (date=date0; date 11, X=code table 4.6 Y=pert num Z=num ens members -1=No Change
-set_ensm_derived_fcst misc X Y convert PDT 0,1,2 -> 2, 8,11,12 -> 12, X=code table 4.7 Y=num ens members
-set_ftime misc X either set_ftime1 or set_ftime2 dep on version_ftime
-set_ftime1 misc X set ftime
-set_ftime2 misc X set ftime2 .. will be replace -set_ftime/ave in the future -- TESTING ---
-set_grib_max_bits misc X sets scaling so number of bits does not exceed N in (new) grib output
-set_grib_type misc X set grib type = jpeg, simple, ieee, complex(1|2|3), aec, same
-set_hex misc X Y Z set bytes in Section X, Octet Y, bytes Z (a|a:b:c|abc) in hexadecimal
-set_ieee misc X Y Z set ieee float in Section X, Octet Y, floats Z (a|a:b:c)
-set_ijval misc X Y Z sets grid point value X=ix Y=iy Z=val
-set_int misc X Y Z set 4-byte ints in Section X, Octet Y, signed integers Z (a|a:b:c)
-set_int2 misc X Y Z set 2-byte ints in Section X, Octet Y, signed integers Z (a|a:b:c)
-set_ival misc X Y sets grid point value X=i1:i2:.. Y=va1:val2:.. grid[i1] = val1,etc
-set_lev misc X changes level code .. not complete
-set_metadata misc X read meta-data for grib writing from file X
-set_metadata_str misc X X = metadata string
-set_pdt misc X makes new pdt, X=(+)PDT_number or X=(+)PDT_number:size of PDT in octets, +=copy metadata
-set_percentile misc X convert PDT 0..6 -> 6, 8..15 -> 10, X=percentile (0..100)
-set_prob misc 5 args X/Y forecasts Z=Code Table 4.9 A=lower limit B=upper limit
-set_radius misc X set radius of Earth X= 0,2,4,5,6,8,9 (Code Table 3.2), 1:radius , 7:major:minor
-set_scaling misc X Y set decimal scaling=X/same binary scaling=Y/same new grib messages
-set_sec_size misc X Y resizes section , X=section number, Y=size in octets, DANGEROUS
-set_ts_dates misc X Y Z changes date code for time series X=YYYYMMDDHH(mmss) Y=dtime Z=#msgs/date
-set_var misc X changes variable name
-start_timer misc starts OpenMP timer
-status misc X X X=file
-submsg misc X process submessage X (0=process all messages)
-sys misc X run system/shell command, X=shell command
-text_col misc X number of columns on text output
-text_fmt misc X format for text output (C)
-udf misc X Y run UDF, X=program+optional_args, Y=return file
-udf_arg misc X Y add grib-data to UDF argument file, X=file Y=name
-undefine misc X Y Z sets grid point values to undefined X=(in-box|out-box) Y=lon0:lon1 Z=lat0:lat1
-undefine_val misc X grid point set to undefined if X=val or X=low:high
-v misc verbose (v=1)
-v0 misc not verbose (v=0)
-v2 misc really verbose (v=2)
-version misc print version
--version misc print version
-AAIG out writes Ascii ArcInfo Grid file, lat-lon grid only (alpha)
-AAIGlong out writes Ascii ArcInfo Grid file, lat-lon grid only long-name *.asc (alpha)
-ave out X Y average X=time step Y=output v2
-ave_var out X Y average/std dev/min/max X=time step, Y=output
-ave0 out X Y average X=time step, Y=output grib file needs file is special order
-bin out X write binary data to X
-cress_lola out X..Z,A lon-lat grid values X=lon0:nlon:dlon Y=lat0:nlat:dlat Z=file A=radius1:radius2:..:radiusN
-csv out X make comma separated file, X=file (WxText enabled)
-csv_long out X make comma separated file, X=file (WxText enabled)
-cubeface2global out X Y write faces X as global cubed grid to Y: X=list of faces to exclude
-ens_processing out X Y ave/min/max/spread X=output Y=future use
-fcst_ave out X Y average X=time step Y=output v2
-fcst_ave0 out X Y average X=time step, Y=output grib file needs file is special order
-fi out depreceated
-grib out X writes GRIB record (one submessage) to X
-GRIB out X writes entire GRIB record (all submessages)
-grib_ieee out X writes data[] to X.grb, X.head, X.tail, and X.h
-grib_out out X writes decoded/modified data in grib-2 format to file X
-grib_out_irr out X Y writes irregular grid grib (GDT=130 not adopted) X=(all|defined) Y=(output file)
-grib_out_irr2 out 5 args writes irregular grid grib GDT 101 X=npnts Y=grid_no Z=grid_ref A=UUID B=(output file)
-gridout out X text file with grid: i j lat lon (1st record)
-ieee out X write (default:big-endian) IEEE data to X
-ijbox out X..Z,A grid values in bounding box X=i1:i2[:di] Y=j1:j2[:dj] Z=file A=[bin|text|spread]
-ijsmall_grib out X Y Z make small domain grib file X=ix0:ix1 Y=iy0:iy1 Z=file
-irr_grid out X Y Z make irregular grid (GDT=130 not adopted), nearest neighbor, X=lon-lat list Y=radius (km) Z=output grib file
-lola out X..Z,A lon-lat grid values X=lon0:nlon:dlon Y=lat0:nlat:dlat Z=file A=[bin|text|spread|grib]
-merge_fcst out X Y merge forecast ave/acc/min/max X=number to intervals to merge (0=every) Y=output grib file
-mysql out 5 args H=[host] U=[user] P=[password] D=[db] T=[table]
-mysql_dump out 7 args H=[host] U=[user] P=[password] D=[db] T=[table] W=[western_lons:0|1] PV=[remove unlikely:0|1]
-mysql_speed out 7 args H=[host] U=[user] P=[password] D=[db] T=[table] W=[western_lons:0|1] PV=[remove unlikely:0|1]
-ncep_norm out X normalize NCEP-type ave/acc X=output grib file
-ncep_uv out X combine U and V fields into one message like NCEP operations
-netcdf out X write netcdf data to X
-new_grid out X..Z,A bilinear interpolate: X=projection Y=x0:nx:dx Z=y0:ny:dy A=grib_file alpha
-reduced_gaussian_grid out X Y Z reduced Gaussian grid, X=outputfile Y=-1 Z=(neighbor|linear)[-extrapolate]
-small_grib out X Y Z make small domain grib file X=lonW:lonE Y=latS:latN Z=file
-spread out X write text - spread sheet format into X (WxText enabled)
-submsg_uv out X combine vector fields into one message
-text out X write text data into X
-time_processing out X..Z,A average X=CodeTable 4.10 Y=CodeTable 4.11 Z=time step A=output
-tosubmsg out X convert GRIB message to submessage and write to file X
-wind_dir out X calculate wind direction, X = output gribfile (direction in degrees, 0=wind from north, 90=wind from east)
-wind_speed out X calculate wind speed, X = output gribfile (U then V in datafile)
-wind_uv out X calculate UGRD/VGRD from speed/dir, X = output gribfile
-alarm init X terminate after X seconds
-append init append mode, write to existing output files
-crlf init make the end of the inventory a crlf (windows) instead of newline (unix)
-d init X dump message X (n or n.m), only 1 -d allowed
-egrep init X egrep X | wgrib2 (X is POSIX regular expression)
-egrep_v init X egrep -v X | wgrib2 (X is POSIX regular expression)
-eof_bin init X Y send (binary) integer to file upon EOF: X=file Y=integer
-eof_string init X Y send string to file upon EOF: X=file Y=string
-err_bin init X Y send (binary) integer to file upon err exit: X=file Y=integer
-err_string init X Y send string to file upon err exit: X=file Y=string
-fgrep init X fgrep X | wgrib2
-fgrep_v init X fgrep -v X | wgrib2
-fix_ncep_2 init ncep bug fix 2, probability observation < -ve number
-fix_ncep_3 init sets flag to fix ncep bug 3 (constant fields)
-fix_ncep_4 init fixes NCEP grib2 files where DX and DY are undefined
-for init X process record numbers in range, X=(start:end:step), only one -for allowed
-for_n init X process inv numbers in range, X=(start:end:step), only one -for allowed
-g2clib init X X=0/1/2 0=WMO std 1=emulate g2clib 2=use g2clib
-i init read Inventory from stdin
-i_file init X read Inventory from file
-inv init X write inventory to X
-match init X process data that matches X (POSIX regular expression)
-match_fs init X process data that matches X (fixed string)
-match_inv_add init X Y Z add new options to match_inventory
-nc_grads init require netcdf file to be grads v1.9b4 compatible (fixed time step only)
-nc_nlev init X netcdf, X = max LEV dimension for {TIME,LEV,LAT,LON} data
-nc_pack init X pack/check limits of all NEW input variables, X=min:max[:byte|short|float]
-nc_table init X X is conversion_to_netcdf_table file name
-nc_time init X netcdf, [[-]yyyymmddhhnnss]:[dt{s[ec]|m[in]|h[our]|d[ay]}], [-] is for time alignment only
-nc3 init use netcdf3 (classic)
-nc4 init use netcdf4 (compressed, controlled endianness etc)
-ncpu init X number of threads, default is environment variable OMP_NUM_THREADS/number of cpus
-no_append init not append mode, write to new output files (default)
-no_nc_grads init netcdf file may be not grads v1.9b4 compatible, variable time step
-no_nc_pack init no packing in netcdf for NEW variables
-no_nc_table init disable previously defined conversion_to_netcdf_table
-no_nc_time init netcdf, disable previously defined initial or relative date and time step
-not init X process data that does not match X (POSIX regular expression)
-not_fs init X process data that does not match X (fixed string)
-one_line init puts all on one line (makes into inventory format)
-order init X decoded data in X (raw|we:sn|we:ns) order, we:sn is default
-persistent init X makes file X persistent if already opened (default on open), CW2
-rewind_init init X rewinds file X on initialization if already opened, CW2
-set_ext_name init X X=0/1 extended name on/off
-set_regex init X set regex mode X = 0:extended regex (default) 1:pattern 2:extended regex & quote metacharacters
-set_version_ftime init X set version of ftime X=1, 2
-tigge init use modified-TIGGE grib table
-transient init X make file X transient, CW2
然后根据网上的两个例子，去使用wgrib2读取grib2数据。
1.2.1按编号读取
import os
# 将当前目录转到wgrib2工具目录
os.chdir("C:\\wgrib2")
# 执行cmd命令，查看各个“小房间”具体信息
#os.system(r"wgrib2 C:\\download\\fnl_20200327_00_00.grib2 -v")
结果：  可以看到数据信息，我们读第17号，温度信息，读到txt文件中。
os.system(r"wgrib2 C:\\download\\fnl_20200327_00_00.grib2 -d 17 -text C:\\download\\test.txt
结果如下：  这样python就可以直接处理txt文件。
1.2.2 按名称匹配读取
这种方法用的较多，因为气象变量的名称是不会变的，就不需要去找仓库号，这次把结果转成csv文件。
# 进行气象变量名称匹配,读取2m处温度这个气象变量到2m_tmp.csv文件中
#os.system("wgrib2 C:\\download\\fnl_20200327_00_00.grib2 -match ':TMP:2 m' -csv C:\\download\\2m_tmp.csv")
结果如下：  这样，我们就通过wgrib2将grib2数据转化为txt或者csv文件，这样就可以方便使用python处理。但是，这样我只能看部分数据，我想要的是能看到整个文件的数据，看网上关于wgrib2也没有能直接读取出来全部数据的。 最后看wgrib2的时候看到了它的一条命令：  这样，就可以通过命令把数据转化为nc文件，这样就可以用python直接读取。 下面看一下它的用法：
wgrib2 ../example/eta.t00z.awphys18.grb2 -netcdf eta.nc
具体的用法，可以详细了解wgrib2官网关于netcdf的用法https://www.cpc.ncep.noaa.gov/products/wesley/wgrib2/netcdf.html
os.system(r"wgrib2 C:\\download\\fnl_20200327_00_00.grib2 -netcdf C:\\download\\test.nc"
这样就将grib文件转化为nc文件，方便处理。
2.python处理nc文件
python处理nc文件，需要下载依赖库：
pip install netCDF4
简单的一些操作如下：
from netCDF4 import Dataset
import numpy as np
nc_obj = Dataset('C:\\wgrib2\\test.nc')
#查看nc文件
print(nc_obj)
print('---------------------------------------')
#查看nc文件中的变量
print(nc_obj.variables.keys())
for i in nc_obj.variables.keys():
print(i)
print('---------------------------------------')
# 查看nc文件中的变量
print(nc_obj.variables.keys())
print('---------------------------------------')
# 查看每个变量的信息
for i in nc_obj.variables.keys ():
print("---------------------------------")
print(nc_obj.variables[i])
print('---------------------------------------')
# 查看每个变量的所有属性
for i in nc_obj.variables.keys ():
print(nc_obj.variables[i].ncattrs())
print('---------------------------------------')
# 查看每个变量的单位
for i in nc_obj.variables.keys():
print(nc_obj.variables[i].units)
print('---------------------------------------')
# 查看每个变量的维数
for i in nc_obj.variables.keys():
print(nc_obj.variables[i].ndim)
print('---------------------------------------')
# 查看每个变量的形状（维度）
for i in nc_obj.variables.keys():
print(nc_obj.variables[i].shape)
print('---------------------------------------')
#取出变量blh 10 metre U wind component(10米U风分量)
u10 = nc_obj['TMP_550mb'][:]
print(u10)
print(np.max(u10), np.min(u10))
# 读取数据值
print("读取数据值")
for i in nc_obj.variables.keys():
print('--------------------------------------------------------------------------------------------------------------')
print(nc_obj.variables[i][:5])
只要转成nc文件了，就可以随便看需要了解的信息，目的就达到了，只要转化为python能处理的数据，就简单多了。
由于最近那个wgrib2网站无法打开，可以参考新的方法使用python读取grib文件，在Windows下使用xarray+cfgrib读取grib文件更加方便，不用转化nc文件，参考Windows下xarray+cfgrib读取grib文件。
|