使用:Python Numpy Ssliceing从2D索引列表(x，y，z)递增计数器的3D数组

Question

我想从 2D 事件数组 (x，y，t( 中递增计数器的 3D 矩阵 (nparray( 下面的代码工作：

TOF_cube=np.zeros((324,324,4095),np.int32) #initialise a 3d array for whole data set
data = np.fromfile(f, dtype='<i2', count=no_I16) #read all events, x,y,t as 1D array
data=data.reshape(events,cols)
xpos=data[:,0]
ypos=data[:,1]
tpos=data[:,2]
i=0
while i < events:              
TOF_cube[xpos[i],ypos[i],tpos[i]] += 1
i+=1

为了使用切片和索引，我将我的 while 循环替换为

TOF_cube[xpos,ypos,tpos] += 1

但是，我没有复制正确的4365520数量的事件(通过while循环并独立检查(，我只记录4365197。

为什么切片方法会丢失事件？

我在 while 循环中使用完全相同的切片并作为索引的"参数"。

Answer 1

如果有重复索引，+=不会加两次。

要以矢量化方式获得等效输出，您需要np.add.at：

np.add.at(TOF_cube, [xpos, ypos, tpos], 1)

Answer 2

由于我们不知道您的数据到底是什么样子的，因此很难猜测实际问题是什么。如果这没有帮助，请举一个例子，我们可以运行我们的elvs(即没有文件f(。

假设你有x_pos = [1,1,2,3,5]

a = np.zeros(10)
for i in range(len(x_pos)):
a[x_pos[i]]+=1
# gives a = array([ 0.,  2.,  1.,  1.,  0.,  1.,  0.,  0.,  0.,  0.])

但是其他代码

a[x_pos]+=1
# gives a = array([ 0.,  1.,  1.,  1.,  0.,  1.,  0.,  0.,  0.,  0.])

因此，如果其中一个索引出现两次，则在短版本中仅更新一次。检查您的 xpos 等是否确实如此。

PS：我做了一个稍微简单的版本，只有一个维度，但规则保持不变。