科学数据存储

科学数据存储

有许多存储数据的方法，特别是以数组形式出现的数据。有很多人在电子表格中存储数据，然后将它们导出到带有逗号或制表符的 ascii 文件中，并希望其他人（或由他们自己编写的其他程序）再次读入这些数据，这样一个过程在以下几个方面的浪费：

  • 一个数字的 ascii 表示法要比内部的二进制表示法占用更多的空间，理想情况下，你希望一个文件和内存中的表示一样紧凑。
  • ascii 的转换很慢；它也可能导致精度的损失。

由于这样的原因，最好能有一种基于二进制存储的文件格式。还有一些对一个有用的文件格式的要求：

  • 由于二进制存储在不同的平台之间会有差异，所以一个好的文件格式是与平台无关的。例如，这将防止混淆大端字节序和小端字节序，以及32位和64位浮点数的约定。
  • 应用程序数据可能是异质的，包括整数、字符和浮点数据，理想情况下，所有这些数据应该被存储在一起。
  • 应用数据也是有结构的，这种结构应该反映在存储的形式中。
  • 一个文件格式最好是自带文档的，如果文件本身能告诉你哪些存储的数字是矩阵，哪些是向量，以及这些对象的大小是什么难道不是一件很好的事吗？

本教程将介绍 HDF5 库，它可以满足这些要求。HDF5 是一个庞大而复杂的库，所以本教程将只涉及到基础知识。欲知更多信息，请查阅http://www.hdfgroup.org/HDF5/ 当你做这个教程的时候，请保持你的浏览器在http://www.hdfgroup.org/HDF5/doc/ 或 http://www.hdfgroup.org/HDF5/RM/RM_H5Front.html 来了解这些例程的确切语法。

HDF5 简介

如上所述，HDF5 是一种与机器无关的文件格式，并且是自我记录的。每个 HDF5 文件的设置就像一个目录树，有子目录和包含实际数据的叶子结点。这意味着可以通过参考文件的名称，而不是文件中的位置来找到数据。本节中，你将学习如何编写向 HDF5 文件写入和从 HDF5 文件中读取的程序。为了检查文件是否符合你的目的，你可以在命令行中使用 h5dump 工具。

HDF5 格式与旧版本的 HDF4 并不兼容，后者已不再开发。由于历史原因，你仍然可以碰到有人使用 HDF4 。本教程是基于 HDF5 的1.6版本。在当前的1.8版本中，一些接口发生了变化；为了在1.8版本的软件中使用1.6版本的API，请在你的编译行中添加一个 -DH5_USE_16_API 标志。

许多 HDF5 例程是关于创建对象的：文件句柄、数据集中的成员等等。一般语法是：

  1.         hid_t h_id;
  2.         h_id = H5Xsomething(...);

创建对象的失败是由一个负的返回参数表示的，所以最好是创建一个包含 myh5defs.h 的文件。

  1.         #include "hdf5.h"
  2.         #define H5REPORT(e) \
  3.             {if (e<0) {printf("\nHDF5 error on line %d\n\n",__LINE__); \ 
  4.             return e;}}

并将此作为：

  1.         #include "myh5defs.h"
  3.         hid_t h_id;
  4.         h_id = H5Xsomething(...); H5REPORT(h_id);

创建一个文件

首先，我们需要创建一个 HDF5 文件。

  1.         hid_t file_id;
  2.         herr_t status;
  4.         file_id = H5Fcreate( filename, ... );
  5.             ...
  6.         status = H5Fclose(file_id);

这个文件将是一些数据项的容器，像目录树一样组织。

练习：通过编译和运行下面的 create.c 例子创建一个 HDF5 文件。

预期结果：应该会创建一个 file.h5 文件。

注意事项：请确保添加 HDF5 的包含目录和库目录。

  1.         cc -c create.c -I. -I/opt/local/include

然后

  1.         cc -o create create.o -L/opt/local/lib -lhdf5.

include 和 lib 目录将取决于系统。

  1.         /*
  2.         * File: create.c
  3.         * Author: Victor Eijkhout
  4.         */
  5.         #include "myh5defs.h"
  6.         #define FILE "file.h5"
  8.         main() {
  10.             hid_t file_id; /* file identifier */
  11.             herr_t status;
  13.             /* Create a new file using default properties. */
  14.             file_id = H5Fcreate(FILE, H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT);
  15.             H5REPORT(file_id);
  17.             /* Terminate access to the file. */
  18.             status = H5Fclose(file_id);
  19.         }

你可以在命令行上显示创建的文件：

  1.         %% h5dump file.h5
  2.         HDF5 "file.h5" {
  3.         GROUP "/" {
  4.         }
  5.         }

请注意，一个空文件只对应于将容纳数据的目录树的根。

数据集

接下来我们创建一个数据集，在这个例子中是一个二维网格。为了描述这个，我们首先需要构建一个数据空间：

  1.         dims[0] = 4; dims[1] = 6;
  2.         dataspace_id = H5Screate_simple(2, dims, NULL);
  3.         dataset_id = H5Dcreate(file_id, "/dset", dataspace_id, .... );
  4.         ....
  5.         status = H5Dclose(dataset_id);        
  6.         status = H5Sclose(dataspace_id);

注意，数据集和数据空间需要像文件一样被关闭。

练习：通过编译和运行下面的 dataset.c 代码创建一个数据集。

预期结果：这将创建一个 dset.h 文件，可以用 h5dump 显示。

  1.         /*
  2.         * File: dataset.c
  3.         * Author: Victor Eijkhout
  4.         */
  5.         #include "myh5defs.h"
  6.         #define FILE "dset.h5"
  8.         main() {
  10.         hid_t file_id, dataset_id, dataspace_id; /* identifiers */
  11.         hsize_t dims[2];
  12.         herr_t status;
  14.         /* Create a new file using default properties. */
  15.         file_id = H5Fcreate(FILE, H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT);
  17.         /* Create the data space for the dataset. */
  18.         dims[0] = 4;
  19.         dims[1] = 6;
  20.         dataspace_id = H5Screate_simple(2, dims, NULL);
  22.         /* Create the dataset. */
  23.         dataset_id = H5Dcreate(file_id, "/dset",H5T_NATIVE_INT,dataspace_id,H5P_DEFAULT);
  25.         /*H5T_STD_I32BE*/
  27.         /* End access to the dataset and release resources used by it. */
  28.         status = H5Dclose(dataset_id);
  30.         /* Terminate access to the data space. */
  31.         status = H5Sclose(dataspace_id);
  32.         /* Close the file. */
  33.         status = H5Fclose(file_id);
  34.     }

我们再次在线查看创建的文件：

  1.         %% h5dump dset.h5
  2.         HDF5 "dset.h5" {
  3.         GROUP "/" {
  4.         DATASET "dset" {
  5.             DATATYPE H5T_STD_I32BE
  6.             DATASPACE SIMPLE { ( 4, 6 ) / ( 4, 6 ) }
  7.             DATA {
  8.             (0,0): 0, 0, 0, 0, 0, 0,
  9.             (1,0): 0, 0, 0, 0, 0, 0,
  10.             (2,0): 0, 0, 0, 0, 0, 0,
  11.             (3,0): 0, 0, 0, 0, 0, 0
  12.             }
  13.         }
  14.     }
  15.     }

数据文件包含诸如你存储的数组的大小等信息，不过你可能还想添加相关的标量信息。例如，如果数组是一个程序的输出，你可以记录用什么输入参数生成的。

  1.         parmspace = H5Screate(H5S_SCALAR);
  2.         parm_id = H5Dcreate(file_id,"/parm",H5T_NATIVE_INT,parmspace,H5P_DEFAULT);

练习：通过编译和运行下面的 parmwrite.c 代码，在 HDF5 文件中添加一个标量数据空间。

预期结果：一个新的 wdset.h5 文件被创建。

  1.         /*
  2.         * File: parmdataset.c
  3.         * Author: Victor Eijkhout
  4.         */
  5.         #include "myh5defs.h"
  6.         #define FILE "pdset.h5"
  8.         main() {
  10.             hid_t file_id, dataset_id, dataspace_id; /* identifiers */
  11.             hid_t parm_id,parmspace;
  12.             hsize_t dims[2];
  13.             herr_t status;
  15.             /* Create a new file using default properties. */
  16.             file_id = H5Fcreate(FILE, H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT);
  18.             /* Create the data space for the dataset. */
  19.             dims[0] = 4;
  20.             dims[1] = 6;
  21.             dataspace_id = H5Screate_simple(2, dims, NULL);
  23.             /* Create the dataset. */
  24.             dataset_id = H5Dcreate(file_id, "/dset", H5T_STD_I32BE, dataspace_id, H5P_DEFAULT);
  26.             /* Add a descriptive parameter */
  27.             parmspace = H5Screate(H5S_SCALAR);
  28.             parm_id = H5Dcreate(file_id,"/parm",H5T_NATIVE_INT,parmspace,H5P_DEFAULT);
  30.             /* End access to the dataset and release resources used by it. */
  31.             status = H5Dclose(dataset_id);
  32.             status = H5Dclose(parm_id);
  34.             /* Terminate access to the data space. */
  35.             status = H5Sclose(dataspace_id);
  36.             status = H5Sclose(parmspace);
  38.             /* Close the file. */
  39.             status = H5Fclose(file_id);
  40.         }
  42.         %% h5dump wdset.h5
  43.         HDF5 "wdset.h5" {
  44.         GROUP "/" {
  45.             DATASET "dset" {
  46.                 DATATYPE H5T_IEEE_F64LE
  47.                 DATASPACE SIMPLE { ( 4, 6 ) / ( 4, 6 ) }
  48.                 DATA {
  49.                 (0,0): 0.5, 1.5, 2.5, 3.5, 4.5, 5.5,
  50.                 (1,0): 6.5, 7.5, 8.5, 9.5, 10.5, 11.5,
  51.                 (2,0): 12.5, 13.5, 14.5, 15.5, 16.5, 17.5,
  52.                 (3,0): 18.5, 19.5, 20.5, 21.5, 22.5, 23.5
  53.                 }
  54.             }
  55.         DATASET "parm" {
  56.             DATATYPE H5T_STD_I32LE
  57.             DATASPACE SCALAR
  58.             DATA {
  59.             (0): 37
  60.             }
  61.         }
  62.     }
  63.     }

写入数据

你创建的数据集分配了 HDF5 文件中的空间，现在你需要把实际的数据放在里面。这就是用 H5Dwrite 调用完成

  1.         /* Write floating point data */
  2.         for (i=0; i<24; i++) data[i] = i+.5;
  3.         status = H5Dwrite(dataset,H5T_NATIVE_DOUBLE,H5S_ALL,H5S_ALL,H5P_DEFAULT,data);
  4.         /* write parameter value */
  5.         parm = 37;
  6.         status = H5Dwrite(parmset,H5T_NATIVE_INT,H5S_ALL,H5S_ALL,H5P_DEFAULT,&parm);
  8.         /*
  9.         * File: parmwrite.c
  10.         * Author: Victor Eijkhout
  11.         */
  12.         #include "myh5defs.h"
  13.         #define FILE "wdset.h5"
  15.         main() {
  17.         hid_t file_id, dataset, dataspace; /* identifiers */
  18.         hid_t parmset,parmspace;
  19.         hsize_t dims[2];
  20.         herr_t status;
  21.         double data[24]; int i,parm;
  23.         /* Create a new file using default properties. */
  24.         file_id = H5Fcreate(FILE, H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT);
  26.         /* Create the dataset. */
  27.         dims[0] = 4; dims[1] = 6;
  28.         dataspace = H5Screate_simple(2, dims, NULL);
  29.         dataset = H5Dcreate(file_id, "/dset", H5T_NATIVE_DOUBLE, dataspace, H5P_DEFAULT);
  31.         /* Add a descriptive parameter */
  32.         parmspace = H5Screate(H5S_SCALAR);
  33.         parmset = H5Dcreate(file_id,"/parm",H5T_NATIVE_INT,parmspace,H5P_DEFAULT);
  35.         /* Write data to file */
  36.         for (i=0; i<24; i++) data[i] = i+.5;
  37.         status = H5Dwrite(dataset,H5T_NATIVE_DOUBLE,H5S_ALL,H5S_ALL,H5P_DEFAULT,data); H5REPORT(status);
  39.         /* write parameter value */
  40.         parm = 37;
  41.         status = H5Dwrite(parmset,H5T_NATIVE_INT,H5S_ALL,H5S_ALL,H5P_DEFAULT,&parm); H5REPORT(status);
  42.         /* End access to the dataset and release resources used by it. */
  43.         status = H5Dclose(dataset);
  44.         status = H5Dclose(parmset);
  46.         /* Terminate access to the data space. */
  47.         status = H5Sclose(dataspace);
  48.         status = H5Sclose(parmspace);
  50.         /* Close the file. */
  51.         status = H5Fclose(file_id);
  52. }
  54.         %% h5dump wdset.h5
  55.         HDF5 "wdset.h5" {
  56.         GROUP "/" {
  57.             DATASET "dset" {
  58.                 DATATYPE H5T_IEEE_F64LE
  59.                 DATASPACE SIMPLE { ( 4, 6 ) / ( 4, 6 ) }
  60.                 DATA {
  61.                 (0,0): 0.5, 1.5, 2.5, 3.5, 4.5, 5.5,
  62.                 (1,0): 6.5, 7.5, 8.5, 9.5, 10.5, 11.5,
  63.                 (2,0): 12.5, 13.5, 14.5, 15.5, 16.5, 17.5,
  64.                 (3,0): 18.5, 19.5, 20.5, 21.5, 22.5, 23.5
  65.                 }
  66.             }
  67.             DATASET "parm" {
  68.                 DATATYPE H5T_STD_I32LE
  69.                 DATASPACE SCALAR
  70.                 DATA {
  71.                 (0): 37
  72.                 }
  73.             }
  74. }
  75. }

如果你仔细看一下源代码和备份文件系统，你会发现数据类型被声明为 “本地”，但是呈现为LE 。“本地”声明使数据类型的行为与内置的 C 或 Fortran 数据类型相似，另外，你也可以明确指出数据是小端字节序还是大端字节序。这些术语描述了一个数据项的字节在内存中的排序方式。正如你在备份文件系统输出中所看到的那样，大多数架构都使用小端字节序（little endian）。但值得注意的是，IBM 使用的是大端字节序。

读取

现在我们有了一个包含一些数据的文件，我们可以从该文件中读取数据。这些基本命令是

  1.         h5file = H5Fopen( .... )
  2.         ....
  3.         H5Dread( dataset, .... data .... )

其中 H5Dread 命令的参数与相应的 H5Dwrite 相同。

练习：通过编译和运行下面的 allread.c 例子，从你在前面的练习中创建的 wdset.h5 文件中读取数据。

预期结果：运行可执行文件 allread 将输出参数的值 37 和数组的(1,2)数据点的值8.5。

注意事项：确保你运行 parmwrite 来创建输入文件。

  1.         /*
  2.         * File: allread.c
  3.         * Author: Victor Eijkhout
  4.         */
  5.         #include "myh5defs.h"
  6.         #define FILE "wdset.h5"
  8.         main() {
  10.             hid_t file_id, dataset, parmset;
  11.             herr_t status;            
  12.             double data[24]; int parm;
  14.             /* Open an existing file */
  15.             file_id = H5Fopen(FILE, H5F_ACC_RDONLY, H5P_DEFAULT);
  16.             H5REPORT(file_id);
  18.             /* Locate the datasets. */
  19.             dataset = H5Dopen(file_id, "/dset"); H5REPORT(dataset);
  20.             parmset = H5Dopen(file_id,"/parm"); H5REPORT(parmset);
  22.             /* Read data back */
  23.             status = H5Dread(parmset,H5T_NATIVE_INT,H5S_ALL,H5S_ALL,H5P_DEFAULT,&parm); H5REPORT(status);
  24.             printf("parameter value: %d\n",parm);
  26.             status = H5Dread(dataset,H5T_NATIVE_DOUBLE,H5S_ALL,H5S_ALL,H5P_DEFAULT,data); H5REPORT(status);
  27.             printf("arbitrary data point [1,2]: %e\n",data[1*6+2]);
  29.             /* Terminate access to the datasets */
  30.             status = H5Dclose(dataset); H5REPORT(status);
  31.             status = H5Dclose(parmset); H5REPORT(status);
  32.             /* Close the file. */
  33.             status = H5Fclose(file_id);
  34. }
  35.         %% ./allread
  36.         parameter value: 37
  37.         arbitrary data point [1,2]: 8.500000e+00c