%**************************************************************************************************************************************************************
%   本程序為主程序,目的是運用Taylor變換對WM1000發動機的前13階質量、剛度及阻尼等“偽”結構參數進行識別,以區別正常信號和故障信號的差異來進行故障診斷
%   
%   中國北方發動機研究所試驗測試中心  蘇志霄   2006年10月
%**************************************************************************************************************************************************************
%
speed=1500;                                                       %發動機轉速
m=13;                                                             %m取值為測點數目
H=1;                                                              %Taylor變換中的常數，該數值對計算結果影響很小，故賦值為1
%**************************************************************************************************************************************************************
%讀入時間加速度數據,各測點在每個測試時間內時間數據列是一致的，因此統一標識
t=acceleration_X;
len=max(size(t));
acc=zeros(m,len);
acc(1,:)=acceleration_Y';acc(2,:)=acceleration2_Y';acc(3,:)=acceleration3_Y';acc(4,:)=acceleration4_Y';acc(5,:)=acceleration5_Y';acc(6,:)=acceleration6_Y';
acc(7,:)=NONE_Y';acc(8,:)=NONE2_Y';acc(9,:)=NONE3_Y';acc(10,:)=NONE4_Y';acc(11,:)=NONE5_Y';acc(12,:)=NONE6_Y';acc(13,:)=NONE7_Y';
%將加速度數據由單精度格式轉換為雙精度格式
acc=double(acc);
%輸入“偽”外力，假設僅在第11號傳感器部位作用了一個轉速頻率的正弦外力
omg=2*pi*speed/60;
%F11=sin(omg*t);
%**************************************************************************************************************************************************************
%取時間范圍從1/5開始為實際有效計算時間
t=t(round(len/5):len);length=max(size(t));acc=acc(:,round(len/5):len);%F11=acc(F11(len/5):len);
%Taylor變換中的兩個系數變量kc(k)和yt(k)，由于使用時間區域均勻分割的識別方法，因此這里只需計算kc(0)和yt(0)
kc=(0+1)*(0+2)/H^2;yt=(0+1)/H;
%將計算時間分為3倍測點數目個小時間區域
dengfen=round(length/(3*m)-3*m);                        %為防止取整后會造成丟點問題而造成點數不夠，故每段區域的點數適當減少(-3*m)
t=t(1:dengfen+1)';t=t-t(1);length=max(size(t));         %由于時間區域均勻分割，因此每子段內的時間點數是相同的，故統一代表
%計算Taylor變換中的系數矩陣
xishujuzhen=zeros(3*m);
for i=1:m;
    i,
    for j=1:3*m;
        for k=1:length;
            denfenxinhao(k)=acc(i,(j-1)*dengfen+k);
        end;
        xishujuzhen(i,j)=kc*taylor(t,denfenxinhao,2,H);
    end;
end;
for i=m+1:2*m;
    i,
    for j=1:3*m;
        for k=1:length;
            denfenxinhao(k)=acc(i-m,(j-1)*dengfen+k);
        end;
        xishujuzhen(i,j)=yt*taylor(t,denfenxinhao,1,H);
    end;
end;
for i=2*m+1:3*m;
    i,
    for j=1:3*m;
        for k=1:length;
            denfenxinhao(k)=acc(i-2*m,(j-1)*dengfen+k);
        end;
        xishujuzhen(i,j)=taylor(t,denfenxinhao,0,H);
    end;
end;
%輸入“偽”外力的Taylor變換項
%FT=zeros()