CN106373157B - 高光谱异常探测方法和装置 - Google Patents

Abstract

本发明涉及一种高光谱异常探测方法和装置，该方法包括：逐行获取目标地物的高光谱数据；计算确定背景信息所需的高光谱数据的最小行数k；若获取的当前行高光谱数据的行标号大于k，则将计算当前行高光谱数据的自相关矩阵，将该自相关矩阵作为第一自相关矩阵，并将当前行之前的所有行高光谱数据的背景自相关矩阵作为第二自相关矩阵；计算当前行高光谱数据和在当前行高光谱数据之前获取的所有行高光谱数据的背景自相关矩阵；对计算得到的背景自相关矩阵进行压制处理，得到压制后的矩阵；利用压制后的矩阵对当前行高光谱数据进行异常探测。本发明逐行进行异常探测，提高处理效率，并可以满足特殊行业对高光谱数据实时或准实时处理的需求。

Description

高光谱异常探测方法和装置

技术领域

本发明涉及遥感技术领域，尤其是涉及一种高光谱异常探测方法和装置。

背景技术

目前，高光谱图像数据具有图谱合一的特征，即可以同时观测目标的空间信息和光谱信息，为目标分类和识别提供了有效的遥感技术手段。

但是，高光谱图像波段众多，而且在对异常探测时一般是在所有高光谱图像数据接收完毕后对整幅图像进行一次性全局处理，因此造成数据处理过程复杂并且效率较低，难以满足特殊行业对高光谱数据实时或准实时处理的需求。

发明内容

针对以上缺陷，本发明提供一种高光谱异常探测方法和装置，可以提高异常探测的效率，也可以满足特殊行业对高光谱数据实时或准实时处理的需求。

第一方面，本发明提供的高光谱异常探测方法包括：逐行获取目标地物的高光谱数据；

根据高光谱数据的行宽度和波段数目，计算确定背景信息所需的高光谱数据的最小行数k；

若获取的当前行高光谱数据的行标号大于k，则将计算当前行高光谱数据的自相关矩阵，将该自相关矩阵作为第一自相关矩阵，并将当前行之前的所有行高光谱数据的背景自相关矩阵作为第二自相关矩阵；

根据所述第一自相关矩阵和所述第二自相关矩阵，计算所述当前行高光谱数据和在所述当前行高光谱数据之前获取的所有行高光谱数据的背景自相关矩阵；

根据矩阵求逆算法，对计算得到的背景自相关矩阵进行压制处理，得到压制后的矩阵；

利用压制后的矩阵对当前行高光谱数据进行异常探测，并将得到的探测结果输出。

可选的，在将得到的探测结果输出之前，所述方法还包括：对所述探测结果进行非线性拉伸处理，以增强异常信号的显示亮度。

可选的，在计算所述当前行高光谱数据和在所述当前行高光谱数据之前获取的所有行高光谱数据的背景自相关矩阵之前，所述方法还包括：

判断所述当前行高光谱数据的行标号是否小于预先设定的用于确定背景信息所需的高光谱数据的最大行数，若是，则采用下式计算所述当前行高光谱数据和在所述当前行高光谱数据之前获取的所有行高光谱数据的背景自相关矩阵R_n：

R_n＝(1-1/n)R_n-1+(1/n)R(y_n)

其中，R_n-1为第二自相关矩阵，R(y_n)为第一自相关矩阵，n为当前行高光谱数据的行标号，y_n为当前行光谱数据。

可选的，所述方法还包括：若所述行标号大于或等于所述最大行数，则采用下式计算所述当前行高光谱数据和在所述当前行高光谱数据之前获取的所有行高光谱数据的背景自相关矩阵R_n：

R_n＝R_n-1+(1/m)[R(y_n)-R(y_n-m)]

其中，m为所述最大行数，R(y_n-m)为第n-m行的高光谱数据的自相关矩阵。

可选的，采用下式计算第k行高光谱数据至第1行高光谱数据的背景自相关矩阵R_k：

其中，y_i为第i行高光谱数据，L为每一行高光谱数据的宽度，i为[0，k]范围内的正整数。

可选的，所述根据矩阵求逆算法，对计算得到的背景自相关矩阵进行压制处理，得到压制后的矩阵，包括：

将该背景自相关矩阵分解为上三角矩阵和下三角矩阵；

根据该背景自相关矩阵的上三角矩阵，确定该背景自相关矩阵的逆矩阵的上三角矩阵；

根据矩阵对称性，由所述逆矩阵的上三角矩阵推导出所述逆矩阵的下三角矩阵；

根据所述逆矩阵的上三角矩阵和下三角矩阵，确定所述逆矩阵；

其中，所述确定该背景自相关矩阵的逆矩阵的上三角矩阵，包括：

构建B×B的对角矩阵，该对角矩阵中的第i'行第j'列的元素为其中，B为高光谱数据的波段数目，c_i'i'为该背景自相关矩阵的上三角矩阵中第i'行第i'列的元素，i'为矩阵元素的行标号，j'为矩阵元素的列标号；

采用下式计算所述逆矩阵的上三角矩阵中的第i'行第i'列的元素y_i'j'：

Cy_j’＝p_j’j'＝B,(B-1),L,1

其中，C为该背景自相关矩阵的上三角矩阵，y_j’为所述逆矩阵的上三角矩阵中的第j'列中各元素形成的列向量，p_j’为所述对角矩阵中第j'列中各元素形成的列向量，k为所述最小行数，c_i'k为上三角矩阵C中第i'行第k列的元素。

可选的，采用下式对当前行高光谱数据进行异常探测处理：

其中，r_i”为当前行高光谱数据的第i”个像元，为压制后的矩阵，δ(n,i”)为当前行高光谱数据中的第i”个像元的探测结果。

可选的，采用下式对所述探测结果进行非线性拉伸处理：

或δ'(n,i”)＝20logδ(n,i”)

其中，δ(n,i”)为非线性拉伸之前的探测结果，δ'(n,i”)为非线性拉伸之后的探测结果。

第二方面，本发明提供的高光谱异常探测装置，包括：

逐行扫描模块，用于逐行获取目标地物的高光谱数据；

第一计算模块，用于根据高光谱数据的行宽度和波段数目，计算确定背景信息所需的高光谱数据的最小行数k；

第二计算模块，用于若获取的当前行高光谱数据的行标号大于k，则将计算当前行高光谱数据的自相关矩阵，将该自相关矩阵作为第一自相关矩阵，并将当前行之前的所有行高光谱数据的背景自相关矩阵作为第二自相关矩阵；

第三计算模块，用于根据所述第一自相关矩阵和所述第二自相关矩阵，计算所述当前行高光谱数据和在所述当前行高光谱数据之前获取的所有行高光谱数据的背景自相关矩阵；

压制处理模块，用于根据矩阵求逆算法，对计算得到的背景自相关矩阵进行压制处理，得到压制后的矩阵；

异常探测模块，用于利用压制后的矩阵对当前行高光谱数据进行异常探测，并将得到的探测结果输出。

本发明提供的高光谱异常探测方法和装置，首先，逐行扫描高光谱数据，针对获取的当前行高光谱数据，计算第一自相关矩阵和第二自相关矩阵，再根据这两个矩阵，计算背景自相关矩阵即背景自相关矩阵。然后，在对该矩阵进行一定的压制，便于突出异常信息。最后，利用压制后的矩阵对当前行高光谱数据进行异常探测，从而完成逐行进行异常探测的目的。可见，本发明相对于传统的在对所有高光谱数据接收完毕后对整幅图像进行一次性全局处理，具有提高处理效率的优点，而且可以满足特殊行业对高光谱数据实时或准实时处理的需求。

附图说明

为了更清楚地说明本公开实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本公开的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些图获得其他的附图。

图1示出了本发明一实施例中高光谱异常探测方法的流程示意图；

图2示出了采用本发明提供的高光谱异常探测方法进行逐行异常探测时的行数与处理时间的曲线图；

图3示出了本发明一实施例中高光谱异常探测装置的结构框图。

具体实施方式

下面将结合本公开实施例中的附图，对本公开实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本公开中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本公开保护的范围。

本发明提供一种高光谱异常探测方法，如图1所示，该方法包括：

S101、逐行获取目标地物的高光谱数据；

在具体实施时，可以采用行扫描成像光谱仪进行逐行扫描，所谓的行扫描成像光谱仪可以是逐行成像的地面成像光谱仪、机载成像光谱仪或星载成像光谱仪。而且还可以采用成像光谱采集软件对行扫描成像光谱仪进行控制，使其逐行扫描，并将扫描数据进行实时显示。

S102、根据高光谱数据的行宽度和波段数目，计算确定背景信息所需的高光谱数据的最小行数k；

S103、若获取的当前行高光谱数据的行标号大于k，则将计算当前行高光谱数据的自相关矩阵，将该自相关矩阵作为第一自相关矩阵，并将当前行之前的所有行高光谱数据的背景自相关矩阵作为第二自相关矩阵；

S104、根据所述第一自相关矩阵和所述第二自相关矩阵，计算所述当前行高光谱数据和在所述当前行高光谱数据之前获取的所有行高光谱数据的背景自相关矩阵；

S105、根据矩阵求逆算法，对计算得到的背景自相关矩阵进行压制处理，得到压制后的矩阵；

可理解的是，矩阵求逆算法是指基于背景统计信息的结构即背景数据的自相关矩阵或协方差矩阵为正定对称矩阵，进行压制处理。

S106、利用压制后的矩阵对当前行高光谱数据进行异常探测，并将得到的探测结果输出。

本发明提供的高光谱异常探测方法，首先，逐行扫描高光谱数据，针对获取的当前行高光谱数据，例如假设最小行数为20，当前获取的是第21行高光谱数据，求第21行高光谱数据的自相关矩阵即第一自相关矩阵，并将第1行至第20行的高光谱数据的背景自相关矩阵作为第二自相关矩阵，再根据这两个矩阵，计算第21行至第1行的高光谱数据的背景自相关矩阵。然后，在对该矩阵进行一定的压制，便于突出异常信息。最后，利用压制后的矩阵对第21行高光谱数据进行异常探测，从而完成逐行进行异常探测的目的。当获取到第22行高光谱数据时，同样计算第22行高光谱数据的自相关矩阵作为第一自相关矩阵，并将上次计算得到的第21行至第1行的高光谱数据的背景自相关矩阵作为第二自相关矩阵，并根据这两个矩阵计算第22行至第1行的高光谱数据的背景自相关矩阵，依次类推。本发明提供的异常探测方法，相对于传统的在对所有高光谱数据接收完毕后对整幅图像进行一次性全局处理，具有提高处理效率的优点，而且可以满足特殊行业对高光谱数据实时或准实时处理的需求。参见图2，可以看出每一行探测过程所需的时间基本保持恒定，不随着扫描行数的增加而增加。

在具体实施时，当前行高光谱数据的行标号的不同，步骤S104中计算背景自相关矩阵的方式可能不同。其中，在计算背景自相关矩阵之前，对当前行高光谱数据的行标号是否小于预先设定的用于计算背景自相关矩阵所需的高光谱数据的最大行数进行判断：

若判断得知当前行高光谱数据的行标号小于上述最大行数，则可以采用下式计算背景自相关矩阵R_n：

R_n＝(1-1/n)R_n-1+(1/n)R(y_n)(1)

其中，R_n-1为第二自相关矩阵，R(y_n)为第一自相关矩阵，n为当前行高光谱数据的行标号，y_n为当前行光谱数据，该光谱数据为B×L的矩阵，B为波段数目，L为每一行高光谱数据的宽度。

若经判断得知当前行高光谱数据的行标号大于或等于上述最大行数，则可以采用下式计算背景自相关矩阵：

R_n＝R_n-1+(1/m)[R(y_n)-R(y_n-m)](2)

其中，R_n为背景自相关矩阵，m为最大行数，R(y_n-m)为第n-m行的高光谱数据的自相关矩阵。

举例来说，设定最大行数为50行，若当前行高光谱数据的行标号为30，则采用式(1)计算背景自相关矩阵；若当前行高光谱数据的行标号为80，则采用式(2)计算背景自相关矩阵。也即是说，在当前行标号小于最大行数时，每增加一行，则背景自相关矩阵增加一行，而当前行标号大于或等于最大行数时，使背景自相关矩阵的行数始终为最大行数，不再增加，便于后续的计算过程。

在具体实施时，可以采用下式计算第一自相关矩阵：

其中，R(y_n)为第一自相关矩阵，L为每一行高光谱数据的宽度，y_n为当前行高光谱数据。

在具体实施时，可以采用下式计算第k行高光谱数据至第1行高光谱数据的背景自相关矩阵R_k：

其中，y_i为第i行高光谱数据，L为每一行高光谱数据的宽度，i为[0，k]范围内的正整数。

其中，最小行数可以采用下式计算：

k＝B/L+1(5)

式中，k为最小行数，B为波段数目，L为每一行高光谱数据的宽度。

举例来说，假设确定背景信息所需的高光谱数据的最小行数为20，扫描的第1-19行中的任一行高光谱数据是无法进行异常探测的，但是当扫描的是第20行高光谱数据时，可以进行第一自相关矩阵和第二自相关矩阵的计算，进而实现对第20行高光谱数据的异常探测。对于当前行的行标号大于20的情况，当然也可以实现对其高光谱数据的异常探测。

在具体实施时，S105中可以采用多种方式对背景自相关矩阵进行压制处理，其中一种方式为求其逆矩阵，具体过程可以包括：

S1051、将背景自相关矩阵分解为上三角矩阵和下三角矩阵；

S1052、根据该背景自相关矩阵的上三角矩阵，确定该背景自相关矩阵的逆矩阵的上三角矩阵；

该S1052的具体过程可以包括：

S1052a、构建B×B的对角矩阵，该对角矩阵中的第i'行第j'列的元素为其中，B为高光谱数据的波段数目，c_i'i'为背景自相关矩阵的上三角矩阵中第i'行第i'列的元素，i'为矩阵元素的行标号，j'为矩阵元素的列标号；

可理解的是，背景自相关矩阵与其上三角矩阵的关系为：R_n＝C^TC，令则C^TCY＝I，其中，I为B×B的单位矩阵。

S1052b、采用下式计算逆矩阵的上三角矩阵中的第i'行第j'列的元素y_i'j'：

Cy_j’＝p_j’ j'＝B,(B-1),L,1(6)

其中，C为背景自相关矩阵的上三角矩阵，y_j’为逆矩阵的上三角矩阵中的第j'列中各元素形成的列向量，p_j’为对角矩阵中第j'列中各元素形成的列向量，k为最小行数，c_i'k为上三角矩阵C中第i'行第k列的元素。

根据上式(6)和(7)便可以得到逆矩阵的上三角矩阵中的各个元素，即得到逆矩阵的上三角矩阵。

S1053、根据矩阵对称性，由逆矩阵的上三角矩阵推导出逆矩阵的下三角矩阵；

S1054、根据逆矩阵的上三角矩阵和下三角矩阵，确定逆矩阵。

根据以上S1051-1054，便可得到背景自相关矩阵的逆矩阵，该逆矩阵为其压制后的矩阵。

在具体实施时，可以采用下式对当前行高光谱数据进行异常探测处理：

其中，r_i”为当前行高光谱数据的第i”个像元，为压制后的矩阵，δ(n,i”)为当前行高光谱数据中的第i”个像元的探测结果。

在具体实施时，为了进一步突出异常信号，还可以在将S106中将得到的探测结果输出之前，对探测结果进行非线性拉伸处理。进行非线性拉伸处理之后，背景信号以较暗的灰度显示，而异常信号则高亮显示，从而增强异常信号的显示亮度，利于视觉观察。

其中，进行非线性拉伸处理的方式有多种，其中一种非线性拉伸处理公式为：

当然，另一种非线性拉伸处理公式为：

δ'(n,i”)＝20logδ(n,i”)(10)

在式(9)和(10)中，δ(n,i”)为非线性拉伸之前的探测结果，δ'(n,i”)为非线性拉伸之后的探测结果。

基于相同的发明构思，本发明还提供一种高光谱异常探测装置，如图2所示，该装置200包括：

逐行扫描模块201，用于逐行获取目标地物的高光谱数据；

第一计算模块202，用于根据高光谱数据的行宽度和波段数目，计算确定背景信息所需的高光谱数据的最小行数k；

第二计算模块203，用于若获取的当前行高光谱数据的行标号大于k，则将计算当前行高光谱数据的自相关矩阵，将该自相关矩阵作为第一自相关矩阵，并将当前行之前的所有行高光谱数据的背景自相关矩阵作为第二自相关矩阵；

第三计算模块204，用于根据所述第一自相关矩阵和所述第二自相关矩阵，计算所述当前行高光谱数据和在所述当前行高光谱数据之前获取的所有行高光谱数据的背景自相关矩阵；

压制处理模块205，用于根据矩阵求逆算法，对计算得到的背景自相关矩阵进行压制处理，得到压制后的矩阵；

异常探测模块206，用于利用压制后的矩阵对当前行高光谱数据进行异常探测，并将得到的探测结果输出。

本发明提供的高光谱异常探测装置为本发明提供的高光谱异常探测方法的功能架构模块，其有关内容的解释、说明、举例、优选方式和有益效果等内容请参考本发明提供的高光谱异常探测方法的相应部分，这里不再赘述。

本领域普通技术人员可以理解：实现上述方法实施例的全部或者部分步骤可以通过程序指令相关的硬件来完成，前述的程序可以存储在计算机可读取的存储介质中，该程序在执行时，执行包括上述方法实施例的步骤。

本发明的说明书中，说明了大量具体细节。然而，能够理解，本发明的实施例可以在没有这些具体细节的情况下实践。在一些实例中，并未详细示出公知的方法、结构和技术，以便不模糊对本说明书的理解。

以上实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解；其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明各实施例技术方案的精神和范围。

Claims (8)

1.一种高光谱异常探测方法，其特征在于，包括：

逐行获取目标地物的高光谱数据；

根据高光谱数据的行宽度和波段数目，计算确定背景信息所需的高光谱数据的最小行数k；

若获取的当前行高光谱数据的行标号大于k，则将计算当前行高光谱数据的自相关矩阵，将该自相关矩阵作为第一自相关矩阵，并将当前行之前的所有行高光谱数据的背景自相关矩阵作为第二自相关矩阵；

根据所述第一自相关矩阵和所述第二自相关矩阵，计算所述当前行高光谱数据和在所述当前行高光谱数据之前获取的所有行高光谱数据的背景自相关矩阵；

根据矩阵求逆算法，对计算得到的背景自相关矩阵进行压制处理，得到压制后的矩阵；

利用压制后的矩阵对当前行高光谱数据进行异常探测，并将得到的探测结果输出；

其中，在计算所述当前行高光谱数据和在所述当前行高光谱数据之前获取的所有行高光谱数据的背景自相关矩阵之前，还包括：判断所述当前行高光谱数据的行标号是否小于预先设定的用于确定背景信息所需的高光谱数据的最大行数；

若所述行标号大于或等于所述最大行数，则采用下式计算所述当前行高光谱数据和在所述当前行高光谱数据之前获取的所有行高光谱数据的背景自相关矩阵R_n：

R_n＝R_n-1+(1/m)[R(y_n)-R(y_n-m)]

其中，m为所述最大行数，R(y_n-m)为第n-m行的高光谱数据的自相关矩阵，R_n-1为第二自相关矩阵，R(y_n)为第一自相关矩阵，n为当前行高光谱数据的行标号，y_n为当前行光谱数据。

2.根据权利要求1所述的方法，其特征在于，在将得到的探测结果输出之前，所述方法还包括：对所述探测结果进行非线性拉伸处理，以增强异常信号的显示亮度。

3.根据权利要求1所述的方法，其特征在于，在计算所述当前行高光谱数据和在所述当前行高光谱数据之前获取的所有行高光谱数据的背景自相关矩阵之前，所述方法还包括：

若所述当前行高光谱数据的行标号小于预先设定的用于确定背景信息所需的高光谱数据的最大行数，则采用下式计算所述当前行高光谱数据和在所述当前行高光谱数据之前获取的所有行高光谱数据的背景自相关矩阵R_n：

R_n＝(1-1/n)R_n-1+(1/n)R(y_n)

其中，R_n-1为第二自相关矩阵，R(y_n)为第一自相关矩阵，n为当前行高光谱数据的行标号，y_n为当前行光谱数据。

4.根据权利要求1所述的方法，其特征在于，采用下式计算第k行高光谱数据至第1行高光谱数据的背景自相关矩阵R_k：

其中，y_i为第i行高光谱数据，L为每一行高光谱数据的宽度，i为[0，k]范围内的正整数。

5.根据权利要求1所述的方法，其特征在于，所述根据矩阵求逆算法，对计算得到的背景自相关矩阵进行压制处理，得到压制后的矩阵，包括：

将该背景自相关矩阵分解为上三角矩阵和下三角矩阵；

根据该背景自相关矩阵的上三角矩阵，确定该背景自相关矩阵的逆矩阵的上三角矩阵；

根据矩阵对称性，由所述逆矩阵的上三角矩阵推导出所述逆矩阵的下三角矩阵；

根据所述逆矩阵的上三角矩阵和下三角矩阵，确定所述逆矩阵；

其中，所述确定该背景自相关矩阵的逆矩阵的上三角矩阵，包括：

构建B×B的对角矩阵，该对角矩阵中的第i'行第j'列的元素为其中，B为高光谱数据的波段数目，c_i'i'为该背景自相关矩阵的上三角矩阵中第i'行第i'列的元素，i'为矩阵元素的行标号，j'为矩阵元素的列标号；

采用下式计算所述逆矩阵的上三角矩阵中的第i'行第j'列的元素y_i'j'：

Cy_j’＝p_j’j'＝B,(B-1),L,1

其中，C为该背景自相关矩阵的上三角矩阵，y_j’为所述逆矩阵的上三角矩阵中的第j'列中各元素形成的列向量，p_j’为所述对角矩阵中第j'列中各元素形成的列向量，k为所述最小行数，c_i'k为上三角矩阵C中第i'行第k列的元素。

6.根据权利要求1所述的方法，其特征在于，采用下式对当前行高光谱数据进行异常探测处理：

其中，r_i”为当前行高光谱数据的第i”个像元，为压制后的矩阵，δ(n,i”)为当前行高光谱数据中的第i”个像元的探测结果。

7.根据权利要求2所述的方法，其特征在于，采用下式对所述探测结果进行非线性拉伸处理：

或δ'(n,i”)＝20logδ(n,i”)

其中，δ(n,i”)为非线性拉伸之前的探测结果，δ'(n,i”)为非线性拉伸之后的探测结果。

8.一种高光谱异常探测装置，其特征在于，包括：

逐行扫描模块，用于逐行获取目标地物的高光谱数据；

第一计算模块，用于根据高光谱数据的行宽度和波段数目，计算确定背景信息所需的高光谱数据的最小行数k；

第二计算模块，用于在获取的当前行高光谱数据的行标号大于k时，将计算当前行高光谱数据的自相关矩阵，将该自相关矩阵作为第一自相关矩阵，并将当前行之前的所有行高光谱数据的背景自相关矩阵作为第二自相关矩阵；

第三计算模块，用于根据所述第一自相关矩阵和所述第二自相关矩阵，计算所述当前行高光谱数据和在所述当前行高光谱数据之前获取的所有行高光谱数据的背景自相关矩阵；

压制处理模块，用于根据矩阵求逆算法，对计算得到的背景自相关矩阵进行压制处理，得到压制后的矩阵；

异常探测模块，用于利用压制后的矩阵对当前行高光谱数据进行异常探测，并将得到的探测结果输出；

其中，所述第二计算模块还用于：判断所述当前行高光谱数据的行标号是否小于预先设定的用于确定背景信息所需的高光谱数据的最大行数；

若所述行标号大于或等于所述最大行数，则采用下式计算所述当前行高光谱数据和在所述当前行高光谱数据之前获取的所有行高光谱数据的背景自相关矩阵R_n：

R_n＝R_n-1+(1/m)[R(y_n)-R(y_n-m)]

其中，m为所述最大行数，R(y_n-m)为第n-m行的高光谱数据的自相关矩阵，R_n-1为第二自相关矩阵，R(y_n)为第一自相关矩阵，n为当前行高光谱数据的行标号，y_n为当前行光谱数据。