CF84 - 4 months ago 35

Python Question

I know this question falls within the category of peak detection and there are answers available, but I think my problem is pretty simplistic and refers to a proof of principle.

Say I generate multiple,

`nxn`

`nxn`

`myarray=`

array([[ 0.82760819, 0.82113999, 0.811576 , 0.80308705, 0.81231903,

0.82296263, 0.78448964, 0.79308308, 0.82160627, 0.83475755,

0.8580934 , 0.8857617 , 0.89901092, 0.92479025, 0.91840606,

0.91029942, 0.88523943, 0.85798491, 0.84190422, 0.83350783,

0.83520675],

[ 0.84971526, 0.84759644, 0.81429419, 0.79936736, 0.81750327,

0.81874686, 0.80666801, 0.82297348, 0.84980788, 0.85698662,

0.87819988, 0.89572185, 0.89009723, 0.90347858, 0.89703473,

0.90092666, 0.88362073, 0.86711197, 0.84791422, 0.83632138,

0.83685225], ...] #you can generate any nxn array of random values

Now let's normalize them:

`peak_value=myarray.max()`

norm_array=myarray/peak_value

And proceed to locate the

`(x,y)`

`from collections import Counter`

longit=[] #x values of the peak

latit=[] #y values of the peak

for x in range(myarray.shape[0]):

for y in range(myarray.shape[1]):

if norm_array[x][y]==1:

longit.append(x)

latit.append(y)

x=numpy.array(longit)

y=numpy.array(latit)

c=zip(x,y)

temp=Counter(elem for elem in c) #Counts the number of peaks in each (x,y) point in the 11x11 grid

d=dict(Counter(temp)) #Has the shape d={(x,y): number of peaks}

Now this is just a single realization of the 2D numpy array. Given multiple arrays, the question is:

Answer

In C/C++ it is considered dangerous to do `==`

with floating point numbers.

If there are no multiple exactly same pikes, you can use `numpy.argmax`

:

```
a = random.rand(13, 11)
idx_unrolled = argmax(a)
the_peak = a[idx_unrolled/11, idx_unrolled%11]
```

If you need all pikes, you can get list of `i`

, `j`

indices using `numpy.where`

:

```
i, j = where(a > 0.99*the_peak)
```

Use required number of `9`

to tune the margin. For single precision floating points it is not than close to `1`

.

The best way may be something like [http://stackoverflow.com/a/19141711/774971 ]:

```
i, j = where(a > (1.0 - 5*np.finfo(a.dtype).eps)*the_peak)
```