clear;
clc; 

%確定局部區域大小: (2*offset+1)*(2*offset+1)
offset=3;
[x,map]=imread('img_in.bmp');

x=double(x);
[ny,nx]=size(x);
y=zeros(ny,nx);

if (offset==1)
    %Epanechnikov kernel 3*3
    wcoef=[0.1899,0.1671,0.0962;0.1671,0.1418,0.0709;0.0962,0.0709,0];
else
    %Epanechnikov kernel 7*7
    wcoef=[0.0322,0.0318,0.0305,0.0284,0.0249,0.0210,0.0163;
           0.0318,0.0314,0.0301,0.0279,0.0245,0.0206,0.0159;
           0.0305,0.0301,0.0288,0.0266,0.0232,0.0193,0.0142;
           0.0284,0.0279,0.0266,0.0241,0.0210,0.0168,0.0103;
           0.0249,0.0245,0.0232,0.0210,0.0180,0.0137,0.0073;
           0.0210,0.0206,0.0193,0.0168,0.0137,0.0099,0.0034;
           0.0163,0.0159,0.0142,0.0103,0.0073,0.0034,0];
end

%用Epanechnikov Kernel 進行背景抑制
for i=1:(ny-2*offset)
 for j=1:(nx-2*offset)		                                
	 arr=x(i:(2*offset+1)+i-1,j:(2*offset+1)+j-1);
     %在局部塊中進行列(或行)均值去除
     for p=1:(2*offset+1)
         arr(:,p)=arr(:,p)-mean(arr(:,p));
     end
     weightedsum=sum(sum(arr.*wcoef));
     yn(i,j)= arr(1,1)-weightedsum;
     yb(i,j)= weightedsum;
     if (j==nx-2*offset)
         for k=1:2*offset
            yn(i,j+k)= arr(1,k)-weightedsum; 
            yb(i,j+k)= weightedsum; 
         end
     end
      if (i==ny-2*offset)
         for k=1:2*offset      
             yn(i+k,j)= arr(k,1)-weightedsum; 
             yb(i+k,j)= weightedsum; 
         end
       end
     if (j==nx-2*offset)&(i==ny-2*offset)  
          for k=1:2*offset
              for l=1:2*offset
                  yn(i+k,j+l)= arr(k,l)-weightedsum; 
                  yb(i+k,j+l)= weightedsum;
              end
         end
     end
 end
end
%從殘留圖像中選取大小為100*100的像素塊
 yn=yn(ny/2:(ny/2+99),nx/2:(nx/2+99));

%計算Friedman統計量及概率alpha
[p,table,stats]=friedman(yn(1:100,1:100),1,'off');
s(1)=table(2,5);
a(1)=table(2,6);
[p,table,stats]=friedman(yn(1:50,1:100),1,'off');
s(2)=table(2,5);
a(2)=table(2,6);
[p,table,stats]=friedman(yn(51:100,1:100),1,'off');
s(3)=table(2,5);
a(3)=table(2,6);
[p,table,stats]=friedman(yn(26:75,1:100),1,'off');
s(4)=table(2,5);
a(4)=table(2,6);
[p,table,stats]=friedman(yn(1:30,1:100),1,'off');
s(5)=table(2,5);
a(5)=table(2,6);
[p,table,stats]=friedman(yn(31:60,1:100),1,'off');
s(6)=table(2,5);
a(6)=table(2,6);
[p,table,stats]=friedman(yn(61:90,1:100),1,'off');
s(7)=table(2,5);
a(7)=table(2,6);
[p,table,stats]=friedman(yn(1:20,1:100),1,'off');
s(8)=table(2,5);
a(8)=table(2,6);
[p,table,stats]=friedman(yn(21:60,1:100),1,'off');
s(9)=table(2,5);
a(9)=table(2,6);
[p,table,stats]=friedman(yn(71:90,1:100),1,'off');
s(10)=table(2,5);
a(10)=table(2,6);

%對高斯白噪聲計算Friedman統計量及概率alpha,與上述結果進行比較
yn=randn(100,100);
[p,table,stats]=friedman(yn(1:100,1:100),1,'off');
sn(1)=table(2,5);
an(1)=table(2,6);
[p,table,stats]=friedman(yn(1:50,1:100),1,'off');
sn(2)=table(2,5);
an(2)=table(2,6);
[p,table,stats]=friedman(yn(51:100,1:100),1,'off');
sn(3)=table(2,5);
an(3)=table(2,6);
[p,table,stats]=friedman(yn(26:75,1:100),1,'off');
sn(4)=table(2,5);
an(4)=table(2,6);
[p,table,stats]=friedman(yn(1:30,1:100),1,'off');
sn(5)=table(2,5);
an(5)=table(2,6);
[p,table,stats]=friedman(yn(31:60,1:100),1,'off');
sn(6)=table(2,5);
an(6)=table(2,6);
[p,table,stats]=friedman(yn(61:90,1:100),1,'off');
sn(7)=table(2,5);
an(7)=table(2,6);
[p,table,stats]=friedman(yn(1:20,1:100),1,'off');
sn(8)=table(2,5);
an(8)=table(2,6);
[p,table,stats]=friedman(yn(21:60,1:100),1,'off');
sn(9)=table(2,5);
an(9)=table(2,6);
[p,table,stats]=friedman(yn(71:90,1:100),1,'off');
sn(10)=table(2,5);
an(10)=table(2,6);