查找2D直方图的峰值

Question

我做一些(x, y)数据的2D直方图和我得到这样一个图像：

我想办法获得分数的(x, y)坐标（s）将最大值存储在H中。例如，在上面的图像的情况下，它将是具有aprox坐标的两点：(1090, 1040)和(1110, 1090)。

这是我的代码：

import numpy as np 
import matplotlib.pyplot as plt 
import matplotlib.cm as cm 
from os import getcwd 
from os.path import join, realpath, dirname 

# Path to dir where this code exists. 
mypath = realpath(join(getcwd(), dirname(__file__))) 
myfile = 'datafile.dat' 

x, y = np.loadtxt(join(mypath,myfile), usecols=(1, 2), unpack=True) 

fig = plt.figure() 
ax = fig.add_subplot(111) 

xmin, xmax = min(x), max(x) 
ymin, ymax = min(y), max(y) 

rang = [[xmin, xmax], [ymin, ymax]] 

binsxy = [int((xmax - xmin)/20), int((ymax - ymin)/20)] 

H, xedges, yedges = np.histogram2d(x, y, range=rang, bins=binsxy) 

extent = [yedges[0], yedges[-1], xedges[0], xedges[-1]] 
cp = ax.imshow(H.transpose()[::-1], interpolation='nearest', extent=extent, cmap=cm.jet) 
fig.colorbar(cp) 

plt.show() 

---

编辑

我已经试过张贴马立克和qarma试图获得分区中的坐标，而不是指数的解决方案他们，就像这样：

# Marek's answer 
x_cent, y_cent = unravel_index(H.argmax(), H.shape) 
print('Marek') 
print(x_cent, y_cent) 
print(xedges[x_cent], yedges[y_cent]) 

# qarma's answer 
idx = list(H.flatten()).index(H.max()) 
x_cent2, y_cent2 = idx/H.shape[1], idx % H.shape[1] 
local_maxs = np.argwhere(H == H.max()) 
print('\nqarma') 
print(x_cent2, y_cent2) 
print(xedges[x_cent2], yedges[y_cent2]) 
print(xedges[local_maxs[0,0]], yedges[local_maxs[0,1]], xedges[local_maxs[1,0]], yedges[local_maxs[1,1]]) 

哪个结果在：

Marek 
(53, 50) 
(1072.7838144329899, 1005.0837113402063) 

qarma 
(53, 50) 
(1072.7838144329899, 1005.0837113402063) 
(1072.7838144329899, 1005.0837113402063, 1092.8257731958763, 1065.3611340206187) 

所以最大的坐标是相同的，这是很好的！现在我有一个小问题，因为如果我放大的2D情节，我看到的坐标是有点偏离中心同时为全球最大和本地最大：

这是为什么？

Answer 1

这里是你如何能找到第一个全球最大

idx = list(H.flatten()).index(H.max()) 
x, y = idx/H.shape[1], idx % H.shape[1] 

找到所有最大值被留给读者做练习的协调...

numpy.argwhere(H == H.max()) 

编辑

您的编号：

H, xedges, yedges = np.histogram2d(x, y, range=rang, bins=binsxy) 

这里的H包含直方图的直方图值和xedges, yedges边界。请注意，edges阵列的尺寸比相应尺寸的H的尺寸大一个。因此：

for x, y in numpy.argwhere(H == H.max()): 
    # center is between x and x+1 
    print numpy.average(xedges[x:x + 2]), numpy.average(yedges[y:y + 2]) 

Answer 2

这个问题应该可以帮助您：Python: get the position of the biggest item in a numpy array

您可以使用H.max()得到最大的价值，然后将其与H比较和使用numpy.nonzero找到所有最大值的位置：numpy.nonzero(H.max() == H)。这将比H.argmax()更昂贵，但您将获得所有最大值。