function?pca?=?f_pcafusion（mulfilepanfile）
mul?=?imread（mulfile）;
%?mul=imresize（mul[512512]）;
pan?=?imread（panfile）;
%顯示原來多光譜圖與全色圖像
subplot（2?2?1）imshow（mul）;title（‘原多光譜圖像‘）;
subplot（2?22）imshow（pan）;title（‘原全色譜圖像‘）;
%預處理數據在0-1間
mul=double（mul）/255;
pan=rgb2gray（pan）;
pan=double（pan）/255;

%求相關矩陣
[r?c?bands]=size（mul）;%得到行、列、波段數
pixels?=?r*c;
%?reshape每個波段成一行
mul?=?reshape（mul?[pixels?bands]）;
correlation?=?sqrt（（mul‘*mul）/pixels）;%correlation?=?（mul‘*mul）/pixels;??%協方差矩陣/相關系數
%求特征向量矩陣與特征值對角矩陣
[vector?value]=eig（correlation）;

%計算最大特征值對應的特征向量（波段）
d=diag（value）;%將特征值取出，構成一個列向量
col=0;??%記錄特征值最大的波段
for?j=1:size（d）-1
????if（d（j）>d（j+1））?%當貢獻率大于95%時循環結束并記下取多少個特征值
????????col=j;
????else
????????col=j+1;
????end
end

?%將特征向量按降序排序
[dummyorder]=sort（diag（-value））;%因為sort函數是升序排列，所以先取負號diag（a）是取出a的對角元素構成
%?這里的dummy是降序排列后的對角線列向量，order是其大小順序
value=value（:order）;%將特征向量按照特征值大小進行降序排列

%求主分量圖像
PC?=?mul*vector;??%?Y=AX（X中列為樣本，若X行為樣本，則Y?=XA）
PC?=?reshape（PC[rcbands]）;
subplot（223）imshow（PC）title（‘多光譜主分量圖像‘）;
%根據第一主分量【直方圖配準】pan后代替第一主分量
%?PC（::bands）/max（max（PC（::bands）））就是使第一主成分的元素值在0-1之間

[countsX]?=?imhist（PC（::col）/max（max（PC（::col））））;%返回直方圖數據向量counts或相應的色彩值向量x
%?直方圖匹配
pan?=?histeq（pancounts）;%使全色圖的直方圖分布同第一主成分的直方圖一致
mm=max（max（PC（::col）））;
%PC（::col）?=?double（pan*max（max（PC（::col））））;??%線性變換

%PCA逆變換重構融合圖象
PC?=?reshape（PC[pixels?bands]）;
fusion?=?PC*vector‘;???%PC為行
fusion?=?reshape（fusion[rcbands]）;
%顯示融合圖象
subplot（224）imshow（fusion（::1:bands））;title（‘PCA變換融合多光譜與全色圖像‘）;
%f_pcafusion（‘bldr_tm_6.tif‘‘bldr_sp_1.tif‘）
end

?屬性????????????大小?????日期????時間???名稱
-----------?---------??----------?-----??----

?????文件????4499821??2018-11-04?18:48??pca_fusion\bldr_sp_1.tif

?????文件?????630859??2018-11-04?18:45??pca_fusion\bldr_tm_6.tif

?????文件???????1998??2018-11-04?20:28??pca_fusion\f_pcafusion.m

?????目錄??????????0??2018-11-04?20:36??pca_fusion

-----------?---------??----------?-----??----

??????????????5132678????????????????????4