Если у вас возникли сложности с курсовой, контрольной, дипломной, рефератом, отчетом по практике, научно-исследовательской и любой другой работой - мы готовы помочь.
Предоплата всего
Подписываем
Уральский технический институт связи и информатики (филиал)
Сибирского государственного университета
Телекоммуникаций и информатики
(УрТИСИ ГОУ ВПО «СибГУТИ»)
Курсовая работа.
По дисциплине
«Языки программирования».
Тема: Спектральный анализ дискретных сигналов
Выполнил: Саляхутдинов Р. Г.
Группа: МЕ-51
Екатеринбург, 2008
Задание для курсовой работы.
Написать программу на языке программирования Паскаль для решения следующей задачи (вариант задания индивидуальный). Результаты расчетов должны выводиться на экран и в файл. Оформление графиков и таблиц выполнять средствами математических пакетов (Maple, MathCad). Демонстрационный вариант программы подготовить в среде визуального программирования Delphi.
Проверить решение промежуточных задач средствами математических пакетов. Построить блок-схемы задачи и вспомогательных частей алгоритма. Оценить погрешность выполненных расчетов.
Напряжение U=U(t) на входе транзистора как функция времени описывается дифференциальным уравнением
с начальными условиями (1), где n – последняя цифра номера зачетной книжки, k – коэффициент усиления (см. ниже), fs(t) – периодический сигнал
(рис. 1), m – коэффициент обратной связи.
Указания и пояснения.
Курсовая работа выполняется в ЧЕТЫРЕ этапа.
Оформление:
Часть 1.
> restart;
> with(linalg):with(plots):with(plottools):
> fun:= proc(t) local z ;z:=piecewise(t<0,0,t<1/2,cos(Pi*t),t<=1,0,0); evalf(z);end;
> fun(t/tau);
> p(x):=x^4-x^3-2;
> Koeff:=fsolve(p(x)=0,x,0..3);
> tau:=Koeff:
> plot(p(x),x=Koeff-0.5..Koeff+0.5,thickness=2,color=black);
> R1:=plot(fun(t),t=0..2.5,thickness=2,linestyle=3,color=blue):
> R11:=plot(fun(t/tau),t=0..2.5,thickness=2,color=black):
> display(R1,R11);
> Fourier_T:=proc(F,T0,TN,k::evaln) local T;
global A0,Ak,Bk;
T:=TN-T0;
A0:=2/T*Int(F(x),x=T0..TN);
Ak:=2/T*int(F(x)*cos(k*x*2*Pi/T),x=T0..TN):
Bk:=2/T*int(F(x)*sin(k*x*2*Pi/T),x=T0..TN):
end proc:
> Trig_polynom:=proc(N,T,a0,ak,bk,k::evaln) local n,Pol,A0,A,B;
global a,b,RisTrig;
a:=array(0..N);b:=array(0..N);
A0:=evalf(a0);a[0]:=A0;b[0]:=0;
A:=seq(evalf(subs(k=n,ak)),n=1..N);
B:=seq(evalf(subs(k=n,bk)),n=1..N);
for n from 1 to N do
a[n]:=A[n];b[n]:=B[n];
end do;
Pol:=A0/2+sum(A[k]*cos(2*Pi*k*x/T)+B[k]*sin(2*Pi*k*x/T),k=1..N):
RisTrig:=plot(Pol,x=-T/2..3*T/2,color=blue,thickness=2):
RETURN(Pol);
end proc:
> ARR:=proc(n::integer,c) local L,H,ma,mi,k::integer,
Sim::array(0..n);
ma:=c[0];mi:=c[0];
L:=line([0,c[0]],[n,c[n]],thickness=2,color=red);
for k from 1 to n do
if c[k]>ma then ma:=c[k];end if;
if c[k]<mi then mi:=c[k];end if;
end do;
H:=ma-mi;
if H=0 then RETURN(L) end;
for k from 0 to n do
if abs(c[k])<H/1000 then
Sim[k]:=ellipse([k,c[k]],0.2,0.01*H,color=blue);
else
Sim[k]:=plottools[arrow]([k,0],[k,c[k]],0.2,0.2,0,color=black);
end if;
end do;
convert(Sim,list);
end:
> Spectr:=proc(n,a,b,c,Risphi) local k,R,phi;
for k from 0 to n do
c[k]:=evalf(abs(I*a[k]+b[k])):
# print(k,c[k]);
phi:=evalf(argument(I*a[k]+b[k]));
R[k]:=[eval(k),eval(phi)];
end:;
Risphi:=plot(convert(R,list)):
end:
> T:=3;# величина периода
> F_for_all:=proc(t) global tau;fun(t/tau);end proc:;
> Ris1:=plot(F_for_all(t),t=0..T,color=brown,thickness=2,discont=true):display(Ris1);
> Fourier_T(F_for_all,0,T,k):
> a0:=evalf(A0);
> Nk:=40;
> Trig_polynom(Nk,T,A0,Ak,Bk,k):
> display(RisTrig,Ris1);
> Spectr(Nk,a,b,c,'Risphi1');
> display(ARR(Nk,c));
>
> Ampl:=display(ARR(Nk,c)):;
> 2: DTF:=proc (y,N,Y) local n,k,j,p,h;
n:=N-1;
h:=2*Pi/N;
2.1: for k from 0 to N do
p:=0;
for j from 0 to n do
p:=p+evalf(y[j]*exp(-I*k*j*h));
end;
Y[k]:=evalf(1/N*p);
end:
end:;
> 3: CDTF:=proc(N,Y,y) local n,k,h,p,j;
n:=N-1;
h:=2*Pi/N;
for k from 0 to n do
p:=0;
for j from 0 to n do
p:=p+Y[j]*exp(I*k*j*h);
end;
y[k]:=evalf(Re(p));
end:
end:
> Setka_DTF:=proc(Nt,T,F,Y::array) local h,j,x,R,RL;
global GrafF;
h:=T/Nt;
for j from 0 to Nt do
x:= evalf(j*h);
Y[j]:= F(x);
R[j]:=[j,eval(Y[j])];
end:
5.1: R[Nt]:=[j,eval(Y[0])];
GrafF:=plot(convert(R,list),0..Nt-1,color=brown,
style=point,symbol=circle):
end:
> Spectr_DTF:=proc(n,C,A,phi) local k,R;global Risphi;
6.1: for k from 0 to n do
A[k]:=evalf(abs(C[k])):
phi[k]:=evalf(argument(C[k]));
R[k]:=[eval(k),eval(phi[k])];
end:;
Risphi:=plot(convert(R,list),thickness=2,color=blue,style=point,symbol=box):
end:
Параметры задачи
> Nt:=33:`число дискретных отсчетов `:
> n:=Nt;N:=Nt-1;# параметры ДПФ
> C:=array(0..n):phi:=array(0..n):A:=array(0..n):;
Y:=array(0..N):
> Setka_DTF(N,T,F_for_all,Y);
> DTF(Y,Nt,C,n):
> Spectr_DTF(n,C,A,phi):
Для четных N
> display(ARR(n,A));
> display(ARR((n-1)/2,A));
> CDTF(Nt,C,F):
> display(GrafF,ARR(n-1,F));
> Setka:=proc(Nt,T,F,Y::array) local h,j,x,R,RL;
global GrafF;
h:=T/Nt;
for j from 0 to Nt do
x:= evalf(j*h);
Y[j]:= F(x);
R[j]:=[x,eval(Y[j])];
end:
5.1: R[Nt]:=[x,eval(Y[0])];
GrafF:=plot(convert(R,list),0..T,color=brown,
style=point,symbol=circle):
end:
> F_Discret:=proc (Y,N,a,b,n) local k,j,p,q,h;
h:=2*Pi/N;
for k from 0 to n do
p:=0;q:=0;
for j from 0 to N do
p:=p+evalf(Y[j]*cos(k*j*h));
q:=q+evalf(Y[j]*sin(k*j*h));
end;
a[k]:=2/N*p;b[k]:=2/N*q;
# print(k,a[k],b[k]);
end:
if 2*n=N then b[n]:=0; end;
RETURN(n);
end:;
> 3: Trig:=proc(t,n,T,a,b) local z,k;
z:=a[0]/2+sum(a[k]*cos(k*t*2*Pi/T)+b[k]*sin(k*t*2*Pi/T),k=1..n);
end:
> M:=Nk:
> a:=array(0..M):b:=array(0..M):c:=array(0..M);
> Setka(N,T,F_for_all,Y):
> F_Discret(Y,N,a,b,M):
> Cl:=blue,red,brown:;
> 15: RT:=seq(plot(Trig(t,8*k,T,a,b),t=-0.1..T+0.1,
numpoints=500,color=Cl[k]),k=1..3):
> 16: display(RT,GrafF);
> Spectr(M,a,b,c,'Grafphi');:
> display(ARR(M,c));:
> display(Ampl);
unit koren1;
interface
uses
Windows, Messages, SysUtils, Variants, Classes, Graphics, Controls, Forms,
Dialogs, StdCtrls, ExtCtrls;
type
TForm1 = class(TForm)
Button1: TButton;
Label1: TLabel;
Button2: TButton;
Label2: TLabel;
Button3: TButton;
ListBox1: TListBox;
Button4: TButton;
Button5: TButton;
Image1: TImage;
procedure Button1Click(Sender: TObject);
procedure Button2Click(Sender: TObject);
procedure Button3Click(Sender: TObject);
procedure Button4Click(Sender: TObject);
procedure FormCreate(Sender: TObject);
procedure Button5Click(Sender: TObject);
private
{ Private declarations }
public
{ Public declarations }
end;
var
Form1: TForm1;
implementation
{$R *.dfm}
{============================}
const
N=30;
m=N div 2;
x0=0;
T=3;
hx=T/N;
type
koeff=array[0..m] of real;
dann=array[0..N] of real;
var
Tau:real;
Y:dann;
a,b:koeff;
h:real;
eps:real;
La:real;
Rb:real;
Nkoeff:integer;
function pol(t:real):real;
begin
pol:=sqr(sqr(t))-sqr(t)*t-2;
end;
function derive(t:real):real;
begin
derive:=sqr(t)*t*4-sqr(t)*3;
end;
function root(a,b:real):real;
var
X0,X1,delta:real;
begin
X0:=(a+b)/2;
repeat
X1:=X0-pol(X0)/derive(X0);
delta:=abs(X1-X0);
X0:=X1;
until delta<0.00001;
root:=X0;
end;
function riter(a,b:real):real;
var
X0,X1,delta:real;
const
lambda=0.0001;
begin
X0:=(a+b)/2;
repeat
X1:=X0-pol(X0)*lambda;
delta:=abs(X1-X0);
X0:=X1;
until delta<0.00001;
riter:=X0;
end;
function signal(t:real):real;
var
z:real;
begin
if t<0 then
z:=0
else
if t<=1/2 then
z:=cos(Pi*t)
else
if t<=1 then
z:=0
else z:=0;
signal:=z;
end;
procedure Trig(m,N:integer;Y:dann;var a,b:koeff);
var
j,k:integer;
p,q:real;
x:real;
h:real;
begin
h:=2*Pi/N;
for k := 0 to m do
begin
p:=0;q:=0;
for j := 1 to N do
begin
x:=j*h;
p:=p+Y[j]*cos(x*k);
q:=p+Y[j]*sin(x*k);
end;
a[k]:=p*2/N;
b[k]:=q*2/N;
end;
end;
function Tpol(m:integer;x:real):real;
var
z:real;
k:integer;
begin
z:=a[0]/2;
for k:=1 to m do
z:=z+(a[k]*cos(k*2*Pi/T*x)+b[k]*sin(k*2*Pi/T*x));
Tpol:=z;
end;
procedure grafik(numvar:integer);
type
dann= array[0..N] of real;
var
L,R,W,H: integer;
X: dann;
Y: dann;
k:integer;
ymin,ymax:real;
Mx,My:real;
x0,y0: integer;
posx,posy:integer;
Nkf:string;
tx:real;
ypol:real;
procedure min_max(N:integer;Y:dann; var min, max:real);
var
k: integer;
begin
min:=Y[0];max:=Y[0];
for k := 1 to N do
if Y[k]> max then
max:=Y[k]
else if Y[k]< min then
min:=Y[k];
{увеличим диапазон}
max:=max+0.1;
min:=min-0.1;
end;
begin
L:=20;
R:=form1.image1.clientHeight-20;
W:=form1.image1.Width-50;
H:=form1.image1.clientheight-50;
case numvar of
1: begin
for k:=0 to N do
X[k]:=signal(hx*k/Tau);
min_max(N,X,ymin,ymax);
Mx:=W/N;
My:=H/(ymax-ymin);
x0:=L;
y0:=R-abs(Round(ymin*My));
with form1.image1.Canvas do
begin
pen.Color:=clblue;
font.Name:='Tahoma';
font.Size:=8;
font.Color:=claqua;
for k:=0 to N do
begin
posx:=x0+round(k*Mx);
posy:=y0-round(X[k]*My);
textout(posx-2,posy-8,'o');
Pixels[posx,posy]:=clRed;
end;
pen.Width:=2;
Moveto(L,R);lineto(L,R-H);
moveto(x0,y0);lineto(x0+W,y0);
font.Color:=clred;
textout(x0+W,y0+10,'x');
textout(x0+W,y0-20,floattostrF(T,ffFixed,3,0));
textout(x0+round(W*Tau/T), y0-20,'tau='+ floattostrF (Tau,ffFixed, 6, 3));
Nkf:=Inputbox('Число коэффициентов полинома','например 10','20');
Nkoeff:=strtoint(Nkf);
pen.Color:=clNavy;
tx:=0;
ypol:=Tpol(Nkoeff,tx/Tau);
posx:=x0+round(0*Mx/2);
posy:=y0-round(ypol*My);
moveto(posx,posy);
for k:=1 to 2*N do
begin
tx:=hx*k/2;
ypol:=Tpol(Nkoeff,tx/Tau);
posx:=x0+round(k*Mx/2);
posy:=y0-round(ypol*My);
lineto(posx,posy);
end;
end;
end;
2: begin
for k:=0 to m do
Y[k]:=sqrt(sqr(a[k])+sqr(b[k]));
min_max(m,Y,ymin,ymax);
Mx:=W/m;
My:=H/(ymax-ymin);
x0:=L;
y0:=R-abs(Round(ymin*My));
with form1.image1.Canvas do
begin
pen.Width:=2;
pen.Color:=clred;
Moveto(L,R);lineto(L,R-H);
moveto(x0,y0);lineto(x0+W,y0);
pen.Width:=5;
pen.Color:=clblue;
for k:=0 to m do
begin
posx:=x0+round(k*Mx);
posy:=y0-round(Y[k]*My);
moveto(posx,y0);
lineto(posx,posy);
end;
end;
end;
end;
end;
{===============================}
procedure TForm1.Button1Click(Sender: TObject);
const
a=0; b=3;
begin
Tau:=root(a,b);
Label1.Caption:='корень='+floattostr(Tau);
end;
procedure TForm1.Button2Click(Sender: TObject);
const
a=0; b=3;
begin
Tau:=riter(a,b);
Label2.Caption:='корень='+floattostr(Tau);
Button2.Visible:=false;
Button3.Visible:=true;
Button3.SetFocus;
Button4.Visible:=false;
end;
procedure TForm1.Button3Click(Sender: TObject);
var
j: Integer;
s1,s2,s3:string;
begin
for j := 0 to N do
Y[j]:=signal(x0+j*hx);
Y[N]:=(Y[0]+Y[N])/2;
Trig(m,N,Y,a,b);
for j := 0 to m do
begin
str(j:2,s1);
str(a[j]:10:5,s2);
str(b[j]:10:5,s3);
listbox1.Items.Add(s1+s2+s3);
end;
label1.Caption:='Вычислены коэффициенты Фурье';
Button3.Visible:=false;
Button4.Visible:=true;
Button4.SetFocus;
end;
procedure TForm1.Button4Click(Sender: TObject);
begin
Form1.Caption:='Построение графиков';
{здесь поместим алгоритм построения графика сигнала и триг. полинома}
Label1.Caption:='графики сигнала и полинома';
Button4.Visible:=false;
Button5.Visible:=true;
Button5.SetFocus;
Label2.Visible:=false;
{Label3.Visible:=false;
Label4.Visible:=false;}
listbox1.Visible:=false;
grafik(1); { 1 - сигнал точками, полином - линией}
Label1.Caption:='графики сигнала и полинома';
Button4.Visible:=false
end;
procedure TForm1.Button5Click(Sender: TObject);
begin
Form1.Caption:='Спектр амплитуд';
{здесь поместим алгоритм построения спектра амплитуд}
Label1.Caption:='График спектра амплитуд';
Button5.Visible:=false;
form1.image1.Canvas.FillRect(rect(0,0,clientwidth,clientheight));
grafik(2);
end;
procedure TForm1.FormCreate(Sender: TObject);
begin
Button2.Visible:=true;
Button3.Visible:=false;
Button4.Visible:=false;
Button5.Visible:=false;
end;
end.
Часть 2.
> restart;
> with(linalg):with(plots):
pp:=(x,y)->[x,y];
> fun:= proc(t) local z ;z:=piecewise(t<0,0,t<1/2,cos(Pi*t),t<=1,0,0); evalf(z);end;
> plot(fun(t),t=-1..2,thickness=2,color=brown);
> p(x):=x^4-x^3-2;
> Koeff:=fsolve(p(x),x,0..2);
> T:=Koeff;
> tau:=1;
> Period:=proc(t,t0,tau,T,f) local x,z;
x:=evalf(t-t0-floor((t-t0)/T)*T);
z:=fun(x/tau);evalf(z);
end;
> plot(Period(x,0,tau,T,fun),x=-1..3,thickness=2,color=brown):;
> #==============================================================================
>
> Koc:=0.2;Nzac:=4;
> ur:=diff(U(t),t);
> F:=Nzac*(cos((4+Nzac/10)*t+U(t))+Koeff*Period(t,0,tau,T,f)-Koc*U(t));
> RK:=dsolve({ur=F,U(0)=0.15},U(t),type=numeric,output=listprocedure);
> fU:=subs(RK,U(t));
> T0:=5;Nt:=50;h:=T0/Nt;
> Tx:=array(0..Nt):U:=array(0..Nt):U_map:=array(0..Nt);
> for j from 0 to Nt do
x:=j*h;z:=fU(x);Tx[j]:=x;U[j]:=z;U_map[j]:=z;
#print(x,z);
od:
> RisU:=zip(pp,Tx,U):
> RU:=plot(RisU):
> display(RU):;
> #====================================
> RisU:=zip(pp,Tx,U):
> whattype([RisU]);
> RU0:=plot(RisU,style=point,symbol=cross):
> display(RU0):;
Определение периода с помощью функции автокорреляции
> R:=proc(Nt,T::array,U::array)
local k,j,t,z,Af::array,Uf::array,RAf,RisAf;
Af:=array(0..Nt);Uf:=array(0..2*Nt);
for j from 0 to Nt do Uf[j]:=U[j];Uf[j+Nt]:=U[j];
end;
for j from 0 to Nt do
t:=0;
for k from 0 to Nt do
t:=t+evalf(Uf[k]*Uf[k+j]);
end;
Af[j]:=evalf(t/Nt);
end;
RAf:=zip(pp,T,Af);RisAf:=plot(RAf):
display(RisAf);
end:
> R(Nt,Tx,U);
>
+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++=
>
> fn:=`E:\\work\\Саляхутдинов.txt`;
>
>
++++++++++++++++++++++++++++++++++++++++++++++++++++
> L:=readdata(fn,2):;
Nstrok:=vectdim(L);
> U_n:=array(1..Nstrok);:
T_n:=array(1..Nstrok);
> for j from 1 to Nstrok do
T_n[j]:=L[j,1];
U_n[j]:=L[j,2];
#print(j,T_n[j],U_n[j]);
od:
> u1:=zip(pp,T_n,U_n):
> RU1:=plot(u1,style=point,symbol=cross,color=black):
> display(RU,RU1);
>
> #printf("%s",` № t U_map U_pas разн \n`);
for k from 0 to Nt do t:=Tx[k]:del:=U_map[k]-U_n[k+1];
#printf("% 3.0f % 6.2f % 8.4f % 8.4f % 8.4f \n",k,t,U_map[k],U_n[k+1],del):
end:;
unit final;
interface
uses
Windows, Messages, SysUtils, Variants, Classes, Graphics, Controls, Forms,
Dialogs, Menus, StdCtrls;
type
TForm1 = class(TForm)
MainMenu1: TMainMenu;
N1: TMenuItem;
N2: TMenuItem;
N3: TMenuItem;
SaveDialog1: TSaveDialog;
Label1: TLabel;
Label2: TLabel;
ListBox1: TListBox;
N4: TMenuItem;
N5: TMenuItem;
N6: TMenuItem;
Label3: TLabel;
Label4: TLabel;
procedure N3Click(Sender: TObject);
procedure N2Click(Sender: TObject);
procedure N5Click(Sender: TObject);
procedure N6Click(Sender: TObject);
private
{ Private declarations }
public
{ Public declarations }
end;
var
Form1: TForm1;
implementation
{$R *.dfm}
{=====================================================}
const
U0=0.15;
nzac=4;
koc=0.2;
hintegr=0.1;
Npoint=50;
var
tau:real; {период сигнала}
Koeff:real; {равен периоду}
min,sec:byte;
frez:string;
ftxt:text;
{===================корень полинома==============}
{фунция, задающая вычисление полинома в точке}
function polynom(t:real):real;
begin
polynom:=sqr(sqr(t))-sqr(t)*t-2;
end;
{процедуры метода простой итерации}
function derive(x:real):real;
begin
derive:=sqr(x)*x*4-sqr(x)*3;
end;
procedure Iter(a,b:real;var root:real;var K:integer);
const
eps=0.0000001;{погрешность определения корня}
Q0=1.5; {нормирующий множитель для ламбда}
var
lambda:real;
g:real;
x,x0:real;
begin
g:=derive(b);
form1.label2.caption:='max производной = '+ floattostrF (g, ffGeneral, 5, 1);
lambda:=1.99/g;
x:=(a+b)/2;
k:=0;
repeat
x0:=x;
x:=x0-lambda*polynom(x0);
form1. listbox1.Items.Add(floattostrF(x,ffGeneral, 12, 9));
k:=k+1;
until abs(x-x0)<eps;
root:=x;
end;
function signal(t:real):real;{сигнал 11-го варианта}
var
z:real;
begin
if t<0 then
z:=0
else
if t<1/2 then
z:=cos(Pi*t)
else
if t<=1 then
z:=0
else
z:=0;
signal:=z;
end;
function Period(x,T:real):real;
var
z:real;
begin
z:=x-trunc(x/T)*T;{выделение дробной части}
Period:=signal(z);
end;
{====================правая часть диф. уравнения===}
function F(t,U:real):real;
begin
F:=Nzac* (cos((4+nzac/10)*t+U)+Koeff*Period(t,Koeff)-Koc*U);
end;
procedure difur;
var
j:integer;
U,t:real;
U1:real;
strU,strt:string;
begin
form1.listbox1.Clear;
U:=U0;
t:=0;
Koeff:=Tau;
for j :=0 to Npoint+1 do
begin
strt:=floattostrF(t,ffGeneral, 7, 4);
strU:=floattostrF(U,ffGeneral, 7, 4);
form1. listbox1.Items.Add(strt+' '+ stru);
writeln(ftxt,t:8:2,' ',U:8:3);
{Метод Эйлера}
U1:=U+hintegr*F(t,U);
U:=U1;
t:=t+hintegr;
end;
end;
procedure TForm1.N2Click(Sender: TObject);
begin
if savedialog1.Execute then
begin
frez:=savedialog1.FileName;
label1.caption:='эапись в файл '+#10+frez;
end;
assignfile(ftxt,frez);
rewrite(ftxt);
end;
procedure TForm1.N3Click(Sender: TObject);
begin
Form1.Close;
Closefile(ftxt);
end;
procedure TForm1.N5Click(Sender: TObject);
var
a,b:real;
kiter:integer;
begin
Form1.Caption:='Вычисляем корень полинома';
b:=3;
a:=0;
{здесь вызов алгоритма вычисления корня уравнения}
Iter(a,b,Tau,Kiter);{вызов процедуры метода итераций}
Label3.Caption:='корень равен '+floattostr(Tau);
Label4.Caption:='число итераций '+inttostr(Kiter);
end;
procedure TForm1.N6Click(Sender: TObject);
begin
difur;
end;
end.