function [Umat,Cmat,it,valJ] = fcm(Xdat,C,mfuz,epsilon,itmax);
%   function [Umat,Cmat,it,valJ] = fcm(Xdat,K,mfuz,epsilon,itmax);
%   模糊聚類主函數%   
%   Xdat ：原數據集，大小為(M,N)
%   C    ：聚類數目
%   M：數據的個數，N：數據的維數
%   mfuz ：加權指數m
%   Umat ：模糊劃分矩陣，大小為(C,M)
%   Cmat ：聚類中心矩陣，大小為(C,N)
%   Dmat ：距離矩陣，大小為(C,M)   
%   epsilon ：運算中止條件
%   itmax   ：最大迭代次數
%   it      ：迭代次數

load iris.txt;
Xdat = iris(:,3:6);
%Xdat=[0.9637,0.2536,0.7415;0.4659,1.0000,0.5815;1.0000,0.7893,1.0000;0.5437,0.9915,0.0000;0.3062,0.6861,0.2364];
C=3;
mfuz=2.5;
epsilon=0.001;
itmax=100;
[M,N] = size(Xdat);

Cmat = fcm_dataInitC(Xdat,C);%   初始化聚類中心（聚類原型模式）
Ct=Cmat;
it=0;
stab = 0;%  迭代中止標記0為繼續1為中止
%  迭代求解
while (it<itmax & stab ~=1)

Dmat = fcm_calcD(Xdat,Cmat);%   計算樣本與聚類中心（原型模式）間的距離矩陣

Umat = fcm_calcU(mfuz,Dmat);%   更新模糊劃分矩陣

Cmat = fcm_calcC(Xdat,C,mfuz,Umat);%   更新聚類中心（原型模式）

Ct=Ct-Cmat;
valJ = fcm_evalJ(Umat,mfuz,Dmat);
a=it+1;
Jmat(a,:)=valJ;
fnorm = norm(Ct,'fro');%   計算是否滿足迭代中止條件
    if (fnorm<epsilon)  
     stab = 1;  
    else   
     Ct=Cmat;    
    end
    it = it+1;
end
it
Jmat
Cmat

return %% end of function fcm.m