在软件包scipy中,有定义二进制结构的函数(例如出租车(2,1)或棋盘(2,2))。
import numpy
from scipy import ndimage
a = numpy.zeros((6,6), dtype=numpy.int)
a[1:5, 1:5] = 1;a[3,3] = 0 ; a[2,2] = 2
s = ndimage.generate_binary_structure(2,2) # Binary structure
#.... Calculate Sum of
result_array = numpy.zeros_like(a)
我想要的是使用给定的结构 s 迭代该数组的所有单元格。然后我想将一个函数附加到空数组中索引的当前单元格值(示例函数总和),该数组使用二进制结构中所有单元格的值。
例如:
array([[0, 0, 0, 0, 0, 0],
[0, 1, 1, 1, 1, 0],
[0, 1, 2, 1, 1, 0],
[0, 1, 1, 0, 1, 0],
[0, 1, 1, 1, 1, 0],
[0, 0, 0, 0, 0, 0]])
# 数组 a。单元格 1,2 中的值当前为 1。给定结构 s 和示例函数(如 sum),结果数组 (result_array) 中的值变为 7(如果排除当前单元格值,则变为 6)。
有人有想法吗?
对于总和的特定情况,您可以使用 ndimage.convolve:
In [42]: import numpy as np
In [43]: a = np.zeros((6,6), dtype=np.int)
a[1:5, 1:5] = 1;
a[3,3] = 0;
a[2,2] = 2
In [48]: s = ndimage.generate_binary_structure(2,2) # Binary structure
In [49]: ndimage.convolve(a,s)
Out[49]:
array([[1, 2, 3, 3, 2, 1],
[2, 5, 7, 7, 4, 2],
[3, 7, 9, 9, 5, 3],
[3, 7, 9, 9, 5, 3],
[2, 4, 5, 5, 3, 2],
[1, 2, 3, 3, 2, 1]])
对于产品的特殊情况,您可以使用log(a*b) = log(a)+log(b)
将问题转换回涉及金额的问题。例如,如果我们想"产品卷积"b
:
b = a[1:-1, 1:-1]
print(b)
# [[1 1 1 1]
# [1 2 1 1]
# [1 1 0 1]
# [1 1 1 1]]
我们可以计算:
print(np.exp(ndimage.convolve(np.log(b), s, mode = 'constant')))
# [[ 2. 2. 2. 1.]
# [ 2. 0. 0. 0.]
# [ 2. 0. 0. 0.]
# [ 1. 0. 0. 0.]]
如果b
包含负值,情况会变得更加复杂:
b[0,1] = -1
print(b)
# [[ 1 -1 1 1]
# [ 1 2 1 1]
# [ 1 1 0 1]
# [ 1 1 1 1]]
但并非不可能:
logb = np.log(b.astype('complex'))
real, imag = logb.real, logb.imag
print(np.real_if_close(
np.exp(
sum(j * ndimage.convolve(x, s, mode = 'constant')
for x,j in zip((real, imag),(1,1j))))))
# [[-2. -2. -2. 1.]
# [-2. -0. -0. 0.]
# [ 2. 0. 0. 0.]
# [ 1. 0. 0. 0.]]
如果使用 2 深的零墙会更容易:
In [11]: a0
Out[11]:
array([[ 0., 0., 0., 0., 0., 0., 0., 0.],
[ 0., 0., 0., 0., 0., 0., 0., 0.],
[ 0., 0., 1., 1., 1., 1., 0., 0.],
[ 0., 0., 1., 2., 1., 1., 0., 0.],
[ 0., 0., 1., 1., 0., 1., 0., 0.],
[ 0., 0., 1., 1., 1., 1., 0., 0.],
[ 0., 0., 0., 0., 0., 0., 0., 0.],
[ 0., 0., 0., 0., 0., 0., 0., 0.]])
In [12]: b0 = zeros_like(a0)
In [13]: for i in range(1,len(a0)-1):
....: for j in range(1,len(a0)-1):
....: b0[i,j] = sum(a0[i-1:i+2, j-1:j+2] * s)
这使您能够根据需要将两个子矩阵相乘并求和。(你也可以在这里做一些更精细的事情...
In [14]: b0
Out[14]:
array([[ 0., 0., 0., 0., 0., 0., 0., 0.],
[ 0., 1., 2., 3., 3., 2., 1., 0.],
[ 0., 2., 5., 7., 7., 4., 2., 0.],
[ 0., 3., 7., 9., 9., 5., 3., 0.],
[ 0., 3., 7., 9., 9., 5., 3., 0.],
[ 0., 2., 4., 5., 5., 3., 2., 0.],
[ 0., 1., 2., 3., 3., 2., 1., 0.],
[ 0., 0., 0., 0., 0., 0., 0., 0.]])
In [15]: b0[1:len(b0)-1, 1:len(b0)-1]
Out[15]:
array([[ 1., 2., 3., 3., 2., 1.],
[ 2., 5., 7., 7., 4., 2.],
[ 3., 7., 9., 9., 5., 3.],
[ 3., 7., 9., 9., 5., 3.],
[ 2., 4., 5., 5., 3., 2.],
[ 1., 2., 3., 3., 2., 1.]])