Graphics Programs Reference
In-Depth Information
The phase shift at a frequency of
f
=0.2 (period
τ
=5) can be interpolated from
the phase spectrum
interp1(f,phase,0.2)
which produces the output
ans =
1.2568
The phase spectrum is normalized to one full period
=2, therefore a phase
shift of 1.2568 equals (1.2568*5)/(2*) = 1.0001, which is the phase shift of
one that we introduced at the beginning.
We now use two sine waves with different periodicities to illustrate
crossspectral analysis. The both have a periodicity of 5, but with a phase
shift of 1, then they have both one other period, which are different, how-
ever.
τ
clear
t = 0.01 : 0.1 : 1000;
y1 = sin(2*pi*t/15) + 0.5*sin(2*pi*t/5);
y2 = 2*sin(2*pi*t/50) + 0.5*sin(2*pi*t/5+2*pi/5);
plot(t,y1,'b-',t,y2,'r-')
Now we compute the crossspectrum, which clearly shows the common pe-
riod of
τ
=5 or frequency of
f
=0.2.
[Pxy,f] = cpsd(y1,y2,[],0,512,10);
magnitude = abs(Pxy);
plot(f,magnitude);
xlabel('Frequency')
ylabel('Power')
title('Cross PSD Estimate via Welch')
The coherence shows a large value of approximately one at
f
=0.2.
[Cxy,f] = mscohere(y1,y2,[],0,512,10);
plot(f,Cxy)
xlabel('Frequency')
ylabel('Magnitude Squared Coherence')
title('Coherence Estimate via Welch')
The complex part is required for calculating the phase shift between the two
sine waves.
