Lucrarea nr. 9

1

Sisteme discrete în timp (sisteme cu eşantionare)

1 Aspecte teoretice

Teorema eşantionării (Teorema lui Shannon) :

Metode de discretizare: Elementul de eşantionare reţinere de ordinul zero (zoh)

2 Scopul lucrării

Obţinerea sistemelor cu eşantionarea pornind de la sisteme continue în timp; Discretizarea funcţiilor de transfer a sisteme continue utilizând mediul Matlab; Simularea răspunsului în timp al sistemelor discrete: răspunsul la impuls, treaptă, rampă și sinus

3 Alegerea perioadei de eşantionare

În cazul în care nu se respectă teorema eșantionării apare fenomenul de frecvenţe false (aliasing) sau fenomenul de îndoire/înjumătăţire al frecvenţei (folding effect). În figura următoare este evidenţiat fenomenul de apariţiei al unui ‚ alias’ sau frecvenţă falsă.

Figura 1 Fenomenul de „aliasing”

În cazul sistemelor de control cu eşantionare sau digitale, interesul este pe stabilirea perioadei de eşantionare (care va apare ca parametru în structura modelului discret, a funcţiei de transfer discrete).

0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1-1

-0.8

-0.6

-0.4

-0.2

0

0.2

0.4

0.6

0.8

1x(t)=cos(2ft)

timp (sec.)

Am

pli

tud

ine

f=1Hz

f=9Hz

f=11Hz

Lucrarea nr. 9

2

Pentru a detecta regimurile tranzitorii în mod corect, perioada de eşantionare se va alege şi în funcţie de constantele de timp ale sistemului / procesului nu numai în funcţie de spectrul semnalului de intrare.

Dispunând de funcţia de trasfer în care sunt evidenţiaţi termenii de ordinul I şi II de la numitor:

Perioda de eşantionare se va alege (conform teoremei eşantionării):

În practică, perioada de eşantionare se va alege de 10 ori maI mică decât cea mai mică constantă de timp a sistemului, respectiv cea mai mică perioadă a semnalului de intrare.

În cazul în care nu se respectă teorema eşantionării este posibil să nu se detecteze corect regimul de funcţionare (figura 2).

Figura 2 Alegere incorectă a perioadei de eşantionare

4 Reprezentarea semnalelor discrete în Matlab

Pentru a reprezenta grafic semnalele discrete în timp se sugerează utilizarea funcţiei stem din Matlab. Pentru punerea în evidenţă a semnalului eșantionat și reţinut (reconstruit) este utilă folosirea funcţiei

stairs.

Pentru familiarizarea cu cele două funcţii, să se realizeze în Matlab, un script care să genereze următorul grafic din Figura 3.

0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 50

0.2

0.4

0.6

0.8

1

1.2

1.4

1.6

1.8

Raspuns la treapta

timp (sec.) (seconds)

Am

plitu

din

e

el. de ord. II

el de ord. I

Lucrarea nr. 9

3

Figura 3 Funcţia stem şi stairs din Matlab

%% figura 3 f=1; % frecventa semnalului n=10; % numar de esantioane Te1=1/f/n ; % perioada de esantionare t=0:Te1:1; % timpul de simulare esantionat cu Te1 x=cos(2*pi*f*t); % functia cosinus de frecventa f subplot(3,2,[1,3,5]); stem(t,x);hold;stairs(t,x,'r');hold;% comparatie stem vs

stairs legend('stem','stairs'); title('stem vs. stairs'); ylabel('cos(2\pift)');xlabel('timp (sec.)')

Te2=1/f/20;% 20 de esantioane Te3=1/f/100;% 100 de esantioane t2=0:Te2:1;t3=0:Te3:1; x2=cos(2*pi*f*t2);% semnal cu 20 de esantioane pe o perioada x3=cos(2*pi*f*t3);% semnal cu 100 de esantioane pe o perioada

% se reprezinta grafic semnalele subplot(3,2,2);stem(t,x);title('10 esantioane'); subplot(3,2,4);stem(t2,x2);title('20 de esantioane'); subplot(3,2,6);stem(t3,x3);title('100 de esantioane');

0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1-1

-0.8

-0.6

-0.4

-0.2

0

0.2

0.4

0.6

0.8

1

timp (sec.)

Am

pli

tud

ine

stairs vs. stem

stairs

stem

0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1-1

-0.5

0

0.5

110 esantioane

0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1-1

-0.5

0

0.5

120 de esantioane

0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1-1

-0.5

0

0.5

1100 de esantioane

Lucrarea nr. 9

4

5 Discretizarea funcţiilor de transfer în Matlab

Există două funcţii disponibile pentru discretizarea în Matlab a unei funcţii de transfer: a) c2dm (continuous to discrete method, versiunea veche)

Sintaxa funcţiei este următoarea: [nz,dz]=c2dm(nc,dc,Te,’metoda’)

- nz și dz (parametrii de ieșire) reprezintă numărătorul și numitorul funcţiei de transfer în z, sub formă de vectori linie care conţin coeficienţii celor două polinoame în ordine descrescătoare a puterilor lui z;

- nc și dc sunt numărătorul respectiv numitorul funcţiei de transfer a sistemului continuu în timp, precizaţi sub formă de vectori linie;

- Te reprezintă perioada de eșantionare (în secunde); - ‘metoda’ este un string (şir de caractere) și referă metoda de discretizare cu următoarele valori posibile: o 'zoh' element de eșantionare reţinere de ordinul zero ( zero order hold).

o ‘foh' element de eșantionare reţinere de ordinul întâi ( first order hold)

o 'tustin' transformarea biliniară

o 'prewarp' transformarea biliniară, cu frecvenţa de tăiere precizată

{C2DM(num,den,Ts,'prewarp',wt) o 'matched' metoda potrivirii parametrilor.

b) funcţia c2d (continuous to discrete) Este varianta modificată a funcţiei c2dm, pentru a putea fi apelată cu obiect de tip funcţie de transfer.

Sintaxa funcţiei este următoarea: sist_discret=c2d(sist_continuu,Te,’metoda’)

- pe lângă metodele de discretizare de la c2dm, în acest caz mai este predefinită şi:

o ‘imp’ (metoda bazată pe invarianţa răspunsului la impuls)

6 Regimuri de funcţionare

Rezultatele teoretice obţinute în urma mapării semiplanului complex stâng ‘s’ în planul ’z’ prin transformarea , sunt următoarele:

- Axa reală negativă se regăsește pe axa reală pozitivă și subunitară (pentru polii pozitivi și subunitari se va obţine regim aperiodic amortizat);

- Interiorul fâșiei primare se regăsește în interiorul cercului de rază 1 (pentru poli complex conjugaţi cu modulul subunitar se va obţine regim oscilant amortizat)

- Axa imaginară se regăsește pe conturul cercului (pentru poli complex conjugaţi cu modulul egal cu 1, se va obţine regim oscilant întreţinut)

- Pentru poli situaţi în afara cercului de rază 1, sistemul este instabil sau pe răspuns la semnal treaptă se obţine un regim oscilant/aperiodic neamortizat.

Lucrarea nr. 9

5

7 Probleme rezolvate

Pentru sistem continuu determinat de funcţia de transfer

1. Să se obţină echivalentul discret al funcţiei de transfer utilizând metoda de discretizare ‚zoh’.

Hc=tf(500*[1 410],conv([1 210],[1 320])); Te=max([1/210,1/320])/10;% se stabileste perioada de

esantionare Hd=c2d(Hc,Te,'zoh')

2. Să se aducă funcţia de transfer discretă de la forma de raport de două polinoame precizate în ordine descrescătoare a puterilor lui z, la forma în care sunt evidenţiaţi polii, zerourile şi factorul de proporţionalitatea k (de remarcat valorile polilor după discretizare: reale, subunitare şi pozitive)

zpk(Hd)

3. Să se modifice formatul de reprezentare în format ‚long’ şi să se reevalueze coeficienţii funcţiei de transfer discrete.

format long [num,den]=tfdata(Hd,'v')

4. Să se compare răspunsul la treaptă al sistemului continuu cu cel al sistemului discret.

step(Hc,Hd); legend('Hc','Hd')

Lucrarea nr. 9

6

5. Să se reprezinte răspunsul la treaptă al sistemului discret prin valori discrete (nu prin trepte de valoare constantă pe durata unei perioade de eşantionare)

[y,t]=step(Hd); step(Hc);hold;stem(t,y,'r');hold legend('Hc','Hd');

6. Să se compare răspunsul la treaptă al sistemului discret cu coeficienţii rotunjiţii nerotunjiţi, rotunjiţi la primele 4 zecimale, respectiv la 2 zecimale.

Hd4=tf([0 0.2312 -0.1902],[ 1 -1.7634 0.7769],Te); Hd3=tf([0 0.231 -0.190],[ 1 -1.763 0.776],Te);

0 0.005 0.01 0.015 0.02 0.025 0.03 0.035 0.04 0.045 0.050

0.5

1

1.5

2

2.5

3

3.5

Step Response

Time (seconds)

Am

plitu

de

Hc

Hd

0 0.005 0.01 0.015 0.02 0.025 0.030

0.5

1

1.5

2

2.5

3

3.5

Step Response

Time (seconds)

Am

plitu

de

Hc

Hd

Lucrarea nr. 9

7

Hd2=tf([0 0.23 -0.19],[ 1 -1.76 0.77],Te); step(Hd,Hd4,Hd3,Hd2);shg legend('15 zecimale','4 zecimale','3 zecimale','2 zecimale')

7. Să se realizeze un script de tip ‚function’ cu numele ‚reprezentare’care să reprezinte grafic polii şi zerourile sistemului discret în raport cu cercul de raza 1; modul de apelare să fie de forma:

reprezentare(Hd), cu parametrul de intrare Hd,un obiect de tip funcţie de transfer discretă.

function []=reprezentare(Hd) x=0:0.01:2*pi; plot(sin(x),cos(x),'r');% cercul de raza unitate [num,den]=tfdata(Hd,'v'); poli=roots(den); zerouri=roots(num); hold; plot([0 0],[-1.2 1.2],'--k',[-1.2 1.2],[0 0],'--k') plot(real(poli),imag(poli),'x','MarkerSize',6); plot(real(zerouri),imag(zerouri),'o','MarkerSize',6); title('Planul z') axis([-1.2 1.2 -1.2 1.2]) hold;shg

0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.080

0.5

1

1.5

2

2.5

3

3.5

4

Step Response

Time (seconds)

Am

plitu

de

15 zecimale

4 zecimale

3 zecimale

2 zecimale

Lucrarea nr. 9

8

8. Să se compare echivalentul numeric al funcţiilor de transfer :

eşantionat cu

şi

eşantionat cu

clear all;clc Hc1=tf(500*[1 410],conv([1 210],[1 320]));Te1=1/2100; Hc2=tf(50*[1 41],conv([1 21],[1 32]));Te2=1/210;

Hd1=c2d(Hc1,Te1,'zoh') Hd2=c2d(Hc2,Te2,'zoh')

Se observă că ambele funcţii de transfer au aceeaşi coeficienţi la numărătorul şi numitorul funcţiei de transfer, deosebirea între cele două realizându-se la nivel de perioade de eşantionare

-1 -0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8 1

-1

-0.8

-0.6

-0.4

-0.2

0

0.2

0.4

0.6

0.8

1

Planul z

z1 = 0.90484

z2 = 0.85866

_z1 = 0.82258

Lucrarea nr. 9

9

8 Probleme propuse spre rezolvare

1. Să se realizeze un script în Matlab pentru evidenţierea regimurilor de funcţionare menţionate şi prezentate în tabelul următor:

1 ----------------- (z-0.95) (z-0.75) Sample time: 0.1 seconds Discrete-time zero/pole/gain model.

H = 1 ------------------- z^2 - 1.9 z + 0.925 Sample time: 0.01 seconds roots([1 -1.9 0.925]) ans = 0.9500 + 0.1500i 0.9500 - 0.1500i

0 5 100

10

20

30

40

50

60

70

80Regim aperiodic amortizat

Time (seconds)

Am

plit

ud

e

-1 -0.5 0 0.5 1

-1

-0.8

-0.6

-0.4

-0.2

0

0.2

0.4

0.6

0.8

1

Planul z

z1 = 0.95z2 = 0.75

-1 -0.5 0 0.5 1

-1

-0.8

-0.6

-0.4

-0.2

0

0.2

0.4

0.6

0.8

1

Planul z

z1 = 0.95+ 0.15iz2 = 0.95-0.15i

0 0.5 10

10

20

30

40

50

60Regim oscilant amortizat

Time (seconds)

Am

plit

ud

e

Lucrarea nr. 9

10

H = 1 ----------------- z^2 - 1.826 z + 1 Sample time: 0.01 seconds >> roots([1 -1.826 1]) ans = 0.9130 + 0.4080i 0.9130 - 0.4080i

1 ----- z - 1

Sample time: 0.01 seconds

1 ---------------------- z^2 + 0.7925 z + 1.064 Sample time: 0.01 seconds roots([1 0.7925 1.064]) ans = -0.3962 + 0.9524i -0.3962 - 0.9524i

-1 -0.5 0 0.5 1

-1

-0.8

-0.6

-0.4

-0.2

0

0.2

0.4

0.6

0.8

1

Planul z

z1 = 0.91295+ 0.40808iz2 = 0.91295-0.40808i

0 0.1 0.2 0.3 0.4-2

0

2

4

6

8

10

12Regim oscilant intretinut

Time (seconds)

Am

plit

ud

e

-1 -0.5 0 0.5 1

-1

-0.8

-0.6

-0.4

-0.2

0

0.2

0.4

0.6

0.8

1

Planul z

z1 = 1

0 0.05 0.10

1

2

3

4

5

6

7

8

9

10Limita de stabilitate

Time (seconds)

Am

plit

ud

e

0 0.5 1 1.5 2-300

-200

-100

0

100

200

300Regim oscilant neamortizat

Time (seconds)

Am

plit

ud

e

-1 -0.5 0 0.5 1

-1

-0.8

-0.6

-0.4

-0.2

0

0.2

0.4

0.6

0.8

1

Planul z

z1 = -0.39625+ 0.95241iz2 = -0.39625-0.95241i

Lucrarea nr. 9

11

9 Bibliografie

[1]. Hăngănuţ Marius. Teoria Sistemelor. vol. II, Ed. Univ. Tehnice din Cluj-Napoca (1996). [2]. Kuo, Benjamin C. Automatic control systems. Prentice Hall PTR, 1981. [3]. Leonard, Naomi Ehrich, and William S. Levine. Using MATLAB to analyze and design Control Systems.

Benjamin-Cummings Publishing Co., Inc., 1995. [4]. Nise, Norman S. Control systems engineering, (With CD). John Wiley & Sons, 2007. [5]. Ogata, Katsuhiko, and Yanjuan Yang. Modern control engineering. Prentice Hall, 1970.