EMD关于信号的重建,心率提取

news/2024/5/19 0:56:24 标签: signal

关于EMD的俩个假设:

IMF 有两个假设条件:

  • 在整个数据段内,极值点的个数和过零点的个数必须相等或相差最多不能超过一 个;
  • 在任意时刻,由局部极大值点形成的上包络线和由局部极小值点形成的下包络线 的平均值为零,即上、下包络线相对于时间轴局部对称。

先安装pyEMD库 

from pyEMD import EMD  (报错)
执行pip uninstall pyEMD
pip install EMD-signal==1.4.0 -i https://pypi.tuna.tsinghua.edu.cn/simple/

代码部分:

from PyEMD import EMD
import numpy as np
import matplotlib
matplotlib.use('TkAgg')
import matplotlib.pyplot as plt
from scipy.signal import find_peaks
from scipy import signal

#读取信号数据
def readtxt(path):
    with open(path,'r') as f:
        str=f.readline()
        list = str.split(' ')
    list1=[];
    for i,x in enumerate(list):
        if ((i%5 == 0) or (i%5 == 1) or (i%5 == 2)) and x !='':
            list1.append(float(x))
    return list1


def pyem(lis):
    emd = EMD()
    IMFs = emd.emd(np.array(lis))
    print(len(IMFs), len(IMFs[0]), type(IMFs))


    #计算周期和频率
    imfs = emd.imfs

    # 估计瞬时频率和周期
    freqs = []
    periods = []
    for imf in imfs:
        if len(imf) > 1:
            # 计算频率
            sample_rate = 1 / (imf.argmax() / len(imf))
            freq = sample_rate / len(imf)
            print(freq)
            freqs.append(freq)

            # 计算周期
            period = 1 / freq
            print(period)
            periods.append(period)

    fig = plt.figure()
    ax = fig.add_subplot(len(IMFs) + 1, 1, 1)
    ax.plot(np.array(lis))
    for i in range(len(IMFs)):
        ax = fig.add_subplot(len(IMFs) + 1, 1, i + 2)
        ax.plot(IMFs[i])

    plt.show()

    lr = 0
    for i,s in enumerate(IMFs):
        if i>len(IMFs-1)/2+1:
            lr +=s
    return lr


#获取心率
def findPeaks(list):
    # x = electrocardiogram()[2000:4000]
    # jus=[1,2,3,4,10,1,2,3,4,21,1,2,2,3]
    #获取列表最小值,然后减去最小值
    # list_N = list[20000:30000]
    list_N = list

    avg = sum(list_N)/len(list_N);
    list_D=[]
    for i in range(len(list_N)):
        list_D.append(list_N[i]-avg)

    #列表转换数组
    y=np.array(list_D)
    #消除趋势线
    z=signal.detrend(y)
    #结果抽取200点,降频,然后再获取数据的脉率
    pl=200;
    fs=len(list_N)
    #参照值比
    BP = fs/pl;
    #进行趋势拟合
    x=signal.resample(z,pl)
    #获取最小值作为条件限制
    hu=min(x)
    peaks, _ = find_peaks(x, height=hu)
    # print(peaks)
    # 实际的心率值
    # print(len(peaks))
    ##获取相邻俩个峰值之间的点数,然后计算心率值
    for i,d in enumerate(peaks):  #打印查看脉搏波的数值
        print(peaks[i])

    a1 = peaks[0]
    a2 = peaks[1]
    a3= a2-a1
    #计算每个脉搏对应的点数
    R_point = a3*BP

    #以60为节点计算的数值
    rate=60*(500/R_point)
    # print("bass",bass)
    #总的点数除以每一个脉搏对应的点数,然后除以90秒对应的值
    # rate = (len(list)/R_point)/1.5

    plt.plot(x)
    plt.plot(peaks, x[peaks], "x")
    plt.plot(np.zeros_like(x), "--", color="gray")
    plt.show()

    return rate

if __name__ == "__main__":
    path = "../362a7e1de4dd484a9b4a3274a0e5a633_1648249928320.txt" #正常
    # path = "../a7c9bff53f2e4a70af7a9f641552507a_1706541122_1706564288403_887_1.txt"  #异常
    ll = readtxt(path)
    imf = pyem(ll[2000:10000])
    plt.plot(imf)
    plt.show()
    print(findPeaks(imf))

运行结果:

这是IMFS的分解图9个,从低频一直到高频

 因为最后一个是趋势项,我们将IMF[5]、IMF[6]、IMF[7]进行叠加,这几本接近我们的目标信号

然后对目标信号进行峰值提取:

总结:

信号分量的处理

通过经验模态分解(EMD)得到了信号的分量,可以进行许多不同的分析和处理操作,以下是一些常见的对分量的利用方向:

(1)信号重构:将分解得到的各个本征模态函数(IMF)相加,可以重构原始信号。这可以用于验证分解的效果,或者用于信号的重建和恢复。

(2)去噪:对于复杂的信号,可能存在噪声或干扰成分。通过分析各个IMF的频率和振幅,可以识别和去除信号中的噪声成分。

(3)频率分析:分析每个IMF的频率成分,可以帮助理解信号在不同频率上的振荡特性,从而揭示信号的频域特征。

(4)特征提取:每个IMF代表了信号的局部特征和振荡模式,可以用于提取信号的特征,并进一步应用于机器学习或模式识别任务中。

(5)信号预测:通过对分解得到的各个IMF进行分析,可以探索信号的未来趋势和发展模式,从而用于信号的预测和预测建模。

(6)模式识别:分析每个IMF的时域和频域特征,可以帮助对信号进行模式识别和分类,用于识别信号中的不同模式和特征。

(7)异常检测:通过分析每个IMF的振幅和频率特征,可以用于探测信号中的异常或突发事件,从而用于异常检测和故障诊断。

在得到了信号的分量之后,可以根据具体的应用需求选择合适的分析和处理方法,以实现对信号的深入理解、特征提取和应用。


http://www.niftyadmin.cn/n/5462290.html

相关文章

java多数据源几种实现方式以及demo

提示:多数据源实现方式、多数据源的使用场景。AbstractRoutingDataSource、DynamicDataSource框架、mybatisplus的Intercepter插件、java中多数据源的几种实现方式、mybatisPlus的插件实现多数据源 文章目录 前言一、多数据源的几种实现方式二、使用场景三、核心原理…

Taskflow:优先级任务(Prioritized Tasking)

Taskflow 支持给一个Task设置优先级,有助于在特定场景下的优化;总共有三种优先级:tf::TaskPriority::HIGH, tf::TaskPriority::NORMAL, 和 tf::TaskPriority::LOW。 对于一个并行Task集合(一组零依赖的Task,可并行执行…

滑动窗口算法详解及应用示例

引言: 滑动窗口算法是一种用于解决数组/字符串问题的有效技巧。它可以用来解决一系列问题,例如求解子数组/子字符串的最大值、最小值、平均值、和、特定条件下的个数等。在本篇博客中,我们将详细介绍滑动窗口算法的原理,并通过几个…

torch中的nonzero() 方法

其实这里的torch的nonzero() 和numpy的nonzero() 是一样的使用方法。 就是返回所有非0元素的index。一般用于CV中mask操作 首先定义一个5*4的tensor,其中为0的元素个数为3个,不为0的个数为17个 import torchXtorch.tensor([[1, 2, 0, 4],[5, 6, 7, 8…

关系型数据库和非关系型数据库介绍

好的,关系型数据库和非关系型数据库是两种不同类型的数据库,它们在数据存储、查询方式、性能、可扩展性等方面存在差异。以下是一些常见的关系型数据库和非关系型数据库的示例、特点和具体使用场景: 关系型数据库 1. MySQL:MySQL …

SQL109 纠错4(组合查询,order by..)

SELECT cust_name, cust_contact, cust_email FROM Customers WHERE cust_state MI UNION SELECT cust_name, cust_contact, cust_email FROM Customers WHERE cust_state IL ORDER BY cust_name;order by子句,必须位于最后一条select语句之后

2010-2021年各省碳排放测算数据(含原始数据+计算过程+结果)

2010-2021年各省碳排放测算数据(含原始数据计算过程结果) 1、时间:2010-2021年 2、指标:原煤(万吨)、原煤(万吨CO2)、焦炭(万吨)、焦炭(万吨CO2)、汽油(万吨)、汽油(万吨CO2)、煤油(万吨)、煤油(万吨CO2)、柴油(万吨)、柴油(万吨…

【隐私计算实训营006隐语PIR介绍及开发实践】

1. 隐语实现PIR总体介绍 隐匿查询(Private Information Retrieval PIR)定义 按服务器数量分类 单服务器方案(Single Server)多服务器方案(Multi-Server) 按查询类型分类 Index PIRKeyword PIR 隐语目前…