i have results given by instrutor ,however i need MATLAB code for this
Task #1: estimation of correlation features, PSD features, and relationships for LTI methods inputs and outputs
Be aware: this doc makes use of the identical notation as within the course slides for Part 2, e.g., mx , mˆx , fyx( )l , fˆyx( )l , pyx( )l , Fyx(ej?), Pyx(ej?), and PyxW (ej?).
Be aware: equations in purple correspond to extra info which doesn’t seem within the slides.
All through the task, we use the next relationships describing two LTI methods:
, with impulse response g n( ) 5 (= dn- +1) four (dn-2)
ok=-Eight
n
, with impulse response hn( )= ? ?? ?1 un( ).
ok=-Eight ? ?2
Every sign x n( ), d n( ), y n( ) is an remark (sequence, realization) of L samples from a large sense stationary, imply ergodic and correlation ergodic random course of X .
The sign x n( ) has an autocorrelation perform of: f sdxx( )l = x2 ( )l + mx 2 with sx2 =1 and mx = 1 four (though generally we are going to usemx = zero).
Question Assignment 1
a) Discover the theoretical expressions for fyx()l ,fdx()l ,fyy()l .
b) Use the next code to generate an occasion of the sequences x n( ), d n( ), y n( ) with L samples:
L=2^16; % nb. of samples for every realization
g=[0,5,4]; % FIR filter coeffs hb=[1]; % IIR filter numerator ha=[1,-1/2]; % IIR filter denominator
extra_samples_to_remove_transients = max(size(g),5/(1-abs(ha(2))))-1;
sigma_x_sq=1; % variance of enter x mx=1/four; % imply of enter x
x=sqrt(sigma_x_sq)*randn(1,L+extra_samples_to_remove_transients)+mx; d=filter(g,1,x); y=filter(hb,ha,x);
% take away transients
x=x(extra_samples_to_remove_transients+1:finish); d=d(extra_samples_to_remove_transients+1:finish); y=y(extra_samples_to_remove_transients+1:finish);
Compute the unbiased and biased cross-correlation estimates fˆyx( )l and pyx( )l and plot them on the identical determine because the theoretical fyx()l (e.g., three subplots). Produce two figures, utilizing the 2 following scales for the “lag” axis: -10= =l 10 and – +L 1= =l L -1.
Be aware: for all of the plots within the task, you need to embody significant titles, label for horizontal axis, and a grid.
Be aware: in this task you should use both the plot() or stem() features, relying on which one produces the clearest plot for every case.
Be aware: If you happen to use features akin to xcorr() as an alternative of programming the small print of the estimated correlation perform, just remember to use it with the proper parameters such that it offers the identical results as our definitions (i.e., identical scale, identical normalization, and so on.).
Repeat for fˆdx()l ,pdx()l ,fdx()l , and fˆyy( )l ,pyy()l ,fyy()l .
c) Producing N =100 sequences of y n( ), you’ll be able to experimentally approximate the variance of the unbiased estimate fˆyy( )l with the next expression:
the place fˆyyi ( )l represents the autocorrelation estimate measured from the i th sequence of y n( ).
Likewise, the variance of the biased autocorrelation estimate pyy()l might be experimentally measured, contemplating the recognized bias of this estimate:
?? ? l? ?2?2
var = E ??(pyy()l – E p?? yy()l ??)2??? = E ??????pyy()l -???1- L ???fyy()l ??? ???? ˜ N1 ?iN=1 ????pyyi ()l -????1- Ll ????fyy()l ????
?
?
In precept, the vary – +L 1= =l L -1 ought to be used for the plots, however you’ll need to discard the bigger values of l close to – +L 1 and L -1 (e.g. minimize 1% on either side) as a result of the measured variance at these lags turns into too giant and prevents seeing the remainder of the curves.
d) Additionally it is potential to check the measured varfˆyy()l and var (earlier Question Assignment) with corresponding theoretical values:
varth fˆyy()l ˜ (L Ll )2 ??? n=- +?L-L1 1fyy( )n 2 + n=- + +L?- -L1 1ll f*yy(n – l) (fyy n + l)???? = (L -Ll )2fˆtmp()l
– ?
varth ˜ L1 ??? n=- +?L-L1 1fyy( )n 2 + n=- + +L?- -L1 1l l f*yy(n – l) (fyy n + l)???? = L1fˆtmp()l – +L 1= =l L -1
?
For instance a code like this might be included in your code to compute such features: vartmp=ones(1,2*L-1)*sum(abs(phi_yy.*phi_yy)); for lag=zero:L-1
vartmp(L+lag)=vartmp(L+lag)+sum(phi_yy(1:end-2*lag).*phi_yy(1+2*lag:finish)); finish
vartmp(1:L-1)=flip(vartmp(L+1:finish));
var_phi_yy_hat_th=vartmp.*(L./(L-abs(-L+1:1:L-1)).^2); var_p_yy_th=vartmp/L;
Nevertheless, these equations are solely legitimate for mx = zero (or they need to be utilized for autocovariance features as an alternative of autocorrelation features), and for pretty small absolute lag l values. So right here you need to repeat the simulations of c) however utilizing mx =zero
. Plot varfˆyy()l ,var,varth fˆyy()l ,varth on the identical determine (4 subplots) for the next ranges of lags (i.e. two figures):
• Utilizing 15% of lag values (i.e., smaller absolute values round l =zero); • From -200= =l 200.
Question Assignment 2
Be aware: for this Question Assignment on PSD estimation, you need to not use Matlab/Octave features akin to pwelch(), periodogram(), or cpsd() to compute the PSD portions. To raised perceive the computations, use the fundamental fft() perform and the equations for the PSD estimates. Be aware: frequency area alerts in this Question Assignment are written with the discrete time Fourier rework (DTFT) notation (e.g. Fxx(ej?) -p ? p= ), however it’s understood that in apply they’re discretized/computed with quick Fourier transforms (FFTs):
Fxx( )ok zero = =ok Nfft -1, ? pk = k2Nfft , the place Nfft is the dimensions of the FFT (after zeropadding if zero-padding is used). FFT results like Fxx( )ok might be plotted to approximate DTFTs, with correct scaling of frequency axis and “wrap-around” or “fftshift” processing in order that frequencies p ? p= 2 are proven as -p ?= zero.
a) Discover the theoretical expressions for Fyx(ej?), Fyy(ej?), for the circumstances with mx =1 four and mx =zero.
b) For the circumstances with mx = 1 four and mx = zero, generate x n( ) and y n( ) sequences of L = 216 samples, compute the periodogram estimates Pyx(ej?)and Pyy(ej?), and plot Pyx(ej?) ,
Pyy(ej?), Fyx(ej?) , Fyy(ej?) in the identical determine (4 subplots). Be aware that for mx = 1 four the features Fyx(ej?), Fyy(ej?) have infinite Dirac elements which can’t be plotted
(simply ignore them for the plots).
You will note that for the case mx = 1 four the DC element (frequency zero) in Pyx(ej?) and
Pyy(ej?) fully hides the remainder of the plot. However when you zoom the results with mx = 1 four in order that small amplitudes are displayed (use one other determine for that), you’ll then see that besides for the DC element the results are the identical as for mx =zero.
Due to this fact, it’s usually handy in apply to take away the DC element of a sign earlier than computing its Fourier rework or PSD estimate. For the remainder of this task, generate alerts with mx = zero.
c) From the earlier half b), you need to discover that the periodogram estimates Pyx(ej?), Pyy(ej?) have a variety of fluctuations: this is the impact of the big variance of those estimators, which doesn’t decay with the size L of the noticed alerts. In this half, you might be requested to confirm this truth.
Producing N =100 sequences of y n( ), you’ll be able to experimentally measure the variance of the periodogram estimate Pyy(ej?) with the next expression: var ( ) ( ) ( ) 2 ,
=1
the place mˆP = N1 ?iN=1 P eyyi ( j?), and Pyyi (ej?) represents the periodogram estimate measured from the i th sequence of y n( ). Utilizing completely different values of L (ex. L = 212 and L = 216 ), confirm that the variance of the periodogram estimate Pyy(ej?) doesn’t considerably decay because the variety of noticed samples L is elevated.
Additionally evaluate (by plotting on the identical determine) with the theoretical variance of the periodogram estimate:
Var P yy(ej?) = E P[( yy(ej?) – E P[ yy(ej?)]) ]2 ˜Fyy2 (ej?).
d) With L = 216, M = 128, Okay = +1 ??(L M M- )/( / 2)?? = 1023 (case with 50% M /2 samples overlap between home windows), and utilizing a (normalized) Hamming window, discover the Welch estimate PyyW (ej?). Plot PyyW (ej?) and the theoretical Fyy(ej?) on the identical determine (two subplots).
e) Producing N =100 sequences of y n( ) with the identical parameters L = 216, M = 128, Okay = +1 ??(L – M)/(M / 2)?? = 1023 , you’ll be able to experimentally measure the variance of the Welch estimate PyyW (ej?) with the next expression: var ( ) ( ) ( ) 2 ,
=1
the place mˆP = N1 ?iN=1 PyyW i, (ej?), and PyyW i, (ej?) represents the Welch estimate measured from the i th sequence of y n( ).
Additionally evaluate (by plotting on the identical determine) with the theoretical variance of the Welch estimate for 50% overlap:
Var P e yyW ( j?) ˜ 89Okay Fyy2 (ej?).
Be aware the distinction between the variance of the Welch estimate var and the earlier variance of the periodogram estimate var.
f) With L = 216, M = 128, Okay = +1 2(L M M- )/ = 1023 (case with 50% overlap between home windows), and utilizing a (normalized) Hamming window, compute an estimate PyxW (ej?) of the cross-PSD between x n( ) and y n( ), compute estimates of the auto-PSDs PxxW (ej?),
PyyW (ej?), and use them to compute an estimated coherence perform
?Wyy(ej?)= PyxW (ej?) PxxW (ej?)PyyW (ej?). Examine ?Wyy(ej?) with the theoretical coherence magnitude ?yx(ej?) by plotting them on the identical determine.

Published by
Medical
View all posts