Click here for help on the
fft command and the
random command.
The periodogram is a simple and common method for estimating a
power spectrum.
Given a finite duration discrete random sequence x(n)x(n),
(n=0,1,⋯,N-1)(n=0,1,⋯,N-1), the periodogram is defined as
P
x
x
(
ω
)
=
1
N
X
(
ω
)
2
=
1
N
∑
n
=
0
N
-
1
x
(
n
)
e
-
j
ω
n
2
P
x
x
(
ω
)
=
1
N
X
(
ω
)
2
=
1
N
∑
n
=
0
N
-
1
x
(
n
)
e
-
j
ω
n
2
(2)where X(ω)X(ω) is the Discrete Time Fourier Transform (DTFT)
of x(n)x(n).
The periodogram Pxx(ω)Pxx(ω) can be computed using the
Discrete Fourier Transformation (DFT), which in turn can be
efficiently computed by
the Fast Fourier Transformation (FFT). If x(n)x(n) is of length NN,
you can compute an NN-point DFT.
P
x
x
(
ω
k
)
=
1
N
∑
n
=
0
N
-
1
x
(
n
)
e
-
j
ω
k
n
2
P
x
x
(
ω
k
)
=
1
N
∑
n
=
0
N
-
1
x
(
n
)
e
-
j
ω
k
n
2
(3)
ω
k
=
2
π
k
N
,
k
=
0
,
1
,
⋯
,
N
-
1
ω
k
=
2
π
k
N
,
k
=
0
,
1
,
⋯
,
N
-
1
(4)Using Equation 3 and Equation 4,
write a Matlab function called Pgram to calculate
the periodogram. The syntax for this function should be
[P,w] = Pgram(x)
where xx is a discrete random sequence of length NN. The
outputs of this command are PP, the samples of the periodogram,
and ww, the corresponding frequencies of the samples. Both
PP and ww should be vectors of length NN.
Now, let xx be a Gaussian (Normal) random variable with mean 0 and
variance 1. Use Matlab function random or randn to
generate 1024 i.i.d. samples of xx, and denote them as
x0x0, x1x1, ..., x1023x1023. Then filter the samples of xx with
the filter which
obeys the following difference equation
y
(
n
)
=
0
.
9
y
(
n
-
1
)
+
0
.
3
x
(
n
)
+
0
.
24
x
(
n
-
1
)
.
y
(
n
)
=
0
.
9
y
(
n
-
1
)
+
0
.
3
x
(
n
)
+
0
.
24
x
(
n
-
1
)
.
(5)Denote the output of the filter as y0,y1,⋯,y1023y0,y1,⋯,y1023.
The samples of x(n)x(n) and y(n)y(n) will be
used in the following sections. Please save them.
Use your Pgram function
to estimate the power spectrum
of y(n)y(n), Pyy(ωk)Pyy(ωk).
Plot Pyy(ωk)Pyy(ωk) vs. ωkωk.
Next, estimate the power spectrum of yy,
Pyy(ωk)Pyy(ωk) using 1/41/4 of the samples of yy. Do this
only using samples y0y0, y1y1, ⋯⋯, y255y255.
Plot Pyy(ωk)Pyy(ωk) vs. ωkωk.
-
Hand in your labeled plots and your
Pgram code.
- Compare the two plots. The first plot
uses 4 times as many samples as the second one. Is the first
one a better estimation than the second one? Does the first
give you a smoother estimation?
- Judging from the results, when the
number of samples of a discrete random variable becomes larger,
will the estimated power spectrum be smoother?
