I got peak indexes on detailed coefficient5. However, I need to get indexes on signal. How do I need to transform index?
The individual recordings are each 10 hours in duration, and contain two ECG signals each sampled at 250 samples per second with 12-bit resolution over a range of ±10 millivolts. ambulatory ECG recorders with a typical recording bandwidth of approximately 0.1 Hz to 40 Hz.
import numpy as np
import pywt
import wfdb
import os
import matplotlib.pyplot as plt
def soft_threshold(signal,std):
# threshold = std *np.sqrt(2 * np.log(len(signal))) #OK
threshold =0.8*max(signal)
return np.sign(signal) * np.maximum(np.abs(signal) - threshold, 0)
def denoising(signal):
(cA10,cD10,cD9,cD8,cD7,cD6,cD5,cD4,cD3,cD2,cD1)= pywt.wavedec(signal, 'db10', level = 10,mode = 'symmetric') #분해해 bior1.3
std = np.median(np.abs(cD1)) / 0.6745
# 기저선 변동 제거를 위해 cA10,cD10 없앰
# 고주파 오염 노이즈원 및 전력간섭 제거를 위 cD1 ~ cD4 제거
# np.zeros_like는 변수의 사이즈 만큼의 0로 가득찬 numpy array를 배출, np.zeros는 array 사이즈만 넣지만 like은 array 넣음
for coeffs in [ cA10 ,cD10,cD4,cD3,cD2,cD1]:
#print(coeffs)
coeffs = np.zeros_like(coeffs)
denoised_coeffs = [cA10]+ [cD10,cD9,cD8,cD7,cD6,cD5,cD4,cD3,cD2,cD1]
denoised_signal = pywt.waverec(denoised_coeffs, 'db10') #디노이징 된 신호 재건하기
# cD5 thresholding
denoised_coeffs[6] = soft_threshold(denoised_coeffs[6],std)
return denoised_signal,denoised_coeffs[6]
def detect_r_peaks(modified_cD5):
# qrs 값 증폭하고 양수로 만들기 위해 제곱
cD5_s = np.square(modified_cD5)
# cD5를 10단계로 분해
(A10,D10,D9,D8,D7,D6,D5,D4,D3,D2,D1) = pywt.wavedec(cD5_s, 'db10', level=10, mode = 'symmetric')
# D1 ~ D5 제거
for coeffs_d5 in [D5,D4,D3,D2,D1]:
coeffs_d5 = np.zeros_like(coeffs_d5) #
# reference window에서 R 피크 찾기
Rf_window = pywt.waverec((A10,D10,D9,D8,D7,D6,D5,D4,D3,D2,D1), 'db10')
# 임계값 설정
T_window = Rf_window.mean() /4 # /는 소수점 반환, //는 round down 하고 정수 반환
# reference window에서 임계값보다 높은 값 가지면 r피크 영역
# r_peaks = []
# for j in range(1, len(Rf_window) - 1):
# if Rf_window[j] > T_window:
# r_peaks.append(j)
##
r_peaks = []
for j in range(1, len(Rf_window) - 1):
if Rf_window[j] > T_window and Rf_window[j] > Rf_window[j-1] and Rf_window[j] > Rf_window[j+1]:
r_peaks.append(j)
return r_peaks,Rf_window
#sampling_rate = 250
ecg_20s = ecg[0: 5000]
#remove dc
ecg_non_dc = ecg_20s - np.mean(ecg_20s)
#denoising
ecg_denoised,cD5 = denoising(ecg_non_dc)
r_peaks,Rf_window = detect_r_peaks(cD5)
#plot
labels = ['denoised signal','D5 Coefficient after Soft Thresholding T=300','D5 Coefficient after Squaring and D1-D5 Band Cancellation','ECG Signal after R-peak Point Detection']
n =4
#for (i,j) in zip(range(n),[ecg_denoised,cD5,Rf_window,r_peaks]):
# plt.subplot(n,1,i+1)
# plt.plot(j) #5초
# plt.title(labels[i])
# plt.tight_layout()
plt.plot(Rf_window)
plt.plot(r_peaks[:99],Rf_window[r_peaks[:99]],'ro')
What I want to get is the exact point of R peaks.