clc;
clear?all;
img_name=‘jiedao.jpg‘;
%?原始圖像
I=double（imread（img_name））/255;
%?獲取圖像大小
[hwc]=size（I）;
win_size?=?7;
img_size=w*h;
figure?imshow（I）;
win_dark=ones（hw）;
%計算分塊darkchannel
?for?j=1+win_size:w-win_size
????for?i=win_size+1:h-win_size
????????m_pos_min?=?min（I（ij:））;
????????for?n=j-win_size:j+win_size
????????????for?m=i-win_size:i+win_size
????????????????if（win_dark（mn）>m_pos_min）
????????????????????win_dark（mn）=m_pos_min;
????????????????end
????????????end
????????end
????end
?end
?%選定精確dark?value坐標
%?win_b?=?zeros（img_size1）;
?figure?imshow（win_dark）;
?win_t=1-0.95*win_dark;
?win_b=zeros（img_size1）;
for?ci=1:h
????for?cj=1:w
????????if（rem（ci-815）<1）
????????????if（rem（cj-815）<1）
????????????????win_b（ci*w+cj）=win_t（ci*w+cj）;
????????????end
????????end
????end
end
?
%顯示分塊darkchannel
%figure?imshow（win_dark）;
neb_size?=?9;
win_size?=?1;
epsilon?=?0.000001;
%指定矩陣形狀
indsM=reshape（1:img_sizehw）;
%計算矩陣L
tlen?=?img_size*neb_size^2;
row_inds=zeros（tlen?1）;
col_inds=zeros（tlen1）;
vals=zeros（tlen1）;
len=0;
for?j=1+win_size:w-win_size
????for?i=win_size+1:h-win_size
????????if（rem（h-815）<1）
????????????if（rem（w-815）<1）
????????????????continue;
????????????end
????????end
??????win_inds=indsM（i-win_size:i+win_sizej-win_size:j+win_size）;
??????win_inds=win_inds（:）;%列顯示
??????winI=I（i-win_size:i+win_sizej-win_size:j+win_size:）;
??????winI=reshape（winIneb_sizec）;?%三個通道被拉平成為一個二維矩陣?9*3
??????win_mu=mean（winI1）‘;??%求每一列的均值?如果第二個參數為2?則為求每一行的均值??//矩陣變向量
??????win_var=winI‘*winI/neb_size-win_mu*win_mu‘+epsilon/neb_size*eye（c）;?%求方差
??????winI=winI-repmat（win_mu‘neb_size1）;%求離差
??????tvals=（1+winI/win_var*winI‘）/neb_size;%?求論文所指的矩陣L
??????row_inds（1+len:neb_size^2+len）=reshape（repmat（win_inds1neb_size）...
??????