Skip to main content
 首页 » 编程设计

python之使用导数确定二维数组中局部最大值的坐标

2025年05月04日15hnrainll

我有一个 fits 图像,我试图在我的图像中找到 local maxima 的坐标,但到目前为止我还不能完全让它工作。 My image可以在这里找到。 到目前为止我所拥有的是

import numpy as np 
import scipy.nimage as ndimage 
from astropy.wcs import WCS 
from astropy import units as u 
from astropy import coordinates as coord 
from astropy.io import fits 
import scipy.ndimage.filters as filters 
from scipy.ndimage.filters import maximum_filter 
hdulist=fits.open("MapSNR.fits") 
#reading a two dimensional array from fits file 
d=hdulist[0].data 
w=WCS("MapSNR.fits")  
idx,idy=np.where(d==np.max(d)) 
rr,dd=w.all_pix2word(idx,idy,o) 
c=coord.SkyCoord(ra=rr*u.degree, dec=dd*u.degree) 
#The sky coordinate of the image maximum 
print c.ra 
print c.dec 

这就是我如何找到图像的全局最大值,但我想获得局部最大值的坐标,它具有大于三的意义.

我在网上查了一下发现是this following answer在我的情况下不能正常工作。 更新:这个功能我用过

def detect_peaks(data, threshold=1.5, neighborhood_size=5): 
 
  data_max = filters.maximum_filter(data, neighborhood_size) 
  maxima = (data == data_max) 
  data_min = filters.minimum_filter(data, neighborhood_size) 
  diff = ((data_max - data_min) > threshold) 
  maxima[diff == 0] = 0 # sets values <= threshold as background 
  labeled, num_objects = ndimage.label(maxima) 
  slices = ndimage.find_objects(labeled) 
  x,y=[],[] 
  for dy,dx in slices: 
    x_center = (dx.start + dx.stop - 1)/2 
    y_center = (dy.start + dy.stop - 1)/2 
    x.append(x_center) 
    y.append(y_center) 
  return x,y 

我想找到一种使用更好方法的方法,例如数组中的导数或分而治之的方法。我会推荐更好的解决方案。

请您参考如下方法:

所以我有了这个,使用 skimage 自适应阈值。希望对您有所帮助:

原创

代码

from skimage.filters import threshold_adaptive 
import matplotlib.pyplot as plt 
from scipy import misc, ndimage 
import numpy as np 
 
im = misc.imread('\Desktop\MapSNR.jpg') 
 
# Apply a threshold 
binary_adaptive = threshold_adaptive(im, block_size=40, offset=-20).astype(np.int) 
# Label regions and find center of mass 
lbl = ndimage.label(binary_adaptive) 
points = ndimage.measurements.center_of_mass(binary_adaptive, lbl[0], [i+1 for i in range(lbl[1])]) 
 
for i in points: 
    p = [int(j) for j in i] 
    binary_adaptive[i] += 5 
 
plt.figure() 
plt.imshow(im, interpolation='nearest', cmap='gray') 
plt.show() 
 
plt.figure() 
plt.imshow(binary_adaptive, interpolation='nearest', cmap='gray') 
plt.show() 

输出

改变阈值的参数将对在哪里找到局部最大值以及找到多少个最大值有很大影响。