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

December 31, 2013 (最后修改: November 19, 2019)

示例：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消息。
主要的编码函数：

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 有相似的函数。