Міністерство освіти і науки,молоді та спорту України
Національний університет «Львівська політехніка»
Кафедра ІВТ
Розрахункова робота
з курсу :
«Комп’ютерне опрацювання вимірювальної інформації»
на тему:
«Комплексна частотна характеристика. Дослідження методом Монте-Карло»
Варіант – 1
ОСНОВНІ ТЕОРЕТИЧНІ ВІДОМОСТІ
Комплексну частотну характеристику лінійної стаціонарної динамічної системи можна отримати у вигляді прямого перетворення Фур’є вхідного та вихідного сигналів, взявши їх відношення.
Нехай маємо вхідний сигнал x(t) та вихідний сигнал y(t). Потрібно знайти комплексну частотну характеристику , де - уявна одиниця (рис.7.1).
Візьмемо перетворення Фур’є від вхідного та вихідного сигналів:
, (7.1(
, (7.2(
Тоді, комплексна частотна характеристика за означенням дорівнюватиме
. (7.3(
Вираз для комплексної частотної характеристики можна також отримати, використавши перетворення Лапласа.
Зображенням Лапласа функції х(t) дійсного аргумента t називається функція комплексного аргумента p=s+jτ, яка вираховується за формулою:
, (7.4(
Властивості перетворення Лапласа
Якщо дві функції f1(t) та f2(t) мають одне і теж зображення F(p) - f1(t)↔F(p), f2(t)↔F(p) то f1(t) = f2(t).
Якщо f1(t)↔F1(p), a f2(t)↔F2(p) то лінійніq комбінаціїq af1(t)+bf2(t)↔aF1(p)+ bF2(p), що є наслідком лінійності інтегралу.
Таблиці оригіналів та зображень Лапласа можна знайти у математичних довідниках. Якщо знайти зображення Лапласа вхідного та вихідного сигналів, взяти їх відношення та провести заміну комплексного аргумента p на то отримаємо комплексну частотну характеристику.
Комплексну частотну характеристику можна подати в алгебраїчній
, (7.5(
тригонометричній
(7.6(
та полярній
(7.7(
формах.
Дійсна та уявна частина комплексної частотної характеристики виходячи з (7.5( та (7.6( дорівнюють
(7.8(
а геометричне представлення подано на рис.7.2.
Найбільш широко на практиці використовуються дві характеристики – амплітудно-частотна характеристика
, (7.9(
та фазо-частотна характеристика.
. (7.10(
Дані характеристики можна отримати безпосередньо з комплексної частотної характеристики. Наприклад, отримано наступний вираз комплексної частотної характеристики:
, (7.11(
Домножимо чисельник та знаменник у (7.11( на комплексно спряжений вираз і згрупуємо елементи, отримавши дійсну та уявну частини комплексної частотної характеристики
, (7.12(
Тоді, підставивши у вирази (7.9( та (7.10( значення , отримані в (7.12(, знаходимо амплітудно-частотну характеристику
(7.13(
та фазо-частотну характеристику
. (7.14(
Якщо вираз комплексної частотної характеристики простий то отримати амплітудно-частотну та фазо-частотну характеристику легко. Проте бувають випадки, коли маємо складний вираз. У такому випадку можна не шукати математичні залежності для амплітудно-частотної та фазо-частотної характеристик, а скористатися наступним алгоритмом з використанням комплексного числення:
Задаємо крок зміни кругової частоти
. (7.15(
де - максимальна частота, - мінімальна частота, - кількість розбиттів у діапазоні частот.
2. Знаходимо біжуче значення кругової частоти
. (7.16(
3. Шукаємо значення комплексного числа Ki для кругової частоти
, (7.17(
4. Рахуємо дійсну та уявну частини комплексного числа Ki
, (7.18(
5. Шукаємо значення амплітудно-частотної характеристики для кругової частоти
, (7.19(
6. Шукаємо значення фазо-частотної характеристики для кругової частоти
. (7.20(
Повторивши пункти 2-6 разів отримуємо дискретні значення амплітудно-частотної та фазо-частотної характеристик з кроком , які графічно можна вивести на монітор комп’ютера.
Для дослідження похибок перетворення амплітудно-частотної характеристики в залежності від точності виготовлення елементів електронних схем можна використати метод Монте-Карло, який дозволяє досліджувати математичні моделі об’єктів, використовуючи математичний апарат статистичного аналізу випадкових процесів та величин.
Якщо отримано деяку математичну модель реального процесу в залежності від деяких параметрів , які можуть приймати випадкові значення
, (7.21(
то дану систему можна досліджувати методом Монте-Карло. Щоб отримати правильні результати дослідження , потрібно мати попередню інформацію про закон розподілу випадкових величин .
Дослідження амплітудно-частотної характеристики (7.13( методом Монте-Карло можна проводити за наступним алгоритмом:
Задаємо випадкові значення опору резистора та ємності конденсатора на конкретній частоті
. (7.22(
де – поточне значення опору резистора та ємності конденсатора, – номінальні значення опору резистора та ємності конденсатора, – точність виготовлення резистора та конденсатора, - функція, яка генерує випадкові значення у діапазоні від 0 до 32767.
2. Для отриманих випадкових значень опору резистора та ємності конденсатора розраховуємо комплексне значення комплексної частотної характеристики
, (7.23(
3. Рахуємо дійсну та уявну частини комплексного числа Ki
, (7.24(
4. Шукаємо значення амплітудно-частотної характеристики для частоти
, (7.25(
5. Повторити пункти 1 - 4 разів
6. Знаходимо математичне сподівання випадкових значень амплітудно-частотної характеристики
. (7.26(
7. Знаходимо дисперсію випадкових значень амплітудно-частотної характеристики
. (7.27(
8. За відомої дисперсії записуємо довірчий інтервал для амплітудно-частотної характеристики
, (7.28(
де — номінальне значення амплітудно-частотної характеристики на частоті , пораховане для номінальних значень резистора та конденсатора, — довірчий коефіцієнт, значення якого залежить від довірчої ймовірності .
Варіант 1
На мові С це буде виглядати наступним чином: #include<graphics.h>
#include<stdio.h>
#include<stdlib.h>
#include<conio.h>
#include<math.h>
#include<complex.h>
#define k 50
#define s 10
void main()
{
int R1=100,R4=1000,fmin=1,fmax=1000;
float C2=5e-6,L3=5e-3,C4=5e-6;
int xdmin=10,xdmax=639,ydmin=10,ydmax=479,i,o;
float xmin,xmax,ymin,ymax,xd,yd;
float wmin,wmax,w,krok,P,Q,D,R1i,R4i,C2i,C4i,L3i;
float wi[k],Ai[k],Fi[k],Mx[k],da[k],pohubka[s];
complex z,j;
int gdriver = DETECT, gmode, errorcode;
initgraph(&gdriver, &gmode,"");
errorcode = graphresult();
if (errorcode != grOk)
{
printf("Graphics error: %s\n", grapherrormsg(errorcode));
printf("Press any key to halt:");
getch();
exit(1);
}
wmin=2*3.14*fmin;
wmax=2*3.14*fmax;
krok=(wmax-wmin)/k;
j=complex(0,1);
z=complex(1,1);
for(i=0;i<k;i++)
{
wi[i]=wmin+i*krok;
z=(R4/(R4*j*wi[i]*C4+1))/(1/(R1*j*wi[i]*C2+1)+j*wi[i]*L3+R4/(R4*j*wi[i]*C4+1));
P=real(z);
Q=imag(z);
Ai[i]=sqrt(P*P+Q*Q);
Fi[i]=atan(Q/P);
}
/* ach */
xmax=wi[0];
xmin=wi[0];
ymax=Ai[0];
ymin=Ai[0];
for(i=0;i<k;i++)
{
if(wi[i]<xmin) xmin=wi[i];
if(wi[i]>xmax) xmax=wi[i];
if(Ai[i]<ymin) ymin=Ai[i];
if(Ai[i]>ymax) ymax=Ai[i];
}
xd=xdmin+(xdmax-xdmin)/(xmax-xmin)*(wi[0]-xmin);
yd=ydmax-(ydmax-ydmin)/(ymax-ymin)*(Ai[0]-ymin);
moveto(xd,yd);
circle(xd,yd,2);
for(i=1;i<k;i++)
{
xd=xdmin+(xdmax-xdmin)/(xmax-xmin)*(wi[i]-xmin);
yd=ydmax-(ydmax-ydmin)/(ymax-ymin)*(Ai[i]-ymin);
lineto(xd,yd);
circle(xd,yd,2);
}
printf("\n\n\n\n\n ACH");
getch();
cleardevice();
/* fch */
ymax=Fi[0];
ymin=Fi[0];
for(i=0;i<k;i++)
{
if(wi[i]<xmin) xmin=wi[i];
if(wi[i]>xmax) xmax=wi[i];
if(Fi[i]<ymin) ymin=Fi[i];
if(Fi[i]>ymax) ymax=Fi[i];
}
xd=xdmin+(xdmax-xdmin)/(xmax-xmin)*(wi[0]-xmin);
yd=ydmax-(ydmax-ydmin)/(ymax-ymin)*(Fi[0]-ymin);
moveto(xd,yd);
circle(xd,yd,2);
for(i=1;i<k;i++)
{
xd=xdmin+(xdmax-xdmin)/(xmax-xmin)*(wi[i]-xmin);
yd=ydmax-(ydmax-ydmin)/(ymax-ymin)*(Fi[i]-ymin);
lineto(xd,yd);
circle(xd,yd,2);
}
printf("\n\n\n\n\n FCH");
getch();
cleardevice();
/* monte-carlo */
pohubka[0]=0.00001;
randomize();
w=(wmax-wmin)/2;
for(o=0;o<s;o++)
{
pohubka[o]=(float)o/100.0;
for(i=0;i<k;i++)
{
R1i=R1*(1.+(rand()/32767.-0.5)*2.*pohubka[o]);
R4i=R4*(1.+(rand()/32767.-0.5)*2.*pohubka[o]);
C2i=C2*(1.+(rand()/32767.-0.5)*2.*pohubka[o]);
C4i=C4*(1.+(rand()/32767.-0.5)*2.*pohubka[o]);
L3i=L3*(1.+(rand()/32767.-0.5)*2.*pohubka[o]);
z=(R4i/(R4i*j*w*C4i+1))/(1/(R1i*j*w*C2i+1)+j*w*L3i+R4i/(R4i*j*w*C4i+1));
P=real(z);
Q=imag(z);
Ai[i]=sqrt(P*P+Q*Q);
}
Mx[o]=0;
for(i=0;i<k;i++) Mx[o]+=Ai[i];
Mx[o]/=k;
D=0;
for(i=0;i<k;i++) D+=pow((Ai[i]-Mx[o]),2);
D/=(k-1);
da[o]=3.0*sqrt(D)/Mx[o]*100;
}
xmax=pohubka[0];
xmin=pohubka[0];
ymax=da[0];
ymin=da[0];
for(i=0;i<s;i++)
{
if(pohubka[i]<xmin) xmin=pohubka[i];
if(pohubka[i]>xmax) xmax=pohubka[i];
if(da[i]<ymin) ymin=da[i];
if(da[i]>ymax) ymax=da[i];
}
xd=xdmin+(xdmax-xdmin)/(xmax-xmin)*(pohubka[0]-xmin);
yd=ydmax-(ydmax-ydmin)/(ymax-ymin)*(da[0]-ymin);
moveto(xd,yd);
circle(xd,yd,2);
for(i=1;i<20;i++)
{
xd=xdmin+(xdmax-xdmin)/(xmax-xmin)*(pohubka[i]-xmin);
yd=ydmax-(ydmax-ydmin)/(ymax-ymin)*(da[i]-ymin);
lineto(xd,yd);
circle(xd,yd,2);
}
getch();
cleardevice();
for(i=s-1;i>=0;i--)
{
if(Mx[i]/10.0>=da[i]/100*Mx[i])
{
printf("pohubka R1;R2;C;L=%3.f%%\n",pohubka[i]*100);
printf("rezultat =%3.4f+/-%3.4f",Mx[i],da[i]/100*Mx[i]);
break;
}
}
getch();
}