ecCode库根据具体气象参数获取经纬度以及参数值

从github下载的ecCodes源码编译安装后，根据其提供的实例代码，调用接口，根据某一具体的气象参数来获取所有格网点的经纬度以及参数值。

具体代码如下:

#include <stdio.h>
#include <stdlib.h>

#include <eccodes.h>

#define MAX_VAL_LEN 1024 

int main(int argc, char** argv)
{
    int err              = 0;
    size_t i             = 0;
    FILE* in             = NULL;

    const char* filename = "/home/he/grib_data/gfs_data/gfs.t06z.sfluxgrbf103.grib2";

    codes_handle* h      = NULL;
    long numberOfPoints  = 0;
    const double missing = 9999.0;
    double *lats, *lons, *values; /* arrays */


    long paramId;
    double value;
    size_t vlen = MAX_VAL_LEN;
    codes_index* index;

    in = fopen(filename, "rb");
    if (!in) {
        fprintf(stderr, "Error: unable to open input file %s\n", filename);
        return 1;
    }

    /* create new handle from a message in a file */
    index = codes_index_new_from_file(0, filename, "paramId", &err);
    /* paramId is replaceable, such as 172 151... */
    CODES_CHECK(codes_index_select_string(index, "paramId", "156"), 0);
    h = codes_handle_new_from_index(index, &err);//biaoji

    if (h == NULL) {
        fprintf(stderr, "Error: unable to create handle from file %s\n", filename);
        return 1;
    }

    CODES_CHECK(codes_get_long(h, "numberOfPoints", &numberOfPoints), 0);
    CODES_CHECK(codes_set_double(h, "missingValue", missing), 0);

    lats = (double*)malloc(numberOfPoints * sizeof(double));
    if (!lats) {
        fprintf(stderr, "Error: unable to allocate %ld bytes\n", (long)(numberOfPoints * sizeof(double)));
        return 1;
    }
    lons = (double*)malloc(numberOfPoints * sizeof(double));
    if (!lons) {
        fprintf(stderr, "Error: unable to allocate %ld bytes\n", (long)(numberOfPoints * sizeof(double)));
        free(lats);
        return 1;
    }
    values = (double*)malloc(numberOfPoints * sizeof(double));
    if (!values) {
        fprintf(stderr, "Error: unable to allocate %ld bytes\n", (long)(numberOfPoints * sizeof(double)));
        free(lats);
        free(lons);
        return 1;
    }

    CODES_CHECK(codes_grib_get_data(h, lats, lons, values), 0);


    for (i = 0; i < numberOfPoints; ++i) {
        if (values[i] != missing) {
            printf("%f %f %f\n", lats[i], lons[i], values[i]);
            CODES_CHECK(codes_get_long(h, "paramId", &paramId), 0);
            printf("paramId=%ld\n", paramId);//biaoji
        }
    }

    free(lats);
    free(lons);
    free(values);
    codes_handle_delete(h);

    fclose(in);
    return 0;
}

在以上代码中，可以根据自己需要的气象参数idECMWF | Parameter database，此处可以找到各个气象参数在ecCodes中对应的ID，然后改变ID即可获取自己需要的各个格网点的经纬度以及气象参数值，值得注意的是，不同的grib文件，包含的参数不同，如要查看自己的grib文件中包含哪些参数，可以运行以下代码：

#include <stdlib.h>
#include <stdio.h>
#include <string.h>


#include "eccodes.h"

#define MAX_VAL_LEN 1024

static void usage(char* progname);

int main(int argc, char* argv[])
{
    /* to skip read only and computed keys
       unsigned long key_iterator_filter_flags=CODES_KEYS_ITERATOR_SKIP_READ_ONLY |
                                             CODES_KEYS_ITERATOR_SKIP_COMPUTED;
     */
    unsigned long key_iterator_filter_flags = CODES_KEYS_ITERATOR_ALL_KEYS |
                                              CODES_KEYS_ITERATOR_SKIP_DUPLICATES;

    /* choose a namespace. E.g. "ls", "time", "parameter", "geography", "statistics" */
    //const char* name_space = "ls";
    const char* name_space = "parameter";

    /* name_space=NULL to get all the keys */
    /* char* name_space=0; */

    //FILE* f         = NULL;

    FILE* f = NULL;  

    const char* filename = "/home/he/grib_data/gfs_data/gfs.t06z.sfluxgrbf103.grib2";

    f = fopen(filename, "rb");//biaoji

    codes_handle* h = NULL;

    int err       = 0;
    int msg_count = 0;

    char value[MAX_VAL_LEN];
    size_t vlen = MAX_VAL_LEN;


    while ((h = codes_handle_new_from_file(0, f, PRODUCT_GRIB, &err)) != NULL) {
        codes_keys_iterator* kiter = NULL;
        msg_count++;
        printf("-- GRIB N. %d --\n", msg_count);

        kiter = codes_keys_iterator_new(h, key_iterator_filter_flags, name_space);
        if (!kiter) {
            fprintf(stderr, "Error: Unable to create keys iterator\n");
            exit(1);
        }

        while (codes_keys_iterator_next(kiter)) {
            const char* name = codes_keys_iterator_get_name(kiter);
            vlen             = MAX_VAL_LEN;
            memset(value, 0, vlen);
            CODES_CHECK(codes_get_string(h, name, value, &vlen), name);
            printf("%s = %s\n", name, value);

            /* alternative way of getting the string value */
            CODES_CHECK(codes_keys_iterator_get_string(kiter, value, &vlen), 0);
        }

        codes_keys_iterator_delete(kiter);
        codes_handle_delete(h);
    }
    fclose(f);
    return 0;
}

static void usage(char* progname)
{
    printf("\nUsage: %s grib_file\n", progname);
    exit(1);
}

该段代码即github中ecCodes提供的示例代码，稍作修改，源码见

eccodes/grib_keys_iterator.c at develop · ecmwf/eccodes · GitHub

根据这个代码，便可以打印出自己的grib文件中包含的所有气象参数的具体信息。然后根据参数ID，运行最开始的代码即可获取气象数据信息。