添加 numpy 数组时避免溢出答案

【问题标题】：Avoid overflow when adding numpy arrays添加 numpy 数组时避免溢出
【发布时间】：2015-06-19 02:37:31
【问题描述】：

我想用 datatyp uint8 添加 numpy 数组。我知道这些数组中的值可能大到足以发生溢出。所以我得到类似的东西：

a = np.array([100, 200, 250], dtype=np.uint8)
b = np.array([50, 50, 50], dtype=np.uint8)
a += b

现在，a 是[150 250 44]。但是，我希望对于 uint8 来说太大的值而不是溢出，而不是 uint8 允许的最大值。所以我想要的结果是[150 250 255]。

我可以用下面的代码得到这个结果：

a = np.array([100, 200, 250], dtype=np.uint8)
b = np.array([50, 50, 50], dtype=np.uint8)
c = np.zeros((1,3), dtype=np.uint16)
c += a
c += b
c[c>255] = 255
a = np.array(c, dtype=np.uint8)

问题是，我的数组非常大，因此创建具有更大数据类型的第三个数组可能是内存问题。是否有一种快速且内存效率更高的方法来实现所描述的结果？

【问题讨论】：

标签： python numpy image-processing integer-overflow numpy-ndarray

【解决方案1】：

怎么样

>>> a + np.minimum(255 - a, b)
array([150, 250, 255], dtype=uint8)

通常使用

获取您的数据类型的最大值

np.iinfo(np.uint8).max

【讨论】：

@PadraicCunningham，确实如此，但 dtype 为 uint8，而不是 uint16。但是，它会创建三个 temp uint8 数组..
@shx2：我只数了两个。 255 - a 和 np.minimum(255 - a, b)。第三个是什么？
@user2357112，a + ...。如果 OP 想要结果代替 a 数组，则可以避免这种情况。
我发现这个答案对于我的图片 (uint8) 添加来说是最快的

【解决方案2】：

您可以通过创建第三个 dtype uint8 数组和一个 bool 数组~~（它们一起比一个 uint16 数组更节省内存）~~。

np.putmask 对于避免临时数组很有用。

a = np.array([100, 200, 250], dtype=np.uint8)
b = np.array([50, 50, 50], dtype=np.uint8)
c = 255 - b  # a temp uint8 array here
np.putmask(a, c < a, c)  # a temp bool array here
a += b

但是，正如@moarningsun 正确指出的那样，bool 数组与 uint8 数组占用的内存量相同，因此这不一定有用。可以通过避免在任何给定时间拥有多个临时数组来解决此问题：

a = np.array([100, 200, 250], dtype=np.uint8)
b = np.array([50, 50, 50], dtype=np.uint8)
b = 255 - b  # old b is gone shortly after new array is created
np.putmask(a, b < a, b)  # a temp bool array here, then it's gone
a += 255 - b  # a temp array here, then it's gone

这种方法用内存消耗换取 CPU。

另一种方法是预先计算所有可能的结果，即 O(1) 额外内存（即与数组大小无关）：

c = np.clip(np.arange(256) + np.arange(256)[..., np.newaxis], 0, 255).astype(np.uint8)
c
=> array([[  0,   1,   2, ..., 253, 254, 255],
          [  1,   2,   3, ..., 254, 255, 255],
          [  2,   3,   4, ..., 255, 255, 255],
          ..., 
          [253, 254, 255, ..., 255, 255, 255],
          [254, 255, 255, ..., 255, 255, 255],
          [255, 255, 255, ..., 255, 255, 255]], dtype=uint8)

c[a,b]
=> array([150, 250, 255], dtype=uint8)

如果您的数组非常大，这种方法是最节省内存的。同样，它的处理时间很昂贵，因为它用较慢的 2dim 数组索引取代了超快速的整数加法。

解释它是如何工作的

上面c 数组的构造使用了一个numpy 广播技巧。添加形状为(N,) 的数组和形状为(1,N) 的数组都广播为(N,N)-like，因此结果是所有可能和的NxN 数组。然后，我们剪辑它。我们得到一个满足：c[i,j]=min(i+j,255) 对于每个 i,j 的 2dim 数组。

剩下的就是使用花哨的索引来抓取正确的值。使用您提供的输入，我们访问：

c[( [100, 200, 250] , [50, 50, 50] )]

第一个索引数组是指第一个昏暗，第二个是第二个昏暗。因此，结果是一个与索引数组 ((N,)) 形状相同的数组，由值 [ c[100,50] , c[200,50] , c[250,50] ] 组成。

【讨论】：

不知道putmask，谢谢！使用该功能，我认为a += b 后跟np.putmask(a, a<b, 255) 也应该可以工作？
@moarningsun 我认为你是对的。然而它依赖于溢出，我个人觉得不太舒服......
@moarningsun 你为什么删除你的答案？我认为这是一个不错的答案，并且有效
该方法需要putmask 才能有效，这是您回答的要点之一。所以我想我不妨把它作为评论放在这里。
我试过你的第三种方法，效果很好。但是，您能否再解释一下，它的作用是什么？我将这第三种方法（“预计算”）的执行时间与我只使用 uint16 数组作为 a、a+=b 和 a[a>255]=255 的情况进行了比较。 “预先计算”大约需要两倍的时间。但如果内存效率更高，我将能够使用 uint16 方法对内存不足的大型数组进行计算。

【解决方案3】：

这是一种方法：

>>> a = np.array([100, 200, 250], dtype=np.uint8)
>>> b = np.array([50, 50, 50], dtype=np.uint8)
>>> a+=b; a[a<b]=255
>>> a
array([150, 250, 255], dtype=uint8)

【讨论】：

这是可行的，因为溢出就像将数字设为负数一样。 200 + 100 => 44 因为 200 距离 256 56(=0) 所以 -56+100 = 44。所以如果发生溢出，结果总是小于 b。

【解决方案4】：

您可以使用 Numba 真正就地完成，例如：

import numba

@numba.jit('void(u1[:],u1[:])', locals={'temp': numba.uint16})
def add_uint8_inplace_clip(a, b):
    for i in range(a.shape[0]):
        temp = a[i] + b[i]
        a[i] = temp if temp<256 else 255

add_uint8_inplace_clip(a, b)

或者用 Numexpr，例如：

import numexpr

numexpr.evaluate('where((a+b)>255, 255, a+b)', out=a, casting='unsafe')

Numexpr upcasts uint8 到 int32 内部，然后将其放回 uint8 数组中。

【讨论】：

【解决方案5】：

这里有a function in numpy：

numpy.nan_to_num(x)[source]

将 nan 替换为零，将 inf 替换为有限数。

返回一个数组或标量，用零替换非数字 (NaN)，用非常大的数字替换（正）无穷大，用非常小的（或负）数字替换负无穷大。

与 x 形状和 x 中精度最高的元素的 dtype 相同的新数组。

如果 x 不精确，则 NaN 被零替换，无穷大 (-infinity) 被适合输出 dtype 的最大（最小或最大负数）浮点值替换。如果 x 不是不精确的，则返回 x 的副本。

我不确定它是否适用于 uint8，因为在输出中提到了浮点，但对于其他读者来说，它可能有用

【讨论】：

我不明白这对这个问题有何帮助。在任何要添加的数组中都没有 NaN 或无限值。所以也许我错过了你的答案？
@Thomas hmm，也许整数类型不同，但是当我遇到浮点数问题时，溢出显示为 +/- 无穷大
@ToshinouKyouko 是的，整数确实不一样，它们只是像OP的例子一样溢出。

【解决方案6】：

def non_overflowing_sum(a, b)
    c = np.uint16(a)+b
    c[np.where(c>255)] = 255
    return np.uint8( c )

它也交换内存，但我发现更优雅，并且在返回转换后临时 uint16 被释放

【讨论】：

【解决方案7】：

OpenCV有这样一个函数：cv2.addWeighted

【讨论】：