59 lines
1.4 KiB
Python
59 lines
1.4 KiB
Python
import numpy as np
|
|
import matplotlib.pyplot as plt
|
|
from scipy import signal
|
|
|
|
from signalprocess import SignalProcess
|
|
|
|
# Read data from the video processing pipeline
|
|
sig = SignalProcess(path="out.tsv")
|
|
# Uncomment to see fft plot
|
|
# sig.fft()
|
|
|
|
# Applies the butter filter from the Signal Process class
|
|
data = sig.butter_filter()
|
|
|
|
|
|
plt.plot(sig.times, data)
|
|
plt.plot(sig.times, sig.ave_counts)
|
|
# Uncomment to see raw waveform
|
|
# plt.plot(sig.times,sig.counts)
|
|
|
|
# Peak finding algo provided by scipy package
|
|
peaks, descript = signal.find_peaks((-1) * data, prominence=0.2, distance=30)
|
|
|
|
# Plot the calculated peaks
|
|
plt.scatter(
|
|
sig.times[peaks],
|
|
data[peaks],
|
|
)
|
|
|
|
# Array of peak to peak distances
|
|
distances = []
|
|
# Array of the peak heights as calculated in the plot
|
|
heights = []
|
|
for i in range(len(peaks)):
|
|
heights.append(data[peaks[i]])
|
|
|
|
for i in range(1, len(peaks) - 1):
|
|
distance = (peaks[i] - peaks[i - 1]) / 30.0
|
|
distances.append(distance)
|
|
|
|
# Basic stats
|
|
print("distances: median={}, std={}".format(np.median(distances), np.std(distances)))
|
|
print("heights: median={}, std={}".format(np.median(heights), np.std(heights)))
|
|
|
|
# Standard deviation rejection of outliers
|
|
d = np.abs(heights - np.median(heights))
|
|
mdev = np.median(d)
|
|
s = d / mdev if mdev else 0.0
|
|
peaks = peaks[s < 3.0]
|
|
|
|
# plot the non-rejected data on top of the old peak data
|
|
plt.scatter(
|
|
sig.times[peaks],
|
|
data[peaks],
|
|
)
|
|
|
|
# Draw that bitch
|
|
plt.show()
|