spkit.peak_interp

spkit.peak_interp(mX, pX, ploc)

Interpolate peak using parabolic interpolation

Interpolate peak values using parabolic interpolation

refined loction:
\[kp_{new} = kp + 0.5*(X[kp-1] - X[kp+1])/(X[kp-1] -2*X[kp]+X[kp+1])\]
refined value:
\[X_{new} = X[kp] + 0.25*(kp_new - kp)*(X[kp-1] - X[kp+1])\]
Parameters:
mXmagnitude spectrum
pXphase spectrum
ploc: locations of peaks
Returns:
interpolated - refined locations of frequencies and corresponding magnitude and phase:
iplocpeak location
ipmagmagnitude values
ipphase: phase values

Examples

#sp.peak_interp
import numpy as np
import matplotlib.pyplot as plt
import spkit as sp
X, fs, ch_names = sp.data.eeg_sample_14ch()
x = X[:,0]
mX, fr = sp.periodogram(x,fs=fs)
N = (len(fr)-1)*2
mX = 20*np.log10(mX)
ploc = sp.peak_detection(mX, thr=0)
iploc, _, _ = sp.peak_interp(mX, pX=mX, ploc=ploc)
peak_int = fs*iploc/N
print('Peaks Detected     :',fr[ploc],'Hz')
print('Peaks Interpolated :',peak_int.round(1),'Hz')
plt.figure(figsize=(7,4))
plt.plot(fr,mX)
plt.plot(fr[ploc],mX[ploc],'o')
plt.vlines(fr[ploc],ymin=-80,ymax=70,color='C1',ls='--',alpha=0.5)
plt.axhline(0,color='k',lw=1,ls='--')
plt.ylim([-80,70])
plt.xlim([0,fr[-1]])
plt.grid()
plt.ylabel('mX (dB)')
plt.title('Spectral Peaks: EEG Signal')
plt.show()
../../_images/spkit-peak_interp-1_00_00.png
#####################################
#sp.peak_interp
import requests
from scipy.io import wavfile
import numpy as np
import matplotlib.pyplot as plt
import spkit as sp

path = 'https://github.com/MLEndDatasets/samples/raw/main/HWD_HP_hum_1.wav?raw=True'
req = requests.get(path)
with open('downloaded_file.wav', 'wb') as f:
        f.write(req.content)
fs, x = wavfile.read('downloaded_file.wav')
t = np.arange(x.shape[0])/fs
mX, pX, N = sp.dft_analysis(x[int(1*fs):int(1.2*fs)])
fr = (fs/2)*np.arange(len(mX))/(len(mX)-1)
thr = 30
ploc = sp.peak_detection(mX, thr=thr)
iploc, ipmag, ipphase = sp.peak_interp(mX, pX, ploc)
peak_int = fs*iploc/N

print('Peaks detected      : ',fr[ploc], 'Hz')
print('Peaks Interpolated  : ',peak_int, 'Hz')

plt.figure(figsize=(7,4))
plt.plot(fr,mX)
plt.plot(fr[ploc],mX[ploc],'.')
plt.axhline(thr,color='k',lw=1,ls='--')
plt.ylim([-80,80])
plt.xlim([-1,2500])
plt.grid()
plt.ylabel('mX (dB)')
plt.title('Spectral Peaks: Audio')
plt.show()
../../_images/spkit-peak_interp-1_01_00.png