CN104198463A - 拉曼光谱预处理方法及系统 - Google Patents

Abstract

本发明公开了拉曼光谱预处理方法，包括小波去噪阶段、非对称最小二乘法去噪阶段、标准归一化处理阶段和三次平滑样条拟合处理阶段，其中所述非对称最小二乘法去噪阶段，通过设计专门的门电路替代预处理系统中的非对称最小二乘基线校正算法的软件处理方法，提高了拉曼光谱信号预处理系统的处理能力和预处理效果，在同样的处理时间条件下，通过更多次的迭代处理可以使拉曼光谱预处理算法的效果趋向于收敛到更佳的效果。这种思想能够广泛的被借鉴在其他的光谱信号预处理系统中。

Description

拉曼光谱预处理方法及系统

技术领域

本发明涉及信号处理方法，特别是一种拉曼光谱信号的预处理方法，同时还涉及一种拉曼光谱预处理系统。

背景技术

现有设计的拉曼光谱仪中，为了消除各种光谱信号噪声提高光谱信号信噪比，拉曼光谱仪中普遍包含了光谱信号预处理系统，但并没有对其工作方式做过详尽的报道和描述。而且这些研究成果由于都采用了软件方式实现，虽然能够满足拉曼光谱信号的预处理去噪、提高信噪比的需求，能够对一般拉曼光谱系统的设计起到一定的借鉴作用，但是涉及到产品的实时处理性能、更优的去噪效果方面还是有一些欠缺。随着拉曼光谱仪分析技术的广泛使用，需要一种具有实时处理能力、去噪效果更佳、低成本的拉曼光谱信号预处理实现方法。

发明内容

鉴于此，本发明提供一种高效、稳定的拉曼光谱信号预处理方法，同时提供一种拉曼光谱信号预处理系统。

为达到上述目的之一，本发明提供如下技术方案：拉曼光谱预处理方法，包括小波去噪阶段、非对称最小二乘法去噪阶段、标准归一化处理阶段和三次平滑样条拟合处理阶段，所述非对称最小二乘法去噪阶段具体包括以下步骤：

S1求解拉曼信号强度向量维度m；

S2生成m维的稀疏矩阵，并对该稀疏矩阵求2阶导生成矩阵D；

S3生成初始化m行1列的单位矩阵w；

S4以w向量为对角元素，生成维度为m行m列的带状稀疏矩阵W；

S5对(W+D’*D)做Cholesky分解，得到维度为m行m列矩阵C；

S6求C’*C矩阵的逆矩阵E；

S7向量w和向量y做点积生成向量R，其中y为常量；

S8做矩阵乘法E*R生成1行m列向量z；

S9对向量z和拉曼信号按位做’>’作生成新的w向量；

S10更新w向量，重复步骤S4～S10，通过多次迭代找到最优的光谱基线。

为达到上述目的之二，本发明提供如下技术方案：拉曼光谱预处理系统，包括小波去噪模块、非对称最小二乘法去噪模块、标准归一化处理模块和三次平滑样条拟合处理模块，所述非对称最小二乘法去噪模块包括中央控制模块、矩阵处理阵列、Cholesky求解模块和片上内存；每一次迭代过程，中央控制模块通过CHOL_RUN信号触发Cholesky求解模块执行，并且将对片上内存的读写访问总线交给Cholesky求解模块；当Cholesky求解模块执行完毕后通过CHOL_RDY信号反馈给中央控制模块；中央控制模块再将MATRIX_PROC_RUN信号拉高，使能矩阵处理阵列；矩阵处理阵列执行完毕后通过拉高MATRIX_PROC_RDY信号反馈给中央处理单元；中央处理单元执行下一次的基线校正迭代；通过多次迭代找到最优的光谱基线。

进一步，所述Cholesky求解模块用于对待求解矩阵进行最大值遍历搜索和C_ij＝C_ij-^a _ik*C_kj变换，i,j表示被搜索矩阵行列坐标，其中a_ik为k列元素，C_kj为k行元素，k＝1,2,3,4,······,m；并将生成的值更新到相对应坐标位置元素，将相对应最大值元素的行列坐标位置信息存放入FILO的存储器中；所述Cholesky求解模块再从FILO存储器中读取坐标信息，并依据此坐标信息进行矩阵的行列互换，并将生成的结果送到片上内存中。

进一步，所述矩阵处理阵列包括乘法单元、加法器和比较器；所述乘法单元进行两个矩阵的乘法运算，即受控制信号操作从内存中分别读取矩阵A的某一行与矩阵B的某一列，将行列元素分别相乘生成1024个乘积值；所述加法器对1024个乘积值做累加，完成后生成向量z，向量z与向量y通过比较器比较生成向量w，y为常量。

由于采用了以上技术方案，本发明具有以下有益技术效果：本发明提供了一种高效、稳定的拉曼光谱信号预处理方法及系统。通过本发明，能够实现对光谱信号的去噪滤波、消除荧光背景对光谱信号分析的干扰、消除噪声和数据量纲对光谱信号分析的影响提高光谱信号的数据信噪比；通过设计专门的门电路替代预处理系统中的非对称最小二乘基线校正算法的软件处理方法，提高了拉曼光谱信号预处理系统的处理能力和预处理效果，在同样的处理时间条件下，通过更多次的迭代处理可以使拉曼光谱预处理算法的效果趋向于收敛到更佳的效果。

附图说明

为了使本发明的目的、技术方案和优点更加清楚，下面将结合附图对本发明作进一步的详细描述，其中：

图1为拉曼光谱预处理流程图；

图2为拉曼光谱预处理任务流水线；

图3为IEEE-754浮点数表达方式图；

图4为IEEE-754浮点数表达方式精度分析图；

图5为本设计所采用的浮点数描述方式；

图6为本设计浮点表达方式精度分析图；

图7为最小二乘基线校正算法电路系统框图；

图8为片上内存存储区域划分；

图9为矩阵存储方式；

图10为最大值遍历搜索处理电路结构图；

图11为在最大值搜索过程中第K行数据存储图；

图12为矩阵处理阵列基本处理单元电路结构图。

具体实施方式

以下将结合附图，对本发明的优选实施例进行详细的描述；应当理解，优选实施例仅为了说明本发明，而不是为了限制本发明的保护范围。

图1为拉曼光谱预处理流程图；图2为拉曼光谱预处理任务流水线；图3为IEEE-754浮点数表达方式图；图4为IEEE-754浮点数表达方式精度分析图；图5为本设计所采用的浮点数描述方式；图6为本设计浮点表达方式精度分析图；图7为最小二乘基线校正算法电路系统框图；图8为片上内存存储区域划分；图9为矩阵存储方式；图10为最大值遍历搜索处理电路结构图；图11为在最大值搜索过程中第K行数据存储图；图12为矩阵处理阵列基本处理单元电路结构图。

拉曼光谱仪，利用从光信号传感器获取的拉曼光谱对物质进行分析，因而所得到光谱信号质量优劣会直接影响分析结果的准确。通常从光信号传感器获得的原始光谱存在一定的噪声，为了尽量的消减这类噪声、提高光谱信号的信噪比，本发明提供一种拉曼光谱预处理方法，包括小波去噪阶段、非对称最小二乘法去噪阶段、标准归一化处理阶段和三次平滑样条拟合处理阶段，所述小波去噪阶段采用基于小波变换的自适应阈值去噪算法实现拉曼原始光谱信息的高频噪声滤波；所述非对称最小二乘法去噪阶段去除荧光背景对光谱分析的干扰；所述标准归一化处理阶段，消除噪声和数据量纲对光谱分析的影响；所述三次平滑样条拟合处理阶段对归一化处理后的拉曼谱图数据，进行基于三次平滑样条拉曼光谱数据似合。

其中，所述非对称最小二乘法去噪阶段具体包括以下步骤：

S1求解拉曼信号强度向量维度m；

S2生成m维的稀疏矩阵，并对该稀疏矩阵求2阶导生成矩阵D；

S3生成初始化m行1列的单位矩阵w；

S4以w向量为对角元素，生成维度为m行m列的带状稀疏矩阵W；

S5对(W+D’*D)做Cholesky分解，得到维度为m行m列矩阵C；

S6求C’*C矩阵的逆矩阵E；

S7向量w和向量y做点积生成向量R；

S8做矩阵乘法E*R生成1行m列向量z；

S9对向量z和拉曼信号y按位做’>’作生成新的w向量；

S10更新w向量，重复步骤S4～S10，通过多次迭代找到最优的光谱基线。

对拉曼光谱预处理算法处理能力分析，预处理算法中含有大量的大矩阵数据乘法操作和多次的循环迭代，并且其中的部分操作数为浮点格式数据，这些特质在实际的处理中将耗费大量的处理器资源，因此为了满足系统性能要求改善拉曼光谱预处理效果，将整个预处理流程按照功能结构上划分为4级的流水结构，以提高整个系统的数据处理能力。其中小波去噪、标准归一化处理、三次平滑样条拟合处理三部分分别从软件上划分到多核处理器的单个核心上，依靠通用处理器进行处理；非对称最小二乘法由于含有大量的矩阵求逆、乘法、Cholesky操作，而且需要多次迭代，是影响拉曼光谱预处理系统的性能瓶颈，其性能的优良直接决定了拉曼光谱预处理系统的性能指标，甚至直接影响到所设计的拉曼光谱仪的处理性能指标。因此在本设计中通过基本门电路替代传统通用处理器来实现该部分拉曼光谱预处理算法，该方法即可以改善拉曼光谱预处理性能，又能通过更多次的迭代进行基线校正，使校正的基线收敛于最佳的基线位置，从而改善拉曼光谱仪的性能指标。

拉曼光谱预处理流水线如图2示。流水线级与级之间通过ping-pong内存块进行级联。在整个流水处理过程中的某一时刻，上级处理单元将处理的数据放入ping-pong内存块中的ping/pong块，而下一级处理单元从pong/ping块中取出待处理数据进行处理，如此方式使得在每个时刻上，流水线上的4个处理单元均在工作，加快了整个拉曼预处理系统的数据处理速率。

本设计中的最小二乘基线校正算法为了保证拉曼光谱去噪效果，所有算法均在浮点域中做原型验证，而基于门电路的实现系统是以布尔代数为基础的表征方式。若不对浮点域的数值做精度不衰减的定点化处理，将很难在数字系统中实现上述算法，即使勉强实现也会对最小二乘基线校正算法的效果和系统的性能产生较大的衰减。

根据IEEE-754标准所描述，32位浮点数的描述方式如图3示；

IEEE-754标准所能描述的浮点数的精度如图4所示，呈现出在实数轴上0值附近精度较高，离0值越远所描述的浮点数精度衰减较快的规律。

为了便于数字电路处理，同时保证算法的精度，在本设计中采用如图5的方式替代IEEE-754中的浮点数。32位中最高位为符号位，整数部分8位，浮点部分23位；

采用新的浮点数描述方式后，在整个实数域内，该描述方式所代表的浮点数精度是一个恒定的精度，所能代表的实数范围和精度为：

范围[-256,255.999999981]；精度0.000000119；

最小二乘基线校正算法中所有输入参数和光谱数据都采用上述的浮点数描述方式替换后，整个最小二乘基线校正算法的处理就在该数值域下进行。

与此同时本发明还提供一种拉曼光谱预处理系统，包括拉曼光谱预处理系统，包括小波去噪模块、非对称最小二乘法去噪模块、标准归一化处理模块和三次平滑样条拟合处理模块，

所述小波去噪模块采用基于小波变换的自适应阈值去噪算法实现拉曼原始光谱信息的高频噪声滤波；所述非对称最小二乘法去噪模块用于去除荧光背景对光谱分析的干扰；所述标准归一化处理模块用于消除噪声和数据量纲对光谱分析的影响；所述三次平滑样条拟合处理模块用于对归一化处理后的拉曼图谱进行基于三次平滑样条拉曼光谱数据拟合。

图7为最小二乘基线校正算法电路系统框图；所述非对称最小二乘法去噪模块包括中央控制模块、矩阵处理阵列、Cholesky求解模块和片上内存；每一次迭代过程，中央控制模块通过CHOL_RUN信号触发Cholesky求解模块执行，并且将对片上内存的读写访问总线交给Cholesky求解模块；当Cholesky求解模块执行完毕后通过CHOL_RDY信号反馈给中央控制模块；中央控制模块再将MATRIX_PROC_RUN信号拉高，使能矩阵处理阵列；矩阵处理阵列执行完毕后通过拉高MATRIX_PROC_RDY信号反馈给中央处理单元；中央处理单元执行下一次的基线校正迭代；通过多次迭代找到最优的光谱基线。

在片上内存中为最小二乘基线校正算法中常量和中间变量划分对应的存储区域如图8。对于稀疏矩阵或者矩阵中0值占大部分的矩阵采用图9所示存储方式，以节省存储空间、加快矩阵运算速度。

所述Cholesky求解模块用于对待求解矩阵进行最大值遍历搜索和C_ij＝C_ij-^a _ik*C_kj变换，i,j表示被搜索矩阵行列坐标，其中a_ik为k列元素，C_kj为k行元素，k＝1,2,3,4,······,m；并将生成的值更新到相对应坐标位置元素，将相对应最大值元素的行列坐标位置信息存放入FILO的存储器中；所述Cholesky求解模块再从FILO存储器中读取坐标信息，并依据此坐标信息进行矩阵的行列互换，并将生成的结果送到片上内存中。

图10为矩阵最大值搜索电路结构图。在本实施例中最大值遍历搜索的过程为：进行第K次遍历的过程中，首先从内存中读取第K行数据，并将第K行数据分别存储在如图11所示的“寄存器行”中，与MAX寄存器中最大值比较，在最大值比较过程中再预读取第K+1行的数据，如果所正在进行比较操作的行中有数据值大于MAX寄存器中值，那么就更新“含有矩阵最大元素的寄存器行数据”，直到完成整个迭代。

在本实施例中行列互换的方法为，以一个4行4列的矩阵为例，例如在0一次遍历矩阵搜寻最大值的迭代中找到最大值为a₃₂，这里对应元素的行列坐标指的是最大值a₃₂中行左边3，列坐标2，将坐标信息存入FILO中。

交换a₃₂所在的行列元素，是指将第3行元素与第0行元素互换，再将第2列元素与第0列元素互换，整个行互换过程为：

所述矩阵处理阵列包括1024个乘法基本处理单元、一个累加单元和一个比较器，如图12所示；所述乘法单元进行两个矩阵的乘法运算，即从内存中分别读取矩阵A的某一行与矩阵B的某一列，将行列元素分别相乘生成1024个乘积值；所述加法器对1024个乘积值做累加，完成后生成向量z，向量z与向量y通过比较器比较生成向量w，y为常量，若z>y,w＝(1-p)；若z<y,w＝p，其中p为权重系数。

本发明能够实现对光谱信号的去噪滤波、消除荧光背景对光谱信号分析的干扰、消除噪声和数据量纲对光谱信号分析的影响提高光谱信号的数据信噪比；通过设计专门的门电路替代预处理系统中的非对称最小二乘基线校正算法的软件处理方法，提高了拉曼光谱信号预处理系统的处理能力和预处理效果，在同样的处理时间条件下，通过更多次的迭代处理可以使拉曼光谱预处理算法的效果趋向于收敛到更佳的效果。

以上所述仅为本发明的优选实施例，并不用于限制本发明，显然，本领域的技术人员可以对本发明进行各种改动和变型而不脱离本发明的精神和范围。这样，倘若本发明的这些修改和变型属于本发明权利要求及其等同技术的范围之内，则本发明也意图包含这些改动和变型在内。

Claims (4)

1.拉曼光谱预处理方法，包括小波去噪阶段、非对称最小二乘法去噪阶段、标准归一化处理阶段和三次平滑样条拟合处理阶段，其特征在于：所述非对称最小二乘法去噪阶段具体包括以下步骤：

S1求解拉曼信号强度向量维度m；

S2生成m维的稀疏矩阵，并对该稀疏矩阵求2阶导生成矩阵D；

S3生成初始化m行1列的单位矩阵w；

S4以w向量为对角元素，生成维度为m行m列的带状稀疏矩阵W；

S5对(W+D’*D)做Cholesky分解，得到维度为m行m列矩阵C；

S6求C’*C矩阵的逆矩阵E；

S7向量w和向量y做点积生成向量R；

S8做矩阵乘法E*R生成1行m列向量z；

S9对向量z和采集到的拉曼信号y按位做’>’作生成新的w向量；

S10更新w向量，重复步骤S4～S10，通过多次迭代找到最优的光谱基线。

2.拉曼光谱预处理系统，包括小波去噪模块、非对称最小二乘法去噪模块、标准归一化处理模块和三次平滑样条拟合处理模块，其特征在于：所述非对称最小二乘法去噪模块包括中央控制模块、矩阵处理阵列、Cholesky求解模块和片上内存；每一次迭代过程，中央控制模块通过CHOL_RUN信号触发Cholesky求解模块执行，并且将对片上内存的读写访问总线交给Cholesky求解模块；当Cholesky求解模块执行完毕后通过CHOL_RDY信号反馈给中央控制模块；中央控制模块再将MATRIX_PROC_RUN信号拉高，使能矩阵处理阵列；矩阵处理阵列执行完毕后通过拉高MATRIX_PROC_RDY信号反馈给中央处理单元；中央处理单元执行下一次的基线校正迭代；通过多次迭代找到最优的光谱基线。

3.根据权利要求2所述的拉曼光谱预处理系统，其特征在于：所述Cholesky求解模块用于对待求解矩阵进行最大值遍历搜索和C_ij＝C_ij-a_ik*C_kj变换，i,j表示被搜索矩阵行列坐标，其中a_ik为k列元素，C_kj为k行元素，k＝1,2,3,4,······,m；并将生成的值更新到相对应坐标位置元素，将相对应最大值元素的行列坐标位置信息存放入FILO的存储器中；所述Cholesky求解模块再从FILO存储器中读取坐标信息，并依据此坐标信息进行矩阵的行列互换，并将生成的结果送到片上内存中。

4.根据权利要求2所述的拉曼光谱预处理系统，其特征在于：所述矩阵处理阵列包括乘法单元、加法器和比较器；所述乘法单元进行两个矩阵的乘法运算，即受控制信号操作从内存中分别读取矩阵A的某一行与矩阵B的某一列，将行列元素分别点乘生成1024个乘积值；所述加法器对1024个乘积值做累加，完成后生成向量z，向量z与向量y通过比较器比较生成向量w，y为采集的拉曼光谱信号常量。