高光谱异常检测KRX算法

在MATLAB中实现高光谱异常检测KRX算法，主要是参靠KRX英文原文《Kernel RX-Algorithm: A Nonlinear Anomaly Detector for Hyperspectral Imagery》实现。

function result=KRX（imagewinlenthtargetShape）
%image 输入3维高光谱图像
%winlenth 局部窗口大小
%targetShape  目标形状
%result 检测结果
switch targetShape%目标形状
case 0
target_spatial=...
    [0 0 1 0 0;
    0 0 1 0 0;
    1 1 1 1 1;
    0 0 1 0 0;
    0 0 1 0 0];

case 1
target_spatial=...
    [1 1 1;
    1 1 1;
    1 1 1];
case 2 
target_spatial=...
    [1 0 0 0 1; 
    0 1 0 1 0;
    0 0 1 0 0;
    0 1 0 1 0;
    1 0 0 0 1];
case 3
target_spatial=...
    [1 0 0 0 0 0 1;
    0 1 0 0 0 1 0;
    0 0 1 0 1 0 0;
    0 0 0 1 0 0 0;
    0 0 1 0 1 0 0;
    0 1 0 0 0 1 0;
    1 0 0 0 0 0 1];
otherwise
disp（‘Unknown shape.‘）
end
sizeTarget=size（target_spatial）;
% 若为二维图象，将其看作一幅第三维为1的三维图象
SizeA=size（image）;
if length（SizeA）==2 
SizeA（3） = 1;
end
% 若为二维图象，将其看作一幅第三维为1的三维图象
SizeA=size（image）;
if length（SizeA）==2 
SizeA（3） = 1;
end
%将image扩充为imageExpand
imageExpand（（winlenth+1）/2:（winlenth-1）/2+SizeA（1） （winlenth+1）/2:（winlenth-1）/2+SizeA（2）:）=image; %中
imageExpand（（winlenth-1）/2:-1:1（winlenth+1）/2:（winlenth-1）/2+SizeA（2）:）=image（1:（winlenth-1）/2::）; %上
imageExpand（winlenth-1+SizeA（1）:-1:（winlenth+1）/2+SizeA（1）（winlenth+1）/2:（winlenth-1）/2+SizeA（2）:）...%下
=image（SizeA（1）-（winlenth-3）/2:SizeA（1）::）; %
imageExpand（（winlenth+1）/2:（winlenth-1）/2+SizeA（1）（winlenth-1）/2:-1:1:）=image（:1:（winlenth-1）/2:）; %左
imageExpand（（winlenth+1）/2:（winlenth-1）/2+SizeA（1）winlenth-1+SizeA（2）:-1:（winlenth+1）/2+SizeA（2）:）...%右
=image（:SizeA（2）-（winlenth-3）/2:SizeA（2）:）; %
imageExpand（（winlenth-1）/2:-1:1（winlenth-1）/2:-1:1:）=image（1:（winlenth-1）/21:（winlenth-1）/2:）; %左上
imageExpand（winlenth-1+SizeA（1）:-1:（winlenth+1）/2+SizeA（1）winlenth-1+SizeA（2）:-1:（winlenth+1）/2+... %右下
SizeA（2）:）=image（SizeA（1）-（winlenth-3）/2:SizeA（1）SizeA（2）-（winlenth-3）/2:SizeA（2）:）; %
imageExpand（winlenth-1+SizeA（1）:-1:（winlenth+1）/2+SizeA（1）（winlenth-1）/2:-1:1:）... %左下
=image（SizeA（1）-（winlenth-3）/2:SizeA（1）1:（winlenth-1）/2:）; %
imageExpand（（winlenth-1）/2:-1:1winlenth-1+SizeA（2）:-1:（winlenth+1）/2+SizeA（2）:）... %右上
=image（1:（winlenth-1）/2SizeA（2）-（winlenth-3）/2:SizeA（2）:）; %
clear image;

spatial_pattern=zeros（winlenthwinlenth）;%winlenth应该是窗口大小
targetSize=size（target_spatial）;
result=zeros（SizeA（1）SizeA（2））;
for i=1:targetSize（1）
    for j=1:targetSize（2）
        spatial_pattern（（winlenth-targetSize（1））/2+i（winlenth-targetSize（1））/2+j）=1;
    end
end
s=zeros（1winlenth*winlenth）;
    k=1;
for i=1:winlenth
for j=1:winlenth

% s（（i-1）*winlenth+j）=spatial_pattern（