1D陣列中的峯值計數

Question

我想讀取一個.txt文件，並從中獲取數據，然後想出一個方法來找出某個塊中有多少個峯值。我有一臺蓋革計數器，它將運行一段時間，其輸出將記錄在一個文件中，我將稱之爲「NuclearCounting.txt」。你可以看到第一百行here。你可以看到平均值約爲1.7，很少會達到〜3或下降到〜1。

然後，我想要將數據分成幾分鐘或一秒（本例中選擇秒）的數據塊，並分析每個時間間隔內有多少「事件」發生。我通過調用一個「事件」數據的最大值或最小值來做到這一點。由於背景輻射有一些噪音，所以我想在那裏有一條線，表示「峯值」必須與平均值相差3秒。如果可能的話。

這裏是到目前爲止的代碼：

import numpy as np 
import matplotlib.pyplot as plt 
from scipy import signal 
import scipy.stats 
x1 = [] 
y1 = [] 
f = open("NC.txt", "r") 
for line in f: 
    line = line.strip() 
    parts = line.split(",") #Separated by commas and spaces 
    time = float(parts[1]) #This column is the column recording time 
    value = float(parts[2]) #This one records the value the detector reads. 
    x1.append(time) 
    y1.append(value) 
f.close() 

xv = np.array(x1) #Time 
yv = np.array(y1) #Values recorded 

#Statistics 
m = np.mean(yv) 
d = np.std(yv) 
timetot = int(round(np.max(xv))) #This gives the total time in seconds to the nearest second. I use this to break it into chunks of one second of data each. 
#Counting events 
def chunk(l, n): 
    n = max(1, n) 
    return [l[k:k+n] for k in range(0, len(l), n)] 
t = chunk(xv,timetot) #Breaking the time into chunks of a second 
v = chunk(yv,timetot) #Doing the same for the values recorded by the detector 
i = [] 
for x in v: 
    peak = signal.find_peaks_cwt(x, np.arange(2,10)) #This is the problem area 
    i.append(peak) 
print(len(i)) 

陣列i被用於存儲的「峯」我得到的數量。我只想知道每個區間的i的長度。問題是python的輸出告訴我我有25000個峯值。這是沒有意義的。每秒鐘，我應該只能得到最大的數。它也只是給了我一個陣列，其中有條目完全是15.

有人可以請告訴我如何糾正？我不知道爲什麼scipy函數不應該工作，當然不是爲什麼它應該吐出15沒有數據輸入有這個價值。

謝謝。

Answer 1

代碼中的列表i不存儲每個塊的峯值數量，但是存儲峯值索引。

爲了計算峯值，您需要使用i.append(len(peak))。如果您還想限制峯值比標準差大3個標準偏差，您可以執行peak = np.sum(x[signal.find_peaks_cwt(x, np.arange(2,10))] > 3*d)。

最後但並非最不重要print(len(i))將返回峯值列表的長度，即塊的數量。我想你想要的是找出每個塊的峯值數量。所以，總結一下，我認爲這是你想要的。

for x in v: 
    peak_index = signal.find_peaks_cwt(x, np.arange(2,10)) 
    peak_sum = np.sum(x[peak_index] > 3*d) 
    i.append(peak_sum) 
for p in i: 
    print(p)