Python：读取 12 位二进制文件答案

【问题标题】：Python: reading 12-bit binary filesPython：读取 12 位二进制文件
【发布时间】：2017-11-27 21:06:03
【问题描述】：

我正在尝试使用 Python 3 读取包含图像（视频）的 12 位二进制文件。

要读取类似的文件，但以 16 位编码，以下方法非常有效：

import numpy as np
images = np.memmap(filename_video, dtype=np.uint16, mode='r', shape=(nb_frames, height, width))

其中 filename_video 是文件，nb_frames 是可以从另一个文件中读取的视频的高度和宽度特征。 “运行良好”是指快速：在我的计算机上读取 140 帧的 640x256 视频大约需要 1 毫秒。

据我所知，当文件以 12 位编码时，我不能使用它，因为没有 uint12 类型。所以我要做的是读取一个 12 位文件并将其存储在一个 16 位 uint 数组中。以下摘自 (Python: reading 12 bit packed binary image) 的作品：

with open(filename_video, 'rb') as f:
    data=f.read()
images=np.zeros(int(2*len(data)/3),dtype=np.uint16)
ii=0
for jj in range(0,int(len(data))-2,3):
    a=bitstring.Bits(bytes=data[jj:jj+3],length=24)
    images[ii],images[ii+1] = a.unpack('uint:12,uint:12')
    ii=ii+2
images = np.reshape(images,(nb_frames,height,width))

但是，这非常慢：使用我的机器读取只有 5 帧的 640x256 视频大约需要 11.5 秒。理想情况下，我希望能够像使用 memmap 读取 8 或 16 位文件一样有效地读取 12 位文件。或者至少不会慢 10^5 倍。我怎样才能加快速度？

这是一个文件示例： http://s000.tinyupload.com/index.php?file_id=26973488795334213426 （nb_frames=5，高度=256，宽度=640）。

【问题讨论】：

这很酷（+1）。我只知道像 OpenCV 的 cv.cvtColor(bayer,rgb,cv.COLOR_BayerBG2BGR)

标签： python python-3.x binary unpack

【解决方案1】：

我的实现与@max9111 提出的实现略有不同，不需要调用unpackbits。

它通过将中间字节切成两半并使用 numpy 的二进制操作直接从三个连续的 uint8 创建两个 uint12 值。在下文中，data_chunks 被假定为一个二进制字符串，其中包含任意数量的 12 位整数的信息（因此其长度必须是 3 的倍数）。

def read_uint12(data_chunk):
    data = np.frombuffer(data_chunk, dtype=np.uint8)
    fst_uint8, mid_uint8, lst_uint8 = np.reshape(data, (data.shape[0] // 3, 3)).astype(np.uint16).T
    fst_uint12 = (fst_uint8 << 4) + (mid_uint8 >> 4)
    snd_uint12 = ((mid_uint8 % 16) << 8) + lst_uint8
    return np.reshape(np.concatenate((fst_uint12[:, None], snd_uint12[:, None]), axis=1), 2 * fst_uint12.shape[0])

我对其他实现进行了基准测试，事实证明这种方法在 ~5 Mb 输入上的速度提高了 ~4 倍：
read_uint12_unpackbits 65.5 ms ± 1.11 ms 每个循环（平均值 ± 标准偏差，7 次运行，10每个循环） read_uint12 每个循环 14 ms ± 513 µs（平均值 ± 标准偏差，7 次运行，每次 100 个循环）

【讨论】：

【解决方案2】：

加快 numpy 向量化方法的一种方法是避免为临时数据分配昂贵的内存，更有效地使用缓存并利用并行化。这可以使用Numba、Cython 或C 轻松完成。请注意，并行化并不总是有益的。如果要转换的数组太小，请使用单线程版本（parallel=False）

Numba 版本的 Cyril Gaudefroy 使用临时内存分配回答

import numba as nb
import numpy as np
@nb.njit(nb.uint16[::1](nb.uint8[::1]),fastmath=True,parallel=True)
def nb_read_uint12(data_chunk):
  """data_chunk is a contigous 1D array of uint8 data)
  eg.data_chunk = np.frombuffer(data_chunk, dtype=np.uint8)"""
  
  #ensure that the data_chunk has the right length
  assert np.mod(data_chunk.shape[0],3)==0
  
  out=np.empty(data_chunk.shape[0]//3*2,dtype=np.uint16)
  
  for i in nb.prange(data_chunk.shape[0]//3):
    fst_uint8=np.uint16(data_chunk[i*3])
    mid_uint8=np.uint16(data_chunk[i*3+1])
    lst_uint8=np.uint16(data_chunk[i*3+2])
    
    out[i*2] =   (fst_uint8 << 4) + (mid_uint8 >> 4)
    out[i*2+1] = ((mid_uint8 % 16) << 8) + lst_uint8
    
  return out

Cyril Gaudefroy 的 Numba 版本通过内存预分配回答

如果你在类似大小的数据块上多次应用这个函数，你只能预分配输出数组一次。

@nb.njit(nb.uint16[::1](nb.uint8[::1],nb.uint16[::1]),fastmath=True,parallel=True,cache=True)
def nb_read_uint12_prealloc(data_chunk,out):
    """data_chunk is a contigous 1D array of uint8 data)
    eg.data_chunk = np.frombuffer(data_chunk, dtype=np.uint8)"""

    #ensure that the data_chunk has the right length
    assert np.mod(data_chunk.shape[0],3)==0
    assert out.shape[0]==data_chunk.shape[0]//3*2

    for i in nb.prange(data_chunk.shape[0]//3):
        fst_uint8=np.uint16(data_chunk[i*3])
        mid_uint8=np.uint16(data_chunk[i*3+1])
        lst_uint8=np.uint16(data_chunk[i*3+2])

        out[i*2] =   (fst_uint8 << 4) + (mid_uint8 >> 4)
        out[i*2+1] = ((mid_uint8 % 16) << 8) + lst_uint8

    return out

带有临时内存分配的 DGrifffith 答案的 Numba 版本

@nb.njit(nb.uint16[::1](nb.uint8[::1]),fastmath=True,parallel=True,cache=True)
def read_uint12_var_2(data_chunk):
    """data_chunk is a contigous 1D array of uint8 data)
    eg.data_chunk = np.frombuffer(data_chunk, dtype=np.uint8)"""

    #ensure that the data_chunk has the right length
    assert np.mod(data_chunk.shape[0],3)==0

    out=np.empty(data_chunk.shape[0]//3*2,dtype=np.uint16)

    for i in nb.prange(data_chunk.shape[0]//3):
        fst_uint8=np.uint16(data_chunk[i*3])
        mid_uint8=np.uint16(data_chunk[i*3+1])
        lst_uint8=np.uint16(data_chunk[i*3+2])

        out[i*2] =   (fst_uint8 << 4) + (mid_uint8 >> 4)
        out[i*2+1] = (lst_uint8 << 4) + (15 & mid_uint8)

    return out

带有内存预分配的 DGrifffith 答案的 Numba 版本

@nb.njit(nb.uint16[::1](nb.uint8[::1],nb.uint16[::1]),fastmath=True,parallel=True,cache=True)
def read_uint12_var_2_prealloc(data_chunk,out):
    """data_chunk is a contigous 1D array of uint8 data)
    eg.data_chunk = np.frombuffer(data_chunk, dtype=np.uint8)"""

    #ensure that the data_chunk has the right length
    assert np.mod(data_chunk.shape[0],3)==0
    assert out.shape[0]==data_chunk.shape[0]//3*2

    for i in nb.prange(data_chunk.shape[0]//3):
        fst_uint8=np.uint16(data_chunk[i*3])
        mid_uint8=np.uint16(data_chunk[i*3+1])
        lst_uint8=np.uint16(data_chunk[i*3+2])

        out[i*2] =   (fst_uint8 << 4) + (mid_uint8 >> 4)
        out[i*2+1] = (lst_uint8 << 4) + (15 & mid_uint8)

    return out

时间

num_Frames=10
data_chunk=np.random.randint(low=0,high=255,size=np.int(640*256*1.5*num_Frames),dtype=np.uint8)

%timeit read_uint12_gaud(data_chunk)
#11.3 ms ± 53.4 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)
#435 MB/s

%timeit nb_read_uint12(data_chunk)
#939 µs ± 24.3 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)
#5235 MB/s

out=np.empty(data_chunk.shape[0]//3*2,dtype=np.uint16)
%timeit nb_read_uint12_prealloc(data_chunk,out)
#407 µs ± 5.4 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)
#11759 MB/s

%timeit read_uint12_griff(data_chunk)
#10.2 ms ± 55.9 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)
#491 MB/s

%timeit read_uint12_var_2(data_chunk)
#928 µs ± 16.5 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)
#5297 MB/s
%timeit read_uint12_var_2_prealloc(data_chunk,out)
#403 µs ± 13.4 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)
#12227 MB/s

【讨论】：

非常感谢 numba 的提示。我在转换所需的计算时间方面遇到了类似的问题。我用 bitArray 和 bitString 包试过，但循环太慢了。以您最近建议的方式实现它！再次感谢，很好的答案。
嘿@max9111，这是一个很好的答案，您能否建议修改此 numba 操作以获得以下答案（它使用不同的打包/位顺序）。这将非常有用！
@Gittb 你指的是哪一个（用户名）？答案的顺序可能不同。（我认为我的答案是最后一个）是的，但不是今天。周末。
用户名@DGrifffith 的解决方案不着急。我感谢您的帮助。如果您在工作时可以 ping 我，那就太棒了。
谢谢 Max，我真的很感激！

【解决方案3】：

发现@cyrilgaudefroy 的答案很有用。但是，最初，它不适用于我的 12 位打包二进制图像数据。发现在这种特殊情况下包装有点不同。 “中间”字节包含最不重要的半字节。三元组的字节 1 和 3 是 12 位中最高的 8 位。因此修改了@cyrilgaudefroy 的答案：

def read_uint12(data_chunk):
    data = np.frombuffer(data_chunk, dtype=np.uint8)
    fst_uint8, mid_uint8, lst_uint8 = np.reshape(data, (data.shape[0] // 3, 3)).astype(np.uint16).T
    fst_uint12 = (fst_uint8 << 4) + (mid_uint8 >> 4)
    snd_uint12 = (lst_uint8 << 4) + (np.bitwise_and(15, mid_uint8))
    return np.reshape(np.concatenate((fst_uint12[:, None], snd_uint12[:, None]), axis=1), 2 * fst_uint12.shape[0])

【讨论】：

【解决方案4】：

这是另一个变体。我的数据格式是：

第一个 uint12：第二个 uint8 的最低有效 4 位中的最高有效 4 位 + 第一个 uint8 的最低有效 8 位

第二个 uint12：第三个 uint8 的最高有效 8 位 + 第二个 uint8 最高有效 4 位的最低有效 4 位

对应的代码是：

def read_uint12(data_chunk):
    data = np.frombuffer(data_chunk, dtype=np.uint8)
    fst_uint8, mid_uint8, lst_uint8 = np.reshape(data, (data.shape[0] // 3, 3)).astype(np.uint16).T
    fst_uint12 = ((mid_uint8 & 0x0F) << 8) | fst_uint8
    snd_uint12 = (lst_uint8 << 4) | ((mid_uint8 & 0xF0) >> 4)
    return np.reshape(np.concatenate((fst_uint12[:, None], snd_uint12[:, None]), axis=1), 2 * fst_uint12.shape[0])

【讨论】：

这个变体也对我有用。它与 Matlab 的 ubit12 表示相匹配。