从 HDF5 文件读取和写入 numpy 数组

Reading and writing numpy arrays to and from HDF5 files

我正在构建模拟软件,我需要将(数千个)2D numpy 数组写入 HDF5 文件中的 tables,其中数组的一维是可变的。传入的array是float32类型;为了节省磁盘 space,每个数组都存储为 table,并为列提供适当的数据类型(因此不使用数组)。当我读取 tables 时,我想检索一个 float32 类型的 numpy.ndarray,这样我就可以进行很好的计算以供分析。下面是一个包含 A、B 和 C 物种加上时间的数组的示例代码。

我目前正在阅读和写作的方式'works'但是很慢。因此问题是:将 array 快速存储到 table 并将其再次读回 ndarrays 的合适方法是什么?我一直在试验 numpy.recarray,但我无法让它工作(类型错误、尺寸错误、完全错误的数字等)?

代码:

import tables as pt
import numpy as np

# Variable dimension
var_dim=100

# Example array, rows 0 and 3 should be stored as float32, rows 1 and 2 as uint16
array=(np.random.random((4, var_dim)) * 100).astype(dtype=np.float32)

filename='test.hdf5'
hdf=pt.open_file(filename=filename,mode='w')
group=hdf.create_group(hdf.root,"group")

particle={
    'A':pt.Float32Col(),
    'B':pt.UInt16Col(),
    'C':pt.UInt16Col(),
    'time':pt.Float32Col(),
    }
dtypes=np.array([
    np.float32,
    np.uint16,
    np.uint16,
    np.float32
    ])

# This is the table to be stored in
table=hdf.create_table(group,'trajectory', description=particle, expectedrows=var_dim)

# My current way of storing
for i, row in enumerate(array.T):
    table.append([tuple([t(x) for t, x in zip(dtypes, row)])])
table.flush()
hdf.close()


hdf=pt.open_file(filename=filename,mode='r')
array_table=hdf.root.group._f_iter_nodes().__next__()

# My current way of reading
row_list = []
for i, row in enumerate(array_table.read()):
    row_list.append(np.array(list(row)))

#The retreived array
array=np.asarray(row_list).T


# I've tried something with a recarray
rec_array=array_table.read().view(type=np.recarray)

# This gives me errors, or wrong results
rec_array.view(dtype=np.float64)
hdf.close()

我得到的错误:

Traceback (most recent call last):
  File "/home/thomas/anaconda3/lib/python3.6/site-packages/numpy/core/records.py", line 475, in __setattr__
    ret = object.__setattr__(self, attr, val)
ValueError: new type not compatible with array.

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "/home/thomas/Documents/Thesis/SO.py", line 53, in <module>
    rec_array.view(dtype=np.float64)
  File "/home/thomas/anaconda3/lib/python3.6/site-packages/numpy/core/records.py", line 480, in __setattr__
    raise exctype(value)
ValueError: new type not compatible with array.
Closing remaining open files:test.hdf5...done

作为一种快速而肮脏的解决方案,可以通过临时将数组转换为列表来避免循环(如果您可以节省内存)。由于某些原因,记录数组很容易转换为 to/from 列表而不是 to/from 常规数组。

存储:

table.append(array.T.tolist())

正在加载:

loaded_array = np.array(array_table.read().tolist(), dtype=np.float64).T

应该有更多 "Numpythonic" 方法在记录数组和常规数组之间进行转换,但我对前者不够熟悉,不知道如何操作。

我没有使用过 tables,但使用 h5py 查看过它的文件。我猜你的 arrayrecarray 是一个结构化数组,其 dtype 如下:

In [131]: dt=np.dtype('f4,u2,u2,f4')
In [132]: np.array(arr.tolist(), float)
Out[132]: 
array([[ 1.,  1.,  1.,  1.],
       [ 1.,  1.,  1.,  1.],
       [ 1.,  1.,  1.,  1.]])
In [133]: arr
Out[133]: 
array([( 1., 1, 1,  1.), ( 1., 1, 1,  1.), ( 1., 1, 1,  1.)], 
      dtype=[('f0', '<f4'), ('f1', '<u2'), ('f2', '<u2'), ('f3', '<f4')])

使用 @kazemakase's tolist 方法(我在其他帖子中推荐过):

In [134]: np.array(arr.tolist(), float)
Out[134]: 
array([[ 1.,  1.,  1.,  1.],
       [ 1.,  1.,  1.,  1.],
       [ 1.,  1.,  1.,  1.]])

astype 弄错了形状

In [135]: arr.astype(np.float32)
Out[135]: array([ 1.,  1.,  1.], dtype=float32)

view 在组件数据类型统一时工作,例如使用 2 个浮点字段

In [136]: arr[['f0','f3']].copy().view(np.float32)
Out[136]: array([ 1.,  1.,  1.,  1.,  1.,  1.], dtype=float32)

但它确实需要重塑。 view 使用数据缓冲区字节,只是重新解释。

许多 recfunctions 函数使用逐字段复制。这里相当于

In [138]: res = np.empty((3,4),'float32')
In [139]: for i in range(4):
     ...:     res[:,i] = arr[arr.dtype.names[i]]
     ...:     
In [140]: res
Out[140]: 
array([[ 1.,  1.,  1.,  1.],
       [ 1.,  1.,  1.,  1.],
       [ 1.,  1.,  1.,  1.]], dtype=float32)

如果字段数与记录数相比较小,则此迭代并不昂贵。


def foo(arr):
    res = np.empty((arr.shape[0],4), np.float32)
    for i in range(4):
        res[:,i] = arr[arr.dtype.names[i]]
    return res

对于大型 4 字段数组,按字段复制显然更快:

In [143]: arr = np.ones(10000, dtype=dt)
In [149]: timeit x1 = foo(arr)
10000 loops, best of 3: 73.5 µs per loop
In [150]: timeit x2 = np.array(arr.tolist(), np.float32)
100 loops, best of 3: 11.9 ms per loop