0 引言

高光谱图像具有光谱分辨率高、图谱合一的独特优势，其每一层图像对应一个光谱波段，每个像素点对应一条连续的光谱曲线^[1-3]。高光谱图像中蕴含丰富的光谱信息，而光谱特征是不同物质的固有特性，因此利用高光谱图像进行目标检测便成为高光谱图像处理中的热点问题，并且在矿物勘探、军事侦察、气体检测、舰船检测等领域都得到了成功的应用。

舰船检测在海上交通管制、渔业监管、海战形势预测等领域应用广泛^[4]。但是，舰船检测经常受到云雾、阴影、海面风浪、码头、小型岛屿的影响，导致检测率低而虚警率高。如何有效提升检测率、降低虚警率成为舰船检测的热点与难点问题。传统舰船检测技术一般借助海上预警雷达来完成，这种方法面临的主要问题是：对于小型舰船目标，在海面背景杂波较为严重时，背景杂波的雷达回波信号强于目标回波信号，淹没了目标的反射特性，从而造成漏检及误检^[5]。基于红外图像的舰船检测方法主要依靠提取目标的形状特征和空间关系来完成舰船检测，但是，由于成像机制的问题，图像中舰船目标尺寸较小，且很难看出纹理特征，这在一定程度上限制了算法的检测精度^[6-7]。

基于高光谱图像的舰船检测方法主要利用舰船目标与周围海水背景的光谱差异较大的特性，不仅能检测出舰船目标的位置、形状、尺寸等空间信息，又能对舰船的表面材质进行初步分析，便于后续舰船类型识别工作的开展。文献[8-9]提出利用独立分量分析的方法进行舰船检测，该方法借助快速独立分量分析，将舰船目标信息从高维影像投影到低维特征影像，以峰度为度量指标对特征影像进行选择，最后用以偏度为指标的直方图分割法提取舰船目标。代威等^[10]提出一种基于光谱异常的舰船检测方法，提取包含舰船目标的海水区域，构建海水背景模型，通过检测海域中的光谱异常点及虚警剔除实现舰船目标检测。陈梦云等^[11]提出基于视觉显著性进行舰船检测，通过在传统检测算子中引入视觉显著性算法，抑制背景像素，加强目标像素，在海况较为复杂的情况下完成了舰船检测。文献[12]等通过现场测量舰船目标与海水背景的先验光谱，结合多种光谱匹配算法对高光谱图像进行舰船检测，总结了每种算法的准确性、优势及劣势。崔震等^[13]将经过最小噪声分离后的高光谱图像进行基于集合顶点的端元提取，结合多种目标检测算法对飞机目标的自动检测进行了深入研究。王文秀等^[7]通过分析红外图像中较为完整的连通区域的矩特征，将之与特征数据库进行匹配，完成了港口区域在轨船只目标的快速检测。

以上几种方法适用于海面背景较为单一、海况不复杂的远海场景，但在近海区域，尤其在海岸线附近，小岛屿、码头等与舰船尺寸相似的目标众多，基于独立分量分析的方法并不能很好地完成特征影像的选择及舰船目标的提取。阈值分割和形态学处理在分割海陆、消除舰船目标的同时，可能会将海岸附近的非舰船目标也消除掉，从而使后续检测产生大量虚警，基于视觉显著性的检测方法在处理海面岛屿的情况时也容易出现误检及漏检。基于光谱匹配的方法由于光谱库数据匮乏和“同物异谱”效应导致目标的先验光谱往往难以获取，从而使其不能广泛适用。近海场景中的小岛屿、码头等目标的光谱特性与陆地相似，与海水却有很大差别。如果能从检测场景中获取陆地和舰船目标的近似后验光谱，这样就可以利用光谱匹配算法在没有任何先验信息的情况下，去除海域中小岛屿、码头等目标，检测出舰船目标。

基于以上分析，本文提出一种基于近似后验光谱的舰船检测方法。首先对高光谱近、短波红外波段的灰度图像进行海陆预分割；然后通过顶点成分分析^[14-15](vertex component analysis，VCA)从陆地区域获取陆地近似后验光谱，将其作为已知光谱对海域进行约束能量最小化^[16](constrained energy minimization，CEM)目标检测，去除海域中的疑似陆地像元，完成海陆的再次分割；最后再利用VCA从新海域中获取舰船近似后验光谱并作为已知光谱，再次对海域进行CEM检测，提取出舰船目标。

1 基于近似后验信息的高光谱图像海上舰船检测方法

算法主要由海陆预分割、海域疑似陆地像元去除和海域舰船目标检测三部分组成。算法流程图如图 1所示。

1.1 海陆预分割

由于水体和陆地在近红外、短波红外波段反射率差异较大，因此，以近红外和短波红外波段数据为基础，通过对其进行对比度拉伸、灰度阈值分割和形态学处理初步完成海陆分割。首先将近红外、短波红外的r波段影像数据{S_i}₁^r转换为灰度均值图像G，即$\mathit{\boldsymbol{G}}\left( {x, y} \right) = \frac{1}{r}\sum\limits_{i = 1}^r {{\mathit{\boldsymbol{S}}_i}\left( {x, y} \right)} $。

然后通过对比度拉伸^[17]，扩大水体与陆地两者的类间方差，以提升阈值分割的效果。灰度图G进行对比度拉伸后得到拉伸图像G′，G′(x, y)=1/{1+[m/G(x, y)]^E}，式中：m、E为对比度拉伸变换的常参数。

最后采用最大类间方差法^[18]确定阈值，对拉伸图像进行灰度分割。由于舰船目标在近、短波红外波段的反射特性与陆地相似，会被误分为陆地，内陆中的水域也会被分割为海域。因此，采用形态学运算来消除海域中的船只、小岛屿目标和陆地区域中的孔洞，同时达到平滑海陆边界的效果。

1.2 海域疑似陆地像元去除

海陆预分割将海域和陆地初步分离，但此时的海域中仍存在与舰船尺寸相似的小型岛屿、码头、浅海滩涂等目标，如果直接在此时的海域中进行舰船检测会导致大量虚警。由于这些目标与陆地有着相似的光谱特性，而与海水的光谱差异较大，因此，利用VCA对海陆预分割后的陆地区域进行端元提取，获取陆地近似后验光谱。由于CEM算子对近似后验光谱有很好的容差性，将陆地近似后验光谱作为已知光谱在海域进行CEM光谱匹配，去除海域中匹配值较高的疑似陆地像元，完成海陆再次分割，从而降低舰船检测的虚警率。

1.2.1 陆地近似后验光谱获取

VCA无须先验知识，仅从混合光谱数据中依次提取端元，且前后提取的端元互不相关。在陆地区域中，陆地作为背景成分占据大部分信息，因此可通过VCA提取陆地区域端元光谱矩阵，作为该区域的背景光谱矩阵，并对其求均值得出陆地近似后验光谱。

已知高光谱数据S=[s₁, s₂, …, s_N]为b×N维矩阵，其中：b为波段数；N为像素数。首先利用奇异值分解法对观测光谱矩阵降维，变换到p维子空间，即X=U_p^TS，其中：U_p是由奇异值分解得到的左变换矩阵的前p个向量组成的b×p维矩阵。

将X投影到超平面上得到Y，

$ {\left[ \mathit{\boldsymbol{Y}} \right]_{:, j}} = {\left[ \mathit{\boldsymbol{X}} \right]_{:, j}}/\left( {\left[ X \right]_{:, j}^{\rm{T}}\mathit{\boldsymbol{u}}} \right), $

(1)

其中：[M]_{:, j}表示矩阵M的第j列；[M]_{:, i:k}表示矩阵M的第i至k列；u=mean(X)是p维的列向量。

初始给定一个方向f，

$ \mathit{\boldsymbol{f}} = \left( {\left( {\mathit{\boldsymbol{I}} - \mathit{\boldsymbol{A}}{\mathit{\boldsymbol{A}}^\# }} \right)\mathit{\boldsymbol{w}}} \right)/\left( {\left\| {\left( {\mathit{\boldsymbol{I}} - \mathit{\boldsymbol{A}}{\mathit{\boldsymbol{A}}^\# }} \right)\mathit{\boldsymbol{w}}} \right\|} \right), $

(2)

其中：A=[e_u|0|…|0]，A是p×p维矩阵，用来存储估计端元信号的投影；e_u=[0, …, 0, 1]^T为单位向量；A^#=(A^TA)^-1A^T为A的伪逆；w=randn(0, I_p)是均值为0、协方差矩阵为单位阵I_p的随机高斯向量；f是正交于由[A]_{:, 1:i}张成的子空间的向量。

将(1)式得到的超平面上的数据Y投影到(2)式给定的方向f上，得到

$ \mathit{\boldsymbol{v}} = {\mathit{\boldsymbol{f}}^{\rm{T}}}\mathit{\boldsymbol{Y}}, $

(3)

该投影的极值对应的像元位置k即为第一个端元的位置，

$ k = \arg {\max _{j = 1, \cdots , N}}\left\| {{{\left[ \mathit{\boldsymbol{v}} \right]}_{:, j}}} \right\|。$

(4)

将(4)式的结果存储到像元索引数组[indice]_i=k。利用(4)式的结果构造下一个端元的投影方向，

$ {\left[ \mathit{\boldsymbol{A}} \right]_{:, j}} = {\left[ \mathit{\boldsymbol{Y}} \right]_{:, k}}。$

(5)

将(5)式中更新后的A重新代入(2)式，重复执行式(3)和(4)的运算，求出下一个端元所在的位置，依次类推，求出p个端元所对应的位置。最后，用M=U_p[X]_{:, indice}求解各端元的光谱曲线，M即为对陆地区域经VCA后的背景光谱矩阵，对其取均值，得到陆地近似后验光谱Land，Land=$\frac{1}{p}\sum\limits_{i = 1}^p {\left[ \mathit{\boldsymbol{M}} \right]:, i} $。

1.2.2 海域CEM陆地目标检测

将陆地近似后验光谱Land作为已知陆地光谱，对海域像元进行CEM陆地像元光谱匹配。CEM在仅知道目标先验光谱，而背景光谱在未知的条件下进行目标检测，在提取感兴趣信号的同时抑制其他方向的信号，可以达到突出目标的效果。

将近似后验光谱Land作为目标光谱d，高光谱数据的自相关矩阵为R，设计一个线性滤波器w=[w₁, w₂, …, w_b]^T，使其在满足式(6)的条件下输出能量最小。

$ {\mathit{\boldsymbol{d}}^{\rm{T}}}\mathit{\boldsymbol{w = }}1。$

(6)

当输入向量为x_i时，滤波器输出能量为y_i，y_i=$\sum\limits_{j = 1}^b {{\mathit{\boldsymbol{w}}_j}{\mathit{\boldsymbol{x}}_{ij}}} $=w^Tx_i。

高光谱所有数据经过滤波器后的平均输出能量为

$ \mathit{\boldsymbol{E}} = \frac{1}{N}\sum\limits_{i = 1}^N {\mathit{\boldsymbol{y}}_i^2 = {\mathit{\boldsymbol{w}}^{\rm{T}}}\mathit{\boldsymbol{Rw}}} , $

(7)

这样，滤波器的设计就转化为求式(6)在式(7)条件下的最小值问题，其解即为CEM匹配算子

$ {\mathit{\boldsymbol{w}}^*} = \left( {{\mathit{\boldsymbol{R}}^{ - 1}}\mathit{\boldsymbol{d}}} \right)/\left( {{\mathit{\boldsymbol{d}}^{\rm{T}}}\mathit{\boldsymbol{Rd}}} \right)。$

像元x经过CEM检测后的输出能量为e=w^*Tx=(d^TR^-1x)/(d^TR^-1d)。

对预分割后的海域进行CEM检测，即可得出其中疑似陆地像元的分布概率图像。CEM检测的输出能量值e越高，表示该像元属于陆地目标的概率越大，反之越小。

1.2.3 疑似陆地像元判别

岛屿边缘、码头、浅海滩涂等目标在海陆预分割阶段不能被准确归类，这些目标大部分位于海岸线附近，同时，海域中部分舰船目标的CEM能量输出值可能也比较高。为达到最佳的判别效果，疑似陆地像元判别的区域只集中在海岸线附近的海域中。

选取判别海域采用滤波的方法，利用滤波器对海陆预分割获取的陆地掩模图像进行判别滤波，确定海岸线附近的海域像元，将这些像元对应的CEM检测结果与判别阈值η进行比较，大于阈值η的像元标记为疑似陆地像元。

判别过程结束后，将疑似陆地像元标记结果与海陆预分割掩模中的陆地区域进行合并，实现陆地区域的扩张，完成海陆再分割。

1.3 海域舰船目标检测

经疑似陆地像元去除后的海域中，舰船属于小概率目标，而海水占据了海域的大部分信息，属于背景成分。因此在利用VCA获取舰船近似后验光谱之前需要对海域进行海水背景抑制，然后再对海域进行CEM检测，提取出舰船目标。

1.3.1 舰船近似后验光谱获取

首先采用主成分分析的方法抑制海水背景，对海域高光谱影像进行降维，再将前几个主成分投影到海域中，即可抑制海水背景成分。对背景成分抑制后的海域利用VCA进行端元提取，提取过程与1.2.1节类似。同样，对端元提取后的端元光谱矩阵M取均值，得到舰船近似后验光谱Fleet，Fleet=$\frac{1}{p}\sum\limits_{i = 1}^p {\left[ \mathit{\boldsymbol{M}} \right]:, i} $。

1.3.2 海域CEM舰船目标检测

同1.2.2节的方法类似，将舰船近似后验光谱作为已知目标光谱，对海陆再分割后的海域进行CEM检测，将检测结果中大于设定阈值η的像元作为舰船目标，经阈值分割得到最终舰船检测结果。

2 实验结果与分析

实验采用EO-1卫星搭载的Hyperion传感器获取的高光谱影像数据^[19]。该数据为L1Gst级产品，已进行辐射校正、几何校正、投影配准和地形校正。实验影像的空间分辨率为30 m，大小为170×140，地点位于新加坡市的巴威岛和安乐岛海域，共242个波段，去除未定标和受水汽影响较重的波段，剩余163个波段^[20]。图 2(a)为实验数据的舰船标记图像，图中共有7艘舰船。

本文提出的基于近似后验信息的舰船检测方法两次利用后验信息，分别用于海域疑似陆地像元去除和舰船检测，进行了两次海陆分割，将其记为DS-CEM(double split CEM)。同时，为说明海域疑似陆地像元对舰船检测结果的影响，本文提出另一种算法：略过海域疑似陆地像元去除过程，直接从海陆预分割的海域中提取舰船近似后验光谱并进行舰船检测。该算法只进行一次海陆分割，将其记为S-CEM(split CEM)。为验证提出算法的有效性，本文用于作对比的基于传统RX^[21]异常探测的方法分别在两种场景下进行：1)只进行海陆预分割的场景，该场景的海域未去除疑似陆地像元，记为S-RX(split RX)；2)进行两次海陆分割的场景，该场景的海域中疑似陆地像元已被去除，记为DS-RX(double split RX)。

首先，将近、短波红外波段的均值灰度图像进行对比度拉伸和最大类间方差阈值分割。并利用形态学处理去除陆地区域的孔洞和海域中的舰船，完成海陆预分割。此处，应恰当选择膨胀、腐蚀算子的尺寸，在消除舰船目标的同时又能最大限度保留住岛屿目标。如果陆地区域没有孔洞，则无须进行孔洞去除。海陆预分割的结果如图 2(b)所示。

然后进行海域疑似陆地像元去除。先通过VCA对陆地区域进行端元提取，端元数目p设为10；对端元光谱取均值得到陆地近似后验光谱。再将其作为已知光谱对海域进行CEM检测。最后对海岸线附近海域像元进行判别，大于指定阈值η的像元将作为疑似陆地像元从海域中去除。经多次实验，η设为0.01时效果较好。

图 3(a)中的蓝色曲线为VCA从陆地区域提取的端元光谱，可以看出，各光谱的形状基本一致，但受光照、温度等因素的影响，各光谱的幅度不同。对各端元光谱取均值得到陆地近似后验光谱，如图 3(a)中的红色曲线所示。图 3(b)为以陆地近似后验光谱为已知光谱进行CEM检测的结果，CEM输出能量较高的像元大部分位于海岸线附近。图 3(c)为海陆再分割后的掩模图像，经海域疑似陆地像元去除后，陆地区域实现进一步扩张。

最后对再分割后的海域进行背景抑制，并通过VCA对背景抑制后的海域进行端元提取，端元数目p设为10，将端元光谱曲线求均值后得到舰船近似后验光谱，再次作为已知光谱对海域进行CEM检测，最后以经验阈值η=1.0进行分割得到舰船目标分割图像。

通过对比图 4与图 3的端元光谱曲线和近似后验光谱，可以看出，舰船目标与陆地目标在可见光波段(1~35)的光谱特性差异较大，在近、短波红外波段(36~163)光谱特性相似。

图 4(b)为以舰船近似后验光谱为已知光谱，进行CEM检测的结果，可以看出，舰船目标的输出能量明显较高。图 4(c)为CEM检测结果经阈值分割后的舰船分割结果。

本文提出的4种算法检测结果如图 5所示。从图 5(a)的检测结果可以看出，在DS-CEM算法中，7艘舰船目标全部被检测出，虚警数目为2；而在图 5(b)和图 5(d)中，S-CEM与S-RX算法由于都未进行海域疑似陆地像元去除，所以检测结果中含有大量位于海岸线附近的虚警，且S-RX算法未完成舰船检测的任务。图 5(c)所示为DS-RX算法的检测结果，由于去除海域中的疑似陆地像元，大幅减少了海岸线附近的虚警，但仍未全部检测出舰船目标。

为更好地说明算法的检测效果，采用接收机操作特性(receiver operating characteristics, ROC)曲线来对算法的检测性能进行定量分析。ROC用来描述检测概率和虚警概率之间的变化关系，曲线越靠近左上角，说明算法的检测性能越好。4种检测方法的ROC曲线如图 6所示。

从图中可以看出，DS-CEM算法的ROC曲线最陡峭。S-CEM和S-RX算法都未进行疑似陆地像元去除，导致虚警率较高，因此它们的曲线都分别在DS-CEM和DS-RX算法的曲线之下，这也说明了疑似陆地像元去除工作的重要性。DS-CEM算法相比于DS-RX算法，通过提取舰船近似后验光谱进行光谱匹配，其ROC曲线在更低的虚警率时便达到了很高的检测率，进一步证实了基于近似后验信息的舰船检测方法的可行性。

3 结论

本文提出一种基于近似后验信息的海上舰船检测算法。该算法首先进行海陆预分割，初步分离海陆区域；然后对陆地区域进行顶点成分分析，获得陆地近似后验光谱，并将其作为已知光谱对海域进行CEM检测，去除海域中的疑似陆地像元；最后对海域进行顶点成分分析，获得舰船近似后验光谱，再次将其作为已知光谱对海域进行CEM检测，提取出舰船目标。实验证明，该算法能在没有任何先验信息的情况下去除海岸线附近的疑似陆地像元，完成舰船目标的检测，有效抑制虚警率，但在进行海域疑似陆地像元去除时，如何选取合适范围的判别海域以及判别阈值，以达到最佳的疑似陆地像元去除效果，还需进一步研究。