Ідентифікація

Автор: Пользователь скрыл имя, 31 Марта 2013 в 20:31, курсовая работа

Описание работы

Маючи справу з системою, нам необхідна будь-яка схема співвідношення між собою, що характеризують систему змінних. В широкому розумінні ми будемо називати моделлю сукупність зв‘язків між сигналами, що спостерігаються. Очевидно, що моделі можуть приймати різноманітну форму і записуватись з різною степінню математичної деталізації. Вибір того рівня складності, який робить модель корисною, визначається запланованим використанням.

Работа содержит 1 файл

курсач.docx

— 297.13 Кб (Скачать)

 

Продовження таблиці 3.1

Номер точки

Сигнал  Y(t)

Номер точки

Сигнал  Y(t)

t

Y(t)

t

Y(t)

91

130,225724

1,39381767

98

138,4322268

1,126083919

92

131,167454

0,889482

99

139,2394237

1,150989385

93

132,37825

1,22881896

100

140,4502191

0,789860135

94

133,589045

1,36268584

101

141,5264817

1,11674437

95

134,79984

0,95797203

102

143,1408756

1,306648544

96

136,010636

0,8708029

103

144,0826054

1,026462057

97

137,086899

1,4809868

104

145,158868

1,070046622


 

3.2 Побудова моделі за допомогою поліноміальної інтерполяції поліномом 7-го степеня

Побудова моделі за допомогою  поліноміальної інтерполяції поліномом 7-го степеня виконуємо тим же способом як і для вхідного сигналу U(t), що описаний в першому розділі.

Нижче наведено лістинг програмного  модуля в середовищі Matlab:

 

clear;

load wt3.mat

for i=1:8;

    for j=1:8;

        A(i,j)=sum(x.^(i+j-2));

    end;

end;

for i=1:8;

    B(i)=sum(x.^(i-1).*y);

end;

a=B/A;

t=1:150;

f=a(8)*t.^7+a(7)*t.^6+a(6)*t.^5+a(5)*t.^4+a(4)*t.^3+a(3)*t.^2+a(2)*t+a(1);

hPlot=plot(t,f);

hold on;

plot(x,y,'red')

xlabel('t')

ylabel('Y(t)')

axis([0,150,-0.1,2.8]);

set( hPlot, 'LineWidth', 2 );

 

Нижче наведемо графічну модель (рис. 3.1):

Рис. 3.1. Поліноміальна модель вихідного сигналу Y(t) поліномом n-го порядку

 

 

Розв’язавши систему, отримаємо  такі значення коефіцієнтів:

 

a0 = 1,02538836673611;

a1 = 0,0268564689044343;

a2 = -0,00252492654257438;

a3= 9,12939501281913e-05;

a4 = -1,58612302537886e-06;

a5 = 1,41286442406615e-08; 

a6 = -6,21079558699196e-11;

a7= 1,06636334655407e-13.

 

3.3 Побудова моделі за допомогою полінома Чебишева

Побудова моделі за допомогою  полінома Чебишева виконуємо тим же способом як і для вхідного сигналу U(t), що описаний в першому розділі.

Нижче наведено лістинг програмного  модуля в середовищі Matlab:

 

clear;

load wt3.mat;

Cheb_por=7;

n=104;

fi=ones(n,7);

for i=1:n

    fi(i,2)=x(i)-(1/n)*sum(x);

end;

for p=3:Cheb_por; i=1:n;

        beta(p)=-sum(x(i).*fi(i,p-1).^2)/sum(fi(i,p-1).^2);

        gama(p)=-sum(x(i).*fi(i,p-1).*fi(i,p-2))/sum(fi(i,p-1).^2);

        fi(i,p)=(x(i)+beta(p)).*fi(i,p-1)+gama(p).*fi(i,p-2);

end;

A=fi'*fi;

B=fi'*y;

a=A\B;

for i=1:n; j=1:7; 

        Cheb(i)=a(j)'*fi(i,j)';

end;

hPlot=plot(x,Cheb);

hold on;

plot(x,y,'red');

xlabel('t');

ylabel('Y(t)');

set( hPlot, 'LineWidth', 2 );

 

Нижче наведено модель досліджуваного об’єкта (рис. 3.2).

Розв’язавши систему, отримаємо  такі значення коефіцієнтів:

a1 = 1,05675426115593;           a5 = 4,21939951118780e-08;

a2 = -0,00317011396636858;           a6 = -3,81381034627028e-10;

a3 = -6,21837031383335e-05;           a7 = -5,66801060511657e-12.

a4 = 2,40061957256621e-06;

 

Рис. 3.2. Модель вихідного сигналу Y(t), побудована поліномом Чебишева

 

3.4 Побудова моделі за допомогою перетворення Фур’є

Побудова моделі за допомогою  перетворень Фур’є виконуємо  тим же способом як і для вхідного сигналу U(t), що описаний в першому розділі.

Нижче наведено лістинг програмного  модуля в середовищі Matlab:

 

clear;

load wt3.mat;

n=104;

Fur_por=4;

m=Fur_por;

for k=1:m+1

    A=(1/n)*sum(y);

    B(k)=(2/n)*sum(y.*cos(k*x));

    C(k)=(2/n)*sum(y.*sin(k*x));

end

for k=1:m+1

    for i=1:n

        FS(i,1)=A+sum(B(k)*sin(k*x(i))+C(k)*cos(k*x(i)));

    end

end

plot(x,y,'red')

hold on;

hPlot=plot(x,FS);

xlabel('t')

ylabel('Y(t)')

axis([0,150,-0.1,2.1]);

set( hPlot, 'LineWidth', 2 );

 

Розв’язавши систему, отримаємо  такі значення коефіцієнтів:

 

a0 = 1,06570612602025;

 

 

Побудована модель представлена нижче, на рис. 3.3.

Рис. 3.3. Модель вихідного сигналу Y(t), побудована за допомогою перетворення  Фур’є

 

3.5 Статистична обробка даних

Розрахуємо коефіцієнти  дисперсії, кореляції та коваріації.

Розрахунок виконуємо  тим же способом як і для вхідного сигналу U(t), що описаний в першому розділі.

Нижче наведено лістинг програмного  модуля в середовищі Matlab:

 

clear;

load wt3.mat;

n=104;

i=1:n;

tser=sum(x)/n;

fser=sum(y)/n;

cov=sum((x-tser).*(y-fser))/(n-1)

Dt=sum((x-tser).^2)/(n-1)

Du=sum((y-fser).^2)/(n-1)

r=cov/sqrt(Dt*Du)

 

Після виконаного розрахунку отримаємо:

 

cov = -0.0393;

Dt = 1.7398e+003;

Du = 0.0574;

r = -0.0039

3.6 Знаходження періодограми сигналів

Побудову періодограми виконуємо тим же способом як і для вхідного сигналу U(t), що описаний в першому розділі.

Нижче наведено лістинг програмного  модуля в середовищі Matlab:

 

clear;

load wt3.mat;

n=104;

for k=1:n

    w(k,1)=(2*pi*(k))/n;

    U(k,1)=(1/sqrt(n))*sum(y.*exp(-j*w(k).*x));

end

Ui=imag(U);

Ur=real(U);

subplot(1,2,1);

plot(w,Ui);

title('a)');

xlabel('Нормована частота (Гц)');

ylabel('Спектр щільності (Дб/Гц)');

grid on

subplot(1,2,2);

plot(w,Ur);

title('б)');

xlabel('Нрмована частота (Гц)');

ylabel('Спектр щільності (Дб/Гц)');

grid on

 

 

Рис. 3.4 - Періодограма вихідного сигналу Y(t), де а) – уявна частина,

б) – дійсна частина 

3.7 Вибір оптимальної моделі

Вибір оптимальної моделі виконуємо тим же способом як і  для вхідного сигналу U(t), що описаний в першому розділі.

 

Оптимальні значення моделей

Таблиця 3.2

Вид функції

Поліном

Поліном Чебишева

Функція Фур’є

Значення критерію

0.0085

4.2762e-015

0.0639


Нижче наведено лістинг програмного  модуля в середовищі Matlab:

 

for i=1:a;n=104;

    Spol_Yt=sqrt(sum((y(i)-f(i)).^2/(n-2)))

end;

     SCheb_Yt=sqrt((sum(y-Cheb').^2)/(n-2))

 

     SFur_Yt=sqrt((sum(y-FS).^2)/(n-2))

 

Як видно з таблиці  3.2, значення критерію мінімальне у тому випадку, коли апроксимуючою функцією виступає Поліном Чебишева, отже модель побудована даним чином буде оптимальною.

 

Розділ  4. Побудова загальної моделі системи

 

4.1 Табулювання двох вхідних та одного вихідного сигналів

Табулювання графіків функцій (сигналів) наведені в попередніх розділах, тому таблиці та функції табулювання  вказувати не будемо. Лише зазначимо  деякі особливості.

Розіб‘ємо часову вісь t на N проміжків з інтервалом 2 секунди. Для кожного інтервалу визначимо значення ур оптимальних моделей для сигналів U(t), V(t), Y(t).

        i:=1.. n       ti:=2.i

 

для U(t) – оптимальною моделлю є модель побудована за допомогою полінома Чебишева;

V(t) – оптимальною моделлю є модель побудована за допомогою полінома 7-ї степені;

Y(t) – оптимальною моделлю є модель побудована за допомогою полінома Чебишева.

        Спробуємо  визначити, яка з слідуючих моделей буде оптимальною для даної системи:

 

    1. Y = a0+a1 U(t)+a2W(t);
    2. Y = a0+a1 U(t)+a2W(t)+ a3U(t)×W(t);
    3. Y = a0 +b .cos(U(t)) + d .sin(W(t)).

Для цього використаємо метод  найменших квадратів.

 

4.2 Побудова першої моделі системи, виду

Y = a0+a1 U(t)+a2W(t)

За методом найменших квадратів  система рівнянь для першої моделі буде мати вигляд:

 

 

Нижче приведемо лістинг модуля рішення даної системи:

 

clear;

load wt1.mat

for i=1:8;

    for j=1:8;

        A(i,j)=sum(x.^(i+j-2));

    end;

end;

for i=1:8;

    B(i)=sum(x.^(i-1).*y);

end;

a1=B/A;

t=1:150;

f1=a1(8)*t.^7+a1(7)*t.^6+a1(6)*t.^5+a1(5)*t.^4+a1(4)*t.^3+a1(3)*t.^2+a1(2)*t+a1(1);

load wt2.mat

for i=1:8;

    for j=1:8;

        A(i,j)=sum(x.^(i+j-2));

    end;

end;

for i=1:8;

    B(i)=sum(x.^(i-1).*y);

end;

a2=B/A;

f2=a2(8)*t.^7+a2(7)*t.^6+a2(6)*t.^5+a2(5)*t.^4+a2(4)*t.^3+a2(3)*t.^2+a2(2)*t+a2(1);

load wt3.mat

for i=1:8;

    for j=1:8;

        A(i,j)=sum(x.^(i+j-2));

    end;

end;

for i=1:8;

    B(i)=sum(x.^(i-1).*y);

end;

a3=B/A;

f3=a3(8)*t.^7+a3(7)*t.^6+a3(6)*t.^5+a3(5)*t.^4+a3(4)*t.^3+a3(3)*t.^2+a3(2)*t+a3(1);

sf1=0;sf2=0;sf3=0;sf12=0;sf22=0;sf1f2=0;sf1f3=0;sf2f3=0;

n=150;

sf1=sum(f1);

sf2=sum(f2);

sf3=sum(f3);

sf12=sum(f1.^2);

sf22=sum(f2.^2);

sf1f2=sum(f1.*f2);

sf1f3=sum(f1.*f3);

sf2f3=sum(f2.*f3);

A1=[n sf1 sf2;sf1 sf12 sf1f2;sf2 sf1f2 sf22];

B1=[sf3 sf1f3 sf2f3];

a=B1/A1;

sysmod1=a(3).*f2+a(2).*f1+a(1);

hPlot=plot(t,f3,'red');

hold on;

hPlot1=plot(t,sysmod1);

grid on;

xlabel('t')

ylabel('Y(t)')

set( hPlot,'LineWidth', 2 );

set( hPlot1,'LineWidth', 2 );

S1=sqrt((sum(f3-sysmod1)^2)/(n-2));

 

Нижче побудуємо її графічне представлення  і знайдемо коефіцієнти а0, а1 , і а2 , які визначатимуть модель:

 

a0 =0,0288517639895134;

a1 =0,567158190502034;

a2 =0,424033755053556.

 

Рис. 4.1 – Загальна модель системи виду

Y = a0+a1 U(t)+a2W(t)

 

4.3 Побудова другої моделі системи, виду

Y = a0+a1 U(t)+a2W(t)+a3U(t)W(t)

За методом найменших квадратів  система рівнянь для першої моделі буде мати вигляд:

 

 

Нижче приведемо лістинг модуля рішення даної системи:

 

clear;

load wt1.mat

for i=1:8;

    for j=1:8;

        A(i,j)=sum(x.^(i+j-2));

    end;

end;

for i=1:8;

    B(i)=sum(x.^(i-1).*y);

end;

a1=B/A;

t=1:150;

f1=a1(8)*t.^7+a1(7)*t.^6+a1(6)*t.^5+a1(5)*t.^4+a1(4)*t.^3+a1(3)*t.^2+a1(2)*t+a1(1);

load wt2.mat

for i=1:8;

    for j=1:8;

        A(i,j)=sum(x.^(i+j-2));

    end;

end;

for i=1:8;

    B(i)=sum(x.^(i-1).*y);

end;

a2=B/A;

f2=a2(8)*t.^7+a2(7)*t.^6+a2(6)*t.^5+a2(5)*t.^4+a2(4)*t.^3+a2(3)*t.^2+a2(2)*t+a2(1);

load wt3.mat

for i=1:8;

    for j=1:8;

        A(i,j)=sum(x.^(i+j-2));

    end;

end;

for i=1:8;

    B(i)=sum(x.^(i-1).*y);

end;

a3=B/A;

f3=a3(8)*t.^7+a3(7)*t.^6+a3(6)*t.^5+a3(5)*t.^4+a3(4)*t.^3+a3(3)*t.^2+a3(2)*t+a3(1);

n=150;

hold off;

sf1=sum(log(f1));

sf2=sum(log(f2));

sf1f1=sum(log(f1).^2);

sf1f2=sum(log(f1).*log(f2));

sf1sf2=sum(sf1.*sf2);

sf1f1f2=sum((log(f1).^2).*log(f2));

sf2f2=sum(log(f2).^2);

sf2f2f1=sum((log(f2).^2).*log(f1));

sf2f2sf1f1=sum((log(f1).^2).*(log(f2).^2));

sf3=sum(log(f3));

sf1f3=sum(log(f1).*log(f3));

sf2f3=sum(log(f2).*log(f3));

sf3f3=sum(log(f3).^2);

A1=[n sf1 sf2 sf1sf2;sf1 sf1f1 sf1f2 sf1f3;sf2 sf1f2 sf2f2 sf2f3;sf1sf2 sf1f1f2 sf2f2f1 sf2f2sf1f1];

B1=[sf3 sf1f3 sf2f3 sf3f3];

a=B1/A1;

sysmod2=a(4).*f1.*f2+a(3).*f2+a(2).*f1+a(1);

hPlot=plot(t,f3,'red');

hold on;

hPlot1=plot(t,sysmod2);

grid on;

xlabel('t')

ylabel('Y(t)')

set( hPlot,'LineWidth', 2 );

set( hPlot1,'LineWidth', 2 );

S2=sqrt((sum(f3-sysmod2)^2)/(n-2));

 

Розв’язавши її, знайдемо коефіцієнти а0, а1, а2 і а3, які визначатимуть модель:

 

a0 = 0,00317016777990134;

a1 =0,665835193742437;

a2 =0,583053398632476;

a3 =0,0181216832221273.

 

Нижче  побудуємо графічне представлення моделі, на якому видно  її співвідношення до реального об’єкта.

 

 

Рис. 4.2 – Загальна модель системи виду

Y = a0+a1 U(t)+a2W(t)+a3W(t)W(t)

 

4.4 Побудова третьої моделі системи, виду

Y = a0 +b×cos(U(t)) + d×sin(W(t))

За допомогою методу найменших  квадратів знаходимо систему  рівнянь, розв’язавши яку можна  отримати коефіцієнти а0, а1, та а2. За допомогою знайдених коефіцієнтів можна побудувати третю модель.

Система рівнянь буде мати вигляд:

 

Нижче приведемо лістинг модуля рішення даної системи:

 

clear;

load wt1.mat

for i=1:8;

    for j=1:8;

        A(i,j)=sum(x.^(i+j-2));

    end;

end;

for i=1:8;

    B(i)=sum(x.^(i-1).*y);

end;

a1=B/A;

t=1:150;

f1=a1(8)*t.^7+a1(7)*t.^6+a1(6)*t.^5+a1(5)*t.^4+a1(4)*t.^3+a1(3)*t.^2+a1(2)*t+a1(1);

load wt2.mat

for i=1:8;

    for j=1:8;

        A(i,j)=sum(x.^(i+j-2));

    end;

end;

for i=1:8;

    B(i)=sum(x.^(i-1).*y);

end;

a2=B/A;

f2=a2(8)*t.^7+a2(7)*t.^6+a2(6)*t.^5+a2(5)*t.^4+a2(4)*t.^3+a2(3)*t.^2+a2(2)*t+a2(1);

load wt3.mat

for i=1:8;

    for j=1:8;

        A(i,j)=sum(x.^(i+j-2));

    end;

end;

for i=1:8;

    B(i)=sum(x.^(i-1).*y);

end;

a3=B/A;

f3=a3(8)*t.^7+a3(7)*t.^6+a3(6)*t.^5+a3(5)*t.^4+a3(4)*t.^3+a3(3)*t.^2+a3(2)*t+a3(1);

n=150;

Fur_por=3;

m=Fur_por;

for k=1:m+1

    Am=(1/n)*sum(f3);

Информация о работе Ідентифікація