[发明专利]一种大气颗粒物PM2.5数密度检测的方法

摘要

本发明提供了一种大气颗粒物PM2.5数密度的检测方法，其采用激光光源输出的红外光，照射含有PM2.5的气体，因散射效应，形成散射光，用光学部件，将各散射光汇聚，形成干涉光，采用光电探测器等测量干涉光强，为消除噪声污染，用自适应FIR滤波器进行信号处理，并计算干涉光强平均值、方差，依据该平均值、方差与散射光的总光强、形成散射光的PM2.5总数目的相关性原理，按系统参数估计算法，完成对散射光强、PM2.5总数目的估计，根据已知测试室容积和PM2.5总数目的估计值，计算被测气体的PM2.5数密度，实现大气PM2.5的检测，其检测原理新颖、独特，其算法基于递推最小二乘RLS原理，无需矩阵求逆运算，运算量低，性能稳定，其系统结构简洁，成本低，可靠性高。

主权项

一种大气颗粒物PM2.5数密度检测方法，其特征在于，其包括以下步骤：步骤1：含有多种粒径颗粒物的被测空气，经过不同孔径的过滤器分离后，形成一股含有PM2.5及更小颗粒物的气流，并进入气体测试室；步骤2：在气体测试室中，用激光光源输出的红外光照射被测气流，一部分红外光与气流中颗粒物发生散射，形成散射光，另一部分红外光不与颗粒物发生散射，直接穿透被测气流，形成穿透光，采用光学部件，将散射光和穿透光汇聚在一起；步骤3：汇聚在一起的散射光和穿透光，发生干涉效应，并输出干涉光；步骤4：采用光电探测器测量以上干涉光的强度信号，其包含各散射光的自相关强度、穿透光的自相关强度、散射光与穿透光之间的互相关强度、各散射光之间的互相关强度四项内容，该干涉光的强度信号是随机信号，其概率分布密度函数与形成散射光的颗粒子总数目N相关；步骤5：用A/D转换器对以上的干涉光总强度信号，以采样周期为T，进行离散化抽样，得到抽样信号序列y(n)，用y(n)表示在nT时刻的抽样信号，n为正整数；步骤6：为了从被噪声污染的抽样信号中，恢复出原始的干涉光强信号，采用递推最小二乘RLS算法，设计一个M阶的自适应FIR滤波器，对抽样信号y(n)，按下式进行滤波：x(n)=Σi=0M-1wi(n)y(n-i)=WT(n)Y(n)=YT(n)W(n)]]>其中，wi(n)为nT时刻的M阶FIR滤波器滤波系数，i＝0,1,...,M‑1，W(n)为nT时刻滤波系数向量：W(n)＝[w0(n),...,wM‑1(n)]T，而x(n)为滤波后的信号序列，x(n)可视为在nT时刻的原始的干涉光强信号估计值，Y(n)为nT时刻所对应的信号向量，该向量由M个抽样信号组成，其如下：Y(n)＝[y(n),y(n‑1)...,y(n‑M+1)]TWT(n)、YT(n)分别为W(n)和Y(n)的转置：WT(n)＝[w0(n),...,wM‑1(n)]，YT(n)＝[y(n),y(n‑1),...,y(n‑M+1)]；步骤7：设滤波后的信号序列x(n)的序号指针为n0，从n0＝0开始，连续取最前面的二段长度为L的信号，L为正整数，用信号向量X1、X2表示：X1＝[x(1),...,x(L)]TX2＝[x(L+1),...,x(2L)]TX1、X2分别对应在[0‑LT]、[LT‑2LT]时间内的干涉光强度信号，分别是信号X1、X2的转置形式：X1T＝[x(1),...,x(L)]X2T＝[x(L+1),...,x(2L)]按下式，计算X1、X2的平均值：X1‾=(Σk=1Lx(k))/L]]>X2‾=(Σk=1Lx(L+k))/L]]>L,M均为正整数，且L＞M；步骤8：按下式，计算X1、X2的方差：σ12=Σk=1L[x(k)-X1‾]2]]>σ22=Σk=1L[x(L+k)-X2‾]2]]>σ12和σ22分别对应在[0‑LT]和[LT‑2LT]时间内的干涉光强度的涨落:步骤9：根据干涉光强度在[0‑LT]、[LT‑2LT]时间内的平均值和方差：σ12、σ22，建立以下的参数C1、C2的矩阵方程：其中，a2为粒子散射光的总光强，N为形成散射光的粒子总数目，将C1、C2称为系统参数，这里假设，a2和N在[0‑LT]和[LT‑2LT]时间内保持不变，因此，在此时间内，系统参数C1、C2近似是一组常量；步骤10：通过以上方程的最小二乘解，得到[0‑LT]和[LT‑2LT]时间内的系统参数C1和C2的估计值：其中，P0＝(A0TA0)‑1，A0T为矩阵A0的转置，(A0TA0)‑1为A0TA0的逆矩阵，令将在[0‑LT]和[LT‑2LT]时间内的系统参数C1、C2，分别用表示；步骤11：调整信号序列x(n)的序号指针n0，令n0＝n0+2L，使序号指针，指向下一段长度为L的滤波后的信号；步骤12：从当前的序号指针n0开始，取一段长度为L滤波后的信号：X3＝[x(n0+1),...,x(n0+L)]T，X3对应于在[n0T‑(n0+L)T]时间内的干涉光强度信号，按下式，计算X3的平均值及方差σ32：X3‾=(Σk=1Lx(n0+k))/Lσ32=Σk=1L[x(n0+k)-X3‾]2,]]>用新生成的平均值及方差构成信号更新向量dT及b(1)：dT=2X3‾σ32-2σ32,b(1)=1,]]>建立在[n0T‑(n0+L)T]时间内的系统参数C1、C2的矩阵方程：dTC1C2=b(1),b(1)=1;]]>步骤13：根据所构成的信号更新向量dT及b(1)，采用递推最小二乘RLS算法，按下式，完成在[n0T‑(n0+L)T]时间内的对系统参数C1、C2的进行估计，并分别用表示之：P1=P0-P0ddTP01+dTP0dC1(1)C2(1)=C1(0)C2(0)+(b(1)-dTC1(0)C2(0))P1d,]]>其中，向量d与dT互为转置关系，与分别对应于本次迭代运算以前与之后的参数C1、C2；步骤14：根据以上[n0T‑(n0+L)T]时间内的系统参数C1和C2，按下式，计算气体测试室在该时间内，散射光的总光强a2及形成散射光的粒子总数目N:a2=C1,N=C12C2;]]>步骤15：根据已知的气体测试室的容积V和以上求得的粒子总数目N，能计算出在[n0T‑(n0+L)T]时间内，被测气体的颗粒物数密度ρ，即单位体积的气体中PM2.5粒子数量：ρ=NV]]>步骤16：根据递推最小二乘RLS算法要求，按下式，重置迭代运算各单元的初值：P0=P1C1(0)C2(0)=C1(1)C2(1);]]>步骤17：调整信号序列X的序号指针n0，令n0＝n0+L，使序号指针指向下一段长度为L的滤波后的信号；步骤18：重复步骤12‑17，如此循环，完成在起始时刻n0＝0之后的，任意时间[n0T‑(n0+L)T]内的系统参数C1、C2的估计，进而，得到任意时间[n0T‑(n0+L)T]内的散射光的总光强a2及形成散射光的粒子总数目N，实现被测气体在起始时刻n0＝0之后的，任意时间[n0T‑(n0+L)T]内的颗粒物数密度的检测，其中，n0＝0,2L,3L,4L,5L,.......。