GRIB API学习笔记07——GRIB API库2

目录

示例:grib_get
[code]
! Load all the GRIB messages contained in file.grib1
call grib_open_file(ifile, ‘file.grib1′,’r’)
n=1
call grib_new_from_file(ifile,igrib(n), iret)
LOOP: do while (iret /= GRIB_END_OF_FILE)
n=n+1; call grib_new_from_file(ifile,igrib(n), iret)
end do LOOP
! Decode/encode data from the loaded message
read, indx ! Choose one grib loaded GRIB message to decode
call grib_get( igrib(indx) , “dataDate”, date)
call grib_get(igrib(indx), “typeOfLevel”, typeOfLevel)
call grib_get(igrib(indx), “level”, level)
call grib_get_size(igrib(indx), “values”, nb_values); allocate(values(nb_values))
call grib_get(igrib(indx), “values”, values)
print
, date, levelType, level, values(1), values(nb_values)
! Release
do i=1,n
call grib_release(igrib(n))
end do
deallocate(values)
call grib_close_file(ifile)
[/code]
每条消息依次从文件中载入,在do while循环中反复判断是否到达文件结尾。

1.1.  GRIB API索引访问

子程序:
grib_index_create(indexid, filename, keys, status)
创建文件内容的索引(index)
grib_index_get_size(indexid, key, size, status)
获取索引中某个key的维数
grib_index_get(indexid, key, values, status)
获取索引中某key的不同值
to get the different “values” for a key in the index
grib_index_select(indexid, key, value, status)
选择索引中某key的一个值
grib_new_from_index(indexid, igrib, status)

to load the GRIB message corresponding to the selection made.
grib_index_release(indexid, status)
释放索引
grib_release(igrib)
 
使用索引方式随机访问通常比序列方式快很多。
 
例:
[code]
! create an index from a grib file using two keys call
grib_index_create(idx,‘ensemble.grib’,’paramId’)
[/code]
 
‘ensemble.grib’:文件ensemble.grib包含所有集合成员的参数
‘paramId’:可以索引一组key,用逗号分隔,中间没有空格。
 
[code]
! get the number of distinct values of parameters in the index
call grib_index_get_size(idx,’paramId’,paramIdSize)
! allocate the array to contain the list of distinct paramId
allocate(paramId(paramIdSize))
! get the list of distinct parameters from the index
call grib_index_get(idx,’paramId’,paramId)
count=1
do i=1,paramIdSize ! loop on paramId
! select paramId=paramId(i)
call grib_index_select(idx,’paramId’,paramId(i))
call grib_new_from_index(idx,igrib,iret)
[/code]
‘paramId’为建立索引的所有key选择一个值。
grib_new_from_index:加载需要的第一条grib消息到内存。
 
[code]
do while (iret /= GRIB_END_OF_INDEX)
call grib_is_missing(igrib,’number’, is_missing);
if (is_missing /= 1) then
call grib_get(igrib,’number’,onumber)
else
onumber=-9999
end if
call grib_get(igrib,’level’,olevel)
print*,’param:”, paramId(i),’ level:’,olevel, ‘ number:’,onumber
call grib_release(igrib)
call grib_new_from_index(idx,igrib,iret)
end do
[/code]
do while (iret /= GRIB_END_OF_INDEX):
对于一个索引,可能有几条grib消息可用,因此需要使用循环。
 
[code]
end do ! loop on paramId
call grib_release(igrib)
call grib_index_release(idx)
[/code]

1.1.1.   i/o

索引可以保存到文件,重复利用。
grib_index_write(indexid, filename, status)
保存索引到文件
grib_index_read(indexid, filename, status)
从之前由grib_index_write创建的索引文件中载入索引
 
可以增加某个数据文件的内容到某个索引。
grib_index_add_file(indexid, filename, status)
增加数据文件的内容到索引
可以使用命令行工具grib_index_build创建索引。
更多的索引功能可能会被加入。

1.2.  编码一条载入的GRIB消息

要点:尽可能少“编码”。从来不需要“编码”整个GRIB消息。
主要的编码函数:

grib_set(igrib, keyname, values, status)

integer, intent(in) :: igrib

character(len=*), intent(in) :: keyname

,[dimension(:),] intent(in) :: values

integer, optional, intent(out) :: status
type是integer或singgle/double real precision或string
 
写入消息:
grib_write(igrib, output_file)
注意,由grib_write写入的消息保持格式上的正确,但可能在语义上不正确。
 

1.2.1.   创建一条新消息

从样例文件中创建一条新消息:
样例文件是保存在示例目录的示例grib消息。默认的示例目录可以用grib_info命令查询。样例文件以后缀”.tmpl”结尾。用户可以创建自己的样例文件,并修改或添加环境变量GRIB_SAMPLES_PATH指向它们。
 
grib_new_from_samples(igrib, samplename, status)
 
消息可以来自其它消息的克隆:
grib_clone(igrib_src,igrib_dest,status)
 
例grib_set
[code]
! STEP-1: open output file and load a GRIB message from a sample “GRIB1”
call grib_open_file(outfile, ‘out.grib1’,’w‘)
call grib_new_from_samples(igrib, “GRIB1”)
! STEP-2: Get some information from the loaded message
call grib_get_size(igrib, “values”, nb_values)
allocate(values(nb_values)) ! Declared as real, dimension(:), allocatable
call model(values); values(1:100) = 9999.0 ! Compute values and set some missing values
! STEP-3: set the new GRIB message
call grib_set(igrib,’missingValues’, 9999.0) ! Tells the GRIB-API 9999.0 is the missing value
call grib_set(igrib,’bitmapPresent’, 1)
call grib_set(igrib,”values”, values) ! Set values as 1D real array of size nb_values
! STEP-4: write modified message to a file
call grib_write(igrib,outfile)
call grib_release(igrib)
call grib_close_file(outfile)
deallocate(values)
[/code]
 

1.2.2.   使用模板创建GRIB消息

可以将一个模板分配给创建的GRIB消息(使用样例或克隆)。
GRIB1:
下面是grid type和packing type的定义:
<//www.ecmwf.int/publications/manuals/d/gribapi/fm92/grib1/show/grids/>
<//www.ecmwf.int/publications/manuals/d/gribapi/fm92/grib1/show/packing/>
修改dataRepresentationType和packingType应用一个定义:
call grib_set(igrib,’dataRepresentationType’, 5)
将为你的消息定义一个Polar Stereographic Projection Grid。
 
GRIB2:
模板用于定义网格几何形状(section 3,例如lat/lon),产品定义(section 4,例如分析场),数据表示(section 5,例如简单打包),参见:
<//www.ecmwf.int/publications/manuals/d/gribapi/fm92/grib2/show/templates/>
 
可以使用关键词gridDefinitionTemplateNumber、productDefinitionTemplateNumber和dataRepresentationTemplateNumber使用任意模板。例如:
call grib_set(igrib,’dataRepresentationTemplateNumber’, 20)
将为GIRB2消息定义一个Polar Stereographic Projection Grid。
 

1.2.3.   使用不同的打包方式

GRIB数据可以使用不同的方式打包:
参考
https://software.ecmwf.int/wiki/display/GRIB/Grib+API+keys/
packingType
The alghorithm used to pack data into the GRIB message.
For GRIB edition 1:
grid_simple
grid_simple_matrix
grid_simple_matrix_bitmap
grid_second_order
spectral_complex
spectral_simple
grid_unknown
spectral_unknown
For GRIB edition 2:
grid_simple
grid_simple_matrix
grid_simple_matrix_bitmap
grid_complex
grid_complex_spatial_differencing
grid_jpeg
grid_second_order
grid_png
grid_ieee
spectral_simple
spectral_complex
spectral_complex
grid_simple_log_preprocessing
打包类型或者对格点数据可用,或者对spectral field可用。
打包种类会影响GRIB消息的大小,例如second order packing比simple packing小两倍。
同时会影响打包或解包解码GRIB消息所花的时间,例如second order packing比simple packing慢好几倍。
打包不丢失信息,无损。
更多内容参见练习。

1.3.  C API

1.3.1.   索引

不需要使用fopen()或fclose()。
grib_index * grib_index_new_from_file(grib_context *c, char *filename, const char *keys, int *err)
从文件中建立索引
grib_context *c 通常被设置为0。
 
int grib_index_get_size(grib_index *index, const char *key, size_t *size)
获取index的参数中包含的key的不同值的个数
int grib_index_get_double(grib_index *index, const char *key, double *values, size_t *size)
获取index中包含参数的key的值。
必须实现申请大小为grib_index_get_size()的内存空间。
int grib_index_get_string(grib_index *index, const char *key, char **values, size_t *size)
获取index中包含key的值
需要预先声明size大小的char*数组
size包含字符串的实际大小
 
例子:
[cpp]
char** paramId=NULL;
GRIB_CHECK(grib_index_get_size(index,“paramId”,&paramIdSize),0);
paramId=(char**)malloc(sizeof(char*) * paramIdSize);
GRIB_CHECK(grib_index_get_string(index,“paramId”,paramId,&paramIdSize),0);
for (i=0; i<paramIDSize; i++) free(paramId[i]); free(paramId); [/cpp]   int grib_index_select_TYPE(grib_index *index, const char *key, TYPE value)
选择key==value的数据子集
grib_handle * grib_handle_new_from_index(grib_index *index, int *err)
选择键值后,从index中创建新的句柄
句柄使用后,需要调用grib_handle_delete() 释放内存
下一次调用grib_handle_new_from_index() 会创建指向index中的下一个grib消息的句柄。
int grib_index_add_file(grib_index *index, const char *filename)
增加另一个文件到已经存在的index中

1.3.2.   编码

int grib_set_double(grib_handle *h, const char *key, double value)
设置double值,对于long有相同的函数。
int grib_set_string(grib_handle *h, const char *key, const char *mesg, size_t *length)
设置字符串值,对于bytes 有相似的函数。
int grib_set_double_array(grib_handle *h, const char *key, const double *vals, size_t length)
设置双精度数组,对于long array 有相似的函数。

1.3.3.   克隆

  1. 为现有GRIB消息创建句柄
  2. grib_handle * grib_handle_clone(grib_handle *h)

n  使用原有句柄的内容克隆已存在的句柄,消息被复制且重新解析。

  1. 编码(覆写)新grib消息中的keys
  2. int grib_get_message(grib_handle *h_new, const void **message, size_t *message_length)

n  获取句柄关联的grib消息的原始二进制格式。

  1. 打开输出文件: out = fopen(file,”w”)
  2. 写入消息: fwrite(message,1, message_length,out)
  3. 关闭文件: fclose(out)

1.4.  Python API

1.4.1.   索引

iid = grib_index_new_from_file(file, keys)
返回创建index的句柄
使用grib_index_release(iid) 释放
例:
[python]
iid=gribapi.grib_index_new_from_file(“uv500.grib2”,[“shortName”])
[/python]
 
file是文件的路径,keys是key的元组或列表。
 
grib_index_add_file(iid, file)
增加一个文件到索引。
例:
[python]
>>>gribapi.grib_index_add_file(iid,”data/rmf.gra.2013103112001.grb2”)
>>>gribapi.grib_index_get_size(iid,”shortName”)
18
[/python]
 
grib_index_write(iid, file)
将索引写入文件便于重用
 
iid = grib_index_read(file)
从文件中载入由grib_index_write() 写入到文件的索引。
 
例:
[python]
>>> iid=gribapi.grib_index_read(“grib_index_file.01″)
>>> gribapi.grib_index_get_size(iid,”shortName”)
18
[/python]
size = grib_index_get_size(iid, key)
返回索引key的不同值个数
 
values = grib_index_get(iid, key, type=str)
获取索引key的不同值
例:
[python]
>>> gribapi.grib_index_get(iid,”shortName”)
(’10u’, ’10v’, ‘2t’, ‘d’, ‘gh’, ‘kx’, ‘papt’, ‘pli’, ‘prmsl’, ‘q’, ‘r’, ‘sp’, ‘sx’, ‘t’, ‘u’, ‘unknown’, ‘v’, ‘vo’)
[/python]
grib_index_select(iid, key, value)
选择key==value的子集
例:
[python]
>>> gribapi.grib_index_select(iid,”shortName”,”q”)
[/python]
 
gid = grib_new_from_index(iid)
与grib_new_from_file相同
用grib_release(gid) 释放
例:上例中,我们选择q
[python]
>>> gid=gribapi.grib_new_from_index(iid)
>>> gribapi.grib_get(gid,”shortName”)
‘q’
[/python]
 

1.4.2.   编码

grib_set(gid, key, value)
n  设置grib消息标量key的值
grib_set_array(gid, key, value)
n  设置grib消息的数组值
n  输入数组可以是numpy.ndarray或Python序列,如元组、列表、数组等
grib_set_values(gid, values)
n  设置values(key名为values)的值。

1.4.3.   克隆

 
clone_id = grib_clone(gid_src)
n  创建一个消息的拷贝
n  直接用grib_write 写入文件
n  不要忘了使用grib_release

1.4.4.   实用函数

[outlat, outlon, value, distance, index] = grib_find_nearest(gid, inlat, inlon, is_lsm=False, npoints=1)
n  找到指定经纬度点的最近邻
n  npoints=4返回四个最近邻点的列表
 
iter_id = grib_iterator_new(gid,mode)
n  [lat,lon,value] = grib_iterator_next(iterid)
grib_iterator_delete(iter_id)