单参总结——心电信号

小懒懒️

于 2024-08-31 22:25:56 发布

阅读量428

点赞数 10

分类专栏：五参学习文章标签： stm32 学习笔记 c语言

本文链接：https://blog.csdn.net/2301_76901677/article/details/141216816

版权

五参学习专栏收录该内容

1 篇文章 0 订阅

订阅专栏

心电波形的滤波及心率计算

心电测量原理

心脏通过一系列电活动驱动其收缩和放松，这些电活动可以通过电极在身体表面检测到。心脏的电活动主要包括：

去极化：心脏的电信号从心房开始，经过心室，导致心脏的肌肉收缩。
复极化：心脏肌肉放松，准备下一次收缩。

心电波形

心电图（ECG）波形中的各种波代表了心脏电活动的不同阶段。标准的心电图波形包括P波、QRS波群和T波。每种波形反映了心脏电活动的特定部分。
在这里插入图片描述

P波

定义：P波是心电图中的第一个波，代表心房的去极化过程。
特点：
形状：通常为圆润而且小的正向波（在大多数导联中）。
持续时间：一般小于0.12秒（120毫秒）。
振幅：通常在0.1到0.3毫伏（mV）之间。
意义：P波的出现表示心房正在去极化，即心房肌肉准备收缩。如果P波缺失或形状异常，可能指示心房电活动的异常，如心房颤动或心房扑动。

QRS波群

定义：QRS波群代表心室的去极化过程，是心电图中最显著的波形。
组成：
Q波：QRS波群的第一个向下的波（负波），并不是所有心电图中都有Q波。 R波：QRS波群的第一个向上的波（正波），通常最为显著。
S波：R波后的向下波（负波），跟在R波后面。
特点：
形状：QRS波群通常为尖锐而且幅度较高。
持续时间：正常情况下，QRS波群的持续时间小于0.12秒（120毫秒）。
振幅：通常较大，尤其是R波。
意义：QRS波群的出现表示心室肌肉的去极化，即心室收缩。如果QRS波群的形状、持续时间或幅度异常，可能指示心室的电活动异常，如心室肥厚、心室早搏或心室传导阻滞。

T波

定义：T波代表心室的复极化过程，即心室电活动恢复到静息状态。
特点：
形状：通常为圆润的正向波（在大多数导联中），但在某些导联中也可能是负向的。
持续时间：比P波和QRS波群更长，但通常不具体测量。
振幅：通常在0.1到0.5毫伏（mV）之间。
意义：T波的出现表示心室正在复极化。如果T波的形状、方向或幅度异常，可能指示心室复极化过程的异常，如心肌缺血、心肌梗死或电解质失衡。

U波

定义：U波是心电图中有时出现的波，紧接在T波之后。
特点：
形状：通常较小且不总是清晰可见。
振幅：较小，通常不超过0.1毫伏（mV）。
意义：U波的意义尚不完全清楚，但可能与心室复极化过程中的一些细微变化有关。异常的U波可能与电解质失衡或药物效应有关。

综上所述，心电图的各个波形反映了心脏不同部位的电活动情况：

P波：心房去极化。
QRS波群：心室去极化。
T波：心室复极化。
U波（如果存在）：心室复极化的细微变化。

标准双极导联

标准双极导联是心电图（ECG）中常用的配置之一。它使用两个电极来记录心脏的电活动，并且每个导联都是双极的，即通过两个电极之间的电位差来测量心脏信号。以下是常见的标准双极导联配置：

I导联
电极位置：左手腕（LA）和右手腕（RA）。
信号来源：测量从右手腕到左手腕的电位差。
视角：主要监测心脏在水平面上的电活动。
II导联
电极位置：左脚踝（LL）和右手腕（RA）。
信号来源：测量从右手腕到左脚踝的电位差。
视角：提供心脏电活动的另一水平面视角，常用于观察心脏的整体电活动。
III导联
电极位置：左脚踝（LL）和左手腕（LA）。
信号来源：测量从左手腕到左脚踝的电位差。
视角：补充第二导联的电活动视角。

标准双极导联的特点是能够广泛地反映心脏的大致情况，如后壁心肌梗死、心律失常等，在导联Ⅱ或导联Ⅲ中可记录到清晰的波形改变。但标准导联只能说明两个肢体间的电位差，不能记录到单个电极处的电位变化。

心电波形处理

在心电图（ECG）信号的测量和分析中，硬件和软件滤波器用于提高信号的质量，去除噪声，增强信号的可读性。

心电测量电路设计

在这里插入图片描述
将原理图绘制好后进行PCB打板和贴片，通过板子实测模拟器发出的心电信号，打印出来得到原始心电波形：

这是心电信号进行归一化处理后得到的图形，能够清晰地看到不同频率在心电信号中的占比情况。图中可见一些干扰需要对其进行算法处理。

软件滤波

软件滤波是指在信号采集之后，通过计算机算法对信号进行处理。这种处理可以在信号记录后的后期分析中进行，本次实践中，只完成了50Hz陷波器、IIR低通滤波器、IIR高通滤波器、滑动窗口平滑器。

50Hz陷波器

50Hz陷波器的本质在于其能够精准地滤除特定频率（50Hz）的干扰信号，而对其他频率的信号影响最小。它是一种带阻滤波器，设计用来降低或消除中心频率为50Hz的干扰，同时尽可能保留其余信号。
关于50Hz陷波器的理论原理详见陷波器及其算法（基于C语言）
以下为具体的设计步骤：
在找到这里插入图片描述找到滤波器设计器进行滤波器设计
在这里插入图片描述
设计好滤波器后，找到相应的系数

Fs是采样频率，这里是500Hz。上述步骤是作为设计一个50Hz陷波器的步骤参考，最终得到的系数非下面代码所使用的系数，下面代码所使用的Fstop1和Fstop2的频率与上图不一致因此系数不一样。

//50Hz陷波器
static IIR_2OR_DATA_DEF notching_data = {
    Notching_filter_50Hz_GAIN,
    Notching_filter_50Hz_B0,       
    Notching_filter_50Hz_B1,       
    Notching_filter_50Hz_B2,       
    Notching_filter_50Hz_A1,        
    Notching_filter_50Hz_A2,  
    0, 0, 0, 0 
};
//陷波滤波器
#define Notching_filter_50Hz_GAIN  0.888856976280745336715938265115255489945 
#define Notching_filter_50Hz_B0      1
#define Notching_filter_50Hz_B1  -1.61816175199937850592846189101692289114 
#define Notching_filter_50Hz_B2      1
#define Notching_filter_50Hz_A1  -1.438314362015320924115258094388991594315
#define Notching_filter_50Hz_A2     0.777713952561490673431876530230510979891
double iir_50Hz_2or_func(IIR_2OR_DATA_DEF *filter_date, double target)
{
     //
     // w0 = x(0) - A1 * W1 - A2 * W2
     //
    filter_date->filter_W0 = (target) - filter_date->coeff_A1 * filter_date->filter_W1 \
                           - filter_date->coeff_A2 * filter_date->filter_W2;

    //
    // Y(0) = Gain * (B0 * W0 + B1 * W1 + B2 * W2)
    //
    filter_date->filter_out = filter_date->coeff_GAIN * ( filter_date->coeff_B0 * filter_date->filter_W0 \
                            + filter_date->coeff_B1 * filter_date->filter_W1 \
                            + filter_date->coeff_B2 * filter_date->filter_W2 );

    //
    // Update last data
    //
    filter_date->filter_W2 = filter_date->filter_W1;
    filter_date->filter_W1 = filter_date->filter_W0;

    //
    // Return Value
    //
    return(filter_date->filter_out);
}

效果展示

在这里插入图片描述
从归一化图中可以看出，50Hz工频干扰滤除的效果不错，但高于50Hz的频率干扰没有滤除，需要进行IIR低通滤波。

IIR低通滤波器

IIR低通滤波器的作用是允许低频信号通过，同时抑制高频信号。关于原理详细可见嵌入式常用的算法 - 二阶IIR低通滤波器

static u16 lowFilters_IIR(u16 data)
{
  static double x1=0,x2=0,y,y1=0,y2=0;

  y=0.064591*(data +2.0*x1+1.0*x2)+1.401314*y1-0.659679*y2;
  x2=x1;
  x1=data;
  y2=y1;
  y1=y;
  return y;
}

对于滤波器的系数同50Hz陷波器步骤可得，只不过需要改变滤波器类型和截止频率，根据具体要求设计相应功能的滤波器。

效果展示

在这里插入图片描述

由图可见，50Hz后的干扰基本滤除了。

IIR高通滤波器

基线漂移通常表现为信号中的低频偏移或漂移。IIR高通滤波器可以去除基线漂移，因为它能够抑制低频信号成分。高通滤波器允许高频成分通过，同时对低频成分进行衰减，从而去除这些低频漂移，使信号的基线更稳定。
关于IIR高通滤波器这部分具体参考了一篇文献
在这里插入图片描述
这篇文献中具体展示了IIR高通滤波器对基线漂移的抑制效果

static u16 highFilters_IIR(u16 data)
{
  static double x1=0,x2=0,y,y1=0,y2=0;
     
  y=0.992173*(data -1.999993*x1+ 1.0*x2)+1.984303*y1-0.984382*y2;   //IIR4阶二型切比雪夫 0.55   30dB
  x2=x1;
  x1=data;
  y2=y1;
  y1=y;
  return y+2000;  
}

由于处理后打印出来的数据比较小所以加2000可以放大波形，当然不加也没有太大关系。

效果展示

在这里插入图片描述

滑动窗口平滑器

滑动窗口平滑器（也称为滑动平均滤波器）是一种常用的信号处理技术，主要用于减少数据中的噪声和波动，平滑波形。

 u16 sliding_average_filter(u16 value)
{
  static int data_num = 0;
  static u16 output = 0;
	int i;
  float sum = 0;

  if (data_num < 10) //不满窗口，先填充
  {
    buffer[data_num++] = value;
    output = value; //返回相同的值
  }
  else
  {
    memcpy(&buffer[0], &buffer[1], (int)(10 - 1) * 4); //将1之后的数据移到0之后，即移除头部
    buffer[10 - 1] = value;                       //即添加尾部

    for (i = 0; i < 10; i++) //每一次都计算，可以避免累计浮点计算误差
    {    
      sum += buffer[i];
    }
    output = sum / 10;
  }

  return output;
}

效果展示

在这里插入图片描述

心率计算

u16 ECG_Calculate(u16 data[Len] ,u16 len)
{
  u16 heart_rate = 0;   //计算得出的心率
  u16 i,j,k,t;
  u16 Max_data = 0;  //数组中的最大值
  u16 Min_data = 0; //数组中的最小值
  u16 T = 0;    //阈值
  u8 Feak_Flag = 0;  //波峰找到时标志位置1
  u16 Feak_Num = 0;   //波峰个数
  u16 Feak_Index[30];   //Feak是储存波峰的索引数组
  u16 Index_num[30];    //相邻两个峰值之间的数据个数
  u16 n;
  u8 swapped = 0;
  
  //第1步找出求出阈值T
  Max_data=data[0];
  Min_data=data[0];
  for(i = 0;i < len;i++)
  {
    if(data[i] > Max_data)
      Max_data = data[i];
    if(Min_data > data[i])
      Min_data = data[i];
  }
  T = (Max_data - Min_data) * 0.9 + Min_data;  //求出阈值T

  //第2步找出数组中所有大于阈值T的索引，并存入新的数组中
  for(i = 0 , k = 0 ;i < len;i++)
  {
    if(data[i] > T)
    {
      Feak_Flag = 0;
      if((data[i] > data[i-1]) || (data[i] == data[i + 1]))
      {
       
        Feak_Flag = 1;   //初步判定该数据为峰值
        if((i + j) > len)   //超出数组长度退出循环
          break;           
        for(j = 1; j < 6; j++)
        {
          if(data[i] < data[i + j])  //判断前一个数据是否大于后一个数据
          {
            Feak_Flag = 0;          //该数据并非峰值，将标志位置零并跳出循环
            break;
          }
        }
        if(Feak_Flag == 1)  //最终判定该数据是波峰
        {
          Feak_Index[k] = i;
          k += 1;
          Feak_Num += 1;
          i = i + n_min;  //提高代码运行速度，找到一个波峰索引后跳到下一个波峰索引的大概位置
        }
      }
    }
  }
 
  //第3步计算Feak数组两相邻元素之间的差值，即两波峰之间的数据个数
  if(Feak_Num > 1)
    for(i = 1 ,j = 0; i < Feak_Num; i++ ,j++)
    {
      Index_num[j] = Feak_Index[i] - Feak_Index[i - 1]; 
    }
  else
    Index_num[0] = Feak_Index[0];
  
  //第4步将Index_num中的数据按从小到大的顺序进行排序
  if(Feak_Num > 1)
    for (i = 0; i < Feak_Num - 1; i++)
    {
      swapped = 0; // 标记是否发生了交换
      for (j = 0; j < Feak_Num - i - 1; j++)
      {
        if (Index_num[j] > Index_num[j + 1])
        {
          // 交换 Index_num[j] 和 Index_num[j + 1]
          t = Index_num[j];
          Index_num[j] = Index_num[j + 1];
          Index_num[j + 1] = t;
          swapped = 1; // 记录发生了交换
        }
      }
      // 如果在这一轮遍历中没有发生交换，说明数组已经有序
      if (!swapped)
        break;
    }
      
  //第5步取中值作为R-R波之间的数据个数,计算心率
    if(Feak_Num >= 2)
    {
      if((Feak_Num % 2) == 0)  //波峰数量为偶数
        n = Index_num[(Feak_Num/2)-1];
      else                     //波峰数量为奇数
        n = Index_num[((Feak_Num+1)/2)-1];
    }
    else   //波峰数量小于2，为1或0
        n = Index_num[0];

  heart_rate = 30000 / n;

  return heart_rate;
}