如何提高 Python 数组操作性能.

[jeeyong]

最近正在写一个小系统. 涉及到把 CT 影像转换成 PNG.
CT 影像不是不是标准通用格式, 所以需要从读取 bytes 再转成数组..
这个过程中需要两次处理数组数据, 但是性能很慢..有相关经验的朋友可以给点建议嘛?
以下是详述:

我读出来的 bytes 数据处理成 uint8 后. 是这样的形式:
[208, 4, 208, 4, 208, 4...196, 8]
不懂图像处理的知识, 我的理解就是, 一个灰度, 一个 Alpha 通道值(透明值?).
第一次处理数组是要把上面的数组, 改成如下形式: <- 暂且叫 生成数组阶段.
[208, 208, 208, 4, 208, 208, 208, 4.......]
就是把 208 这个灰度值变成 rgb 的形式.
然后再通过一个循环变成 pillow 支持的格式, 如下: <- 暂且叫 数组转换阶段吧.
[
[ [r, g, b, a], [r, g, b, a], [r, g, b, a], ],
[ [r, g, b, a], [r, g, b, a], [r, g, b, a], ],
[ [r, g, b, a], [r, g, b, a], [r, g, b, a], ],
]

然后这个数组要通过 numpy 转换成 uint8 类型才能提交给 pillow 或者类似图片处理的库..

所以大致耗时分为三个阶段:
1. 生成数组
2. 转换数组
3. 通过 numpy 转换数组类型为 uint8.

数组大小为: 5022 * 4200 * 4 = 84,369,600

环境 A 介绍:
Win10, 10900K, 64GB 3600, NVME 2TB SSD
Python 3.10.4, pillow 最新版
阶段一耗时: 9 秒左右
阶段二耗时: 45 秒左右
阶段三耗时: 5 秒左右

觉得太慢, 尝试全部改成 numpy 操作.
但是阶段二慢的出奇, 后来查资料发现, numpy 在大量数组下标赋值操作的性能还不如 Python 原生 list.
阶段二, 上厕所, 洗把脸回来还没跑完, 就直接放弃了. 也尝试过 np.insert 操作, 一样慢.

后来尝试环境 B, 基于 pypy 的.
pypy3.8, numpy, pillow, 同样的硬件.
第一阶段耗时: 0.2 秒
第二阶段耗时: 9.7 秒
第三阶段耗时: 77 秒
查阅资料得知, pypy 的 C 扩展接口性能很差, 不如原生 python.
而 Numpy 就是 C 扩展的库..所以导致 Numpy 性能急剧下降.

好了..剩下的办法超出我的知识范畴了..
有朋友能分享一下经验或者给个可能的方向我去看也行.

需要转换的数据量大约有 23 万个 CT 文件. 按照一个文件 1 分钟, 我即便 4 个进程跑(文件转换的服务器 4 核 8GB 内存, 还要跑一些其他服务.), 5 万分钟 / 60 分钟 / 24 小时 一点意外没有的情况要跑 42 天...

顺便吐槽以下, 朋友给写了个 nodejs 的版本, 转换一个图 2.3 秒, WTF!!

[LnTrx]

所有运算尽量用矩阵实现，然后直接用 pytorch 就行了

[shinsekai]

同样结构的代码，matlab 比 numpy 快 20 倍，楼主可以试试

[Muniesa]

我应该没理解错吧

[inframe]

简单，就是向量化思想方法
resharp 2 列 N 行，[v0,v1],v 是个列向量
然后拼装列向量为[v0,v0,v0,v1]

raw=np.random.randint(0,10, 10**8).reshape((-1,2))
v1=np.vstack([raw[:,0],raw[:,0],raw[:,1]])
output=v1.T

[mizuBai]

楼主可以试试写 cython 扩展，或者可以试试 f2py ，不过感觉 Fortran 这么上古的语言楼主会的概率不大（

[Purelove]

rust pyo3 了解一下

[jeeyong]

@Cy86
@jdhao
@chashao
@paopjian
老哥们, 辛苦啊...
我现在是直接用真实的患者 CT 片子写..
不能外传....

另外, 我看有人已经上 pydicom 了...
另外我请教一下, 我现在处理的 dicom 格式不是通用的 dicom 格式...
是排版好之后的 dcm 文件.
很多信息不能直接获取..最重要的 pydicom.pixel_array 和 pixelData 是没有的..
我是通过 pydicom
dcm = pydicom.read_file(filePath)
imageInfo = dcm.ReferencedImageBoxSequence[0].BasicGrayscaleImageSequence[0]
这种方式获取的 pixelData 的 bytes
有人能科普一下这个问题吗?

[jeeyong]

@Muniesa 老哥, 图完全看不清楚

[lucays]

dicom 是可以写入 pixel array 生成新的 dcm 文件的
只是模板没有的话需要摸索下
画图直接上 numpy+opencv

[lucays]

真实的 ct 片子也可以用 pydicom 读取然后抹掉 patient 信息生成新的不记名的不就完了

你这是完全从零开始啊，没点时间搞不定的，听起来一点 dicom 的基础建设都没

[Nugine0]

cython ，rust pyo3 ，c++，nodejs ，c 语言动态库，都可以试试。
如果只需要一次性转换，还可以开好点的临时云服务器，大力出奇迹。

[Mayye]

@Nugine0 感觉现在不是机器配置的问题而是代码效率的问题了 云服务器应该帮助不大

[Xs0ul]

test_image = [208, 4, 208, 4, 208, 4, 100, 4, 100, 4, 100, 4]
test_image = np.array(test_image).reshape(2,3,2).astype(np.uint8)
image = PIL.Image.fromarray(test_image, mode='LA')
# if need RGBA
image = image.convert('RGBA')

pillow 本身就支持灰度 + alpha 的格式: https://pillow.readthedocs.io/en/stable/handbook/concepts.html

[jeeyong]

@lucays 第三天....T_T

[jeeyong]

@Xs0ul Hey~`
这方法效率很高..
但是有个问题....
他显示出来的图片...如果是 LA 模式, 很浅. 黑色的几乎看不清楚
我尝试用循环把 Alpha 通道值修改为 255 - 原值.
结果图像感觉噪点又很多
LA 模式:
https://imgur.com/QmmmcX9

转为 RGBA 模式
https://imgur.com/GPtekhA

期望得到的结果:
https://imgur.com/undefined

[Xs0ul]

@jeeyong #55 你最后一个链接是 undefined,我看不了期望的是怎么样的。但你开头给的那些数字，alpha 看起来非常的小？那显示出来就几乎是透明的。或许 208 是 alpha ，4 是灰度？具体 ct 影像的格式我不了解

[milkpuff]

上面说的向量化操作都是可以的。不要用循环索引的方式就行了。

[lambdaq]

python 要性能第一个原则就是只要你写 for 就已经输了。

[ferstar]

向量化加太极😏
https://docs.taichi-lang.org/blog/accelerate-python-code-100x

[jeeyong]

@Xs0ul 太着急了, 没粘好.
https://imgur.com/WWKXy9A