Description

Instructions
Write down the answers and gures in a PDF document and submit it in Brightspace. When you’re asked to give a sketch, you can make a drawing by hand, and scan it (or make a photo of it).
The report will be checked and the result of the assignment will be either PASS or FAIL. The grade for the written exam will only count in case you pass this assignment.

Signal Processing for RadioPulsar Navigation
2.1 Introduction
A pulsar (short name for pulsating star), is a rotating neutron start that emits electromagnetic radiation. While rotating, the star emits a beam of radiation. This radiation can be observed on earth as a pulse when the beam is pointed towards earth. You can compare this with a lighthouse, where the light can only be seen when the lighthouse is pointing towards the observer.
Neutron stars have short and very regular rotational periods. As a result, the received pulses can be observed in intervals that are extremely regular. Depending on the star, the intervals (periods) between the pulses range from milliseconds to seconds. For some pulsars, the periodicity of the received pulses is more precise than an atomic clock. That makes pulsars extremely useful for many applications.
One of these (future) applications could be navigation based on pulsars, as the preciseness of the periodic pulses can be used for timeofarrival based navigation. Among the advantages of pulsar based navigation are the fact that it is independent of satellites (like GPS) and the fact that the pulsar signals have an extremely wide frequency range, which makes it di cult to distort by potential enemies.
Besides the advantages of pulsars they carry one disadvantage, and that is the fact that the signaltonoise ratio (SNR) of the received pulses is extremely low (around 90 dB). That makes the pulsar signals very di cult to detect.
This assignment is about the detection of pulsar signals and how this can be done using concepts of random signal processing.
2
2.2 Modeling Pulsar Signals
Usually, it is assumed that the received pulsar signals are degraded by additive noise that is uncorrelated with the pulsar signal. That is,

y(t) = x(t) + n(t);
(1)
where t indicates the (continuous) time index, x(t) the pulsar signal, n(t) the noise and y(t) the received noisy pulsar signal. The pulsar signal is assumed to be deterministic, while the noise n(t) is assumed to be a stochastic stationary and ergodic uncorrelated signal with zero mean, i.e., E[N(t)] = 0, and variance _{N}^{2} . For many pulsars, the period of the pulsar signal x(t), say T_{0}, is known by approximation. This information is obtained by analyzing extremely long recordings of received noisy pulsar signals. By performing speci c signal processing operations on these recordings it is also possible to estimate one period of a pulsar signal x(t). A method that is traditionally often used is called epoch folding. With epoch folding one divides the noisy signal in frames of length T_{0} and averages the data over these periods. As the noise is assumed to be uncorrelated, it will cancel out when enough averaging operations are performed.
Exercise 1
Sketch the autocorrelation function of the noise process N(t)
2.3 Epoch folding
The Epoch folded signal is given by

1
K 1
X_{k}
^{z}e^{(t) =} _{K}
(2)
y(t + kT_{0});
=0
where K denotes the number of frames over which the signal is averaged. Substituting Eq. (1) into (2) we obtain

1
K 1
1
K 1
X
X_{k}
^{z}e^{(t) =} _{K}
K
(3)
x(t + kT_{0}) +
n(t + kT_{0}) ;
k=0
=0

{z
}

{z
}
x_{0}(t)
n(t) with n(t) ! 0
as K!1
3
Figure 1: Example of epoch folding.
where x_{0}(t) indicates one complete period of the pulsar signal and n(t) the averaged noise. As the noise is an uncorrelated zeromean process, n(t) will average out if K is su ciently large. An example of this procedure is given in Fig. 1. To get an impression how well the pulse can be detected in the noise, one can compute the SNR. Let the SNR of the measured noisy signal, computed over one period, be given by

0
R
^{T}_{=0}^{0} x_{0}(t)^{2}dt
1
SN R = 10 log_{10}
@
E
t
T_{0}
n(t)^{2}dt
i^{A}
:
(4)
h_{R}
t=0
The SNR of the epoch folded signal is given by

0
R
^{T}^{0} x_{0}(t)^{2}dt
1
SN R_{e} = 10 log_{10}
@
E
A
:
(5)
t
T_{0}
n(t)^{2}dt
=0
^{hR}t=0
i
Exercise 2
Compute how much the SNR of the epoch folded signal has increased compared to the original SNR if K = 10 averaging operations are performed (that is, compute SN R_{e} SN R).
Exercise 3
How large should K be, in order to obtain an SNR increase of 90 dB using epoch folding? Assume that one period of the pulsar signal takes 0.5 seconds. How many hours of data do we then need to increase the SNR with 90 dB?
4
From the blackboard website you can download a Matlab data le named pulsar data.mat. Download this le and load it using the command
load pulsar data,
This data le contains a synthetically generated example pulsar signal in the variable x, the noise in variable n, the mixed signal in variable y and a variable called template that we will use later on. These signals are sampled at 10 kHz.
Exercise 4
Determine the periodicity T_{0} in samples by observing signal x. Also compute the SNR of the original signal using signals x and n.
Exercise 5
Write an m le [z_{e}]=epoch folding(signal, T_{0}, K);
This m le takes a signal as input, together with the period T_{0} and the number of frames K over which we average. The output is an epoch folded signal, i.e., one period of the signal.
Exercise 6
Use the m le [z_{e}]=epoch folding(signal, T_{0}, K); to perform epoch folding in order to detect the pulse in the pulsar signal. Use for input variable signal the noisy signal y and use for input variable T_{0} the period that you determined in exercise 4. For K, use the values that should give you an SNR increase of 10, 20, and 30 dB. Give the corresponding values for K and plot the three resulting epoch folded signals for the three di erent SNRs. Can you detect the position of the pulse in the signal?
2.4 Matched ltering
The epoch folding approach estimates the signal of interest. However, in order to have a su ciently high output SNR such that detection is easy, K should be very large for these extremely low SNRs. Moreover, we are mainly interested in the position of the pulse in the signal, and not so much in the signal itself. Furthermore, for many pulsars, the signal x is known as researchers performed epoch folding techniques on extremely long data
5
records. This means that templates of one period of the signal x(t), i.e., x_{0}(t), are available. The availability of such a template simpli es detection of the pulse a lot as additional information on the signal is used. The data le that you loaded also contains one period of a template in the variable template.
A wellknown method that makes use of such a template is the matchedltering approach. The matched lter has a very speci c impulse response h(t). Under the signal model y = x + n we can write the output of the matched lter by the following convolution

z_{m}(t) =
Z +1
y(s)h(t
s)ds
(6)
Z +1
=
x(s)h(t
s)ds + ^{Z} ^{+1}
n(s)h(t s)ds
(7)
=
z_{x}(t) + z_{n}(t);
(8)
where z_{x}(t) and z_{n}(t) are the signals that result when the impulse response of the matched lter is applied to signal x(t) and n(t), respectively. To derive the impulse response of the matched lter, the per sample ratio between signal power and noise variance is maximized, that is,

( ) =
jz_{x}( )j^{2}
:
(9)
E [jz_{n}( )j^{2}]
Using the CauchySchwartz inequality, it can be shown that ( ) is maximized by choosing h(t) as h(t) = cx( t + ), that is, a mirrored and shifted version of x(t) (due to t and in the argument of x( t + ), respectively) and where c is an arbitrary constant. This impulse response can also be seen as a template of the signal of interest. By convolving the noisy signal y with this impulse response, a correlation between the noisy signal y and the known template of x is essentially computed, that is,

+1
z_{m}(t) = c y(s)x( t + + s)ds: (10)
For now we consider the case for which = 0, that is,

+1
z_{m}(t) = c y(s)x( t + s)ds: (11)
6
2.5 Relation Between Matched Filtering and Epoch Folding
By limiting the integration boundary of Eq. (11) to an integer multiple of the pulse period T_{0}, say KT_{0}, it is possible to show that the matched ltering can be implemented as the correlation of one pulse period of the template and the estimate z_{e}(t) obtained by epoch folding, that is,

^{KT0}

z_{m}(t) =
y(s)x(
t + s)ds
(12)
=
0
^{K} k=0
y(s + kT_{0})^{!}
x( t + s)ds
(13)
Z_{0}
T_{0}
1
K 1
X

^{T0}

=
z_{e}(s)x( t + s)ds;
(14)
0
where constant c is set to c = _{K}^{1} From this we can conclude that the matched lter can indeed be implemented as the correlation between one pulse period of the template, i.e., x( t + s), and the estimate z_{e}(t) that is obtained by epoch folding.
By substituting Eq. (3) into Eq. (14) this can also be written in terms of x_{0}(t) and n(t),
^{Z}T_{0} ^{Z}T_{0}
z_{m}(t) = x_{0}(s)x( t + s)ds + n(s)x( t + s)ds: (15)
0 0
Exercise 7
Assume that x_{0}(t) consists of a perfect pulse as shown in Fig. 2. Make a sketch of ^{R}_{0}^{T}^{0} x_{0}(s)x( t + s)ds. Explain the shape of the sketch.
The matched lter output z_{m}(t) given by Eq. (15) still contains some noise that is expressed by the second term. The expected value of the statistical process that underlies z_{m}(t), i.e., E [Z_{m}(t)], can be shown to be

E [Z_{m}(t)] = T_{0}R_{X}_{0}_{X} ( t):
(16)
Exercise 8
Derive the expression in Eq. (16) and explain which statistical assumptions you make in order to derive this.
7
Figure 2: Example of one period of x(t).
In Matab we work with data that is sampled at a certain sampling frequency f_{s}. To implement Eq. (15) in Matlab, the integrals in Eq. (15) need to be written as summations and instead of continuous time we use sampled time. This leads to

1
N 1
X
z_{m}(l) =
z_{e}(m)x( l + m)
(17)
f_{s}
m=0

1

=
T_{0}
X
z_{e}(m)x(
l + m)
(18)
N
m=0
=
T_{0}
N 1
x_{0}(m)x(
l + m) +
T_{0}
N 1
n(m)x( l + m); (19)
X
X
N
m=0
N
m=0
where the relation ^{1} = ^{T}0 is used.
f_{s} N
Exercise 9
For large N, i.e., N ! 1, z_{m}(l) converges to

z_{m}(l)
^{T}0^{R}X_{0}X^{(}
l) + T_{0}R_{NX} ( l)
(20)
=
^{T}0^{R}X_{0}X^{(}
l):
(21)
Which property related to estimation of autocorrelation functions should the underlying processes have such that the approximation from Eq. (18) to (20) is valid.
Eq. (18) resembles a correlation estimator between the epoch folded signal z_{e}(l) and the template x( l). Implementing this estimator can be done
8
by writing Eq. (18) in terms of a matrixvector multiplication z_{m} = xz_{e}, i.e.,

2
z_{m}(1)
3
T_{0}
2
x( 1)
x(0)
: : :
x(N
2)^{32}
z_{e}(1)
3
z_{m}(0)
x(0)
x(1)
: : :
x(N
1)
z_{e}(0)
6
_{.}.^{.}
7
=
6
_{.}.^{.}
^{.} . _{.}
_{.}.^{.}
76
_{.}.^{.}
7
N
N +2)
^{6}z_{m}(N
1)^{7}
6
x( N + 1) x(
: : :
x(0)
^{76}z_{e}(N
1)^{7}
6
7
6
76
7
4
5
4
54
5
As x(m) is periodic with N samples, we can replace x(m) by x(m mod N). Leading to z_{m} = x_{mod N} z_{e}, that is,

2
z_{m}(1)
3
^{T}0
2
x(N
1) x(0)
: : :
x(N
2)
32
z_{e}(1)
3
6
z_{m}(0)
7
6
x(0)
x(1)
: : :
x(N
1)
76
z_{e}(0)
^{7}(22):
_{.}.^{.}
=
_{.}.^{.}
^{.} . _{.}
_{.}.^{.}
_{.}.^{.}
N
x(2)
^{6}z_{m}(N
1)^{7}
6
x(1)
: : :
x(0)
76
z_{e}(N
1)^{7}
6
7
6
76
7
4
5
4
54
5
To construct the matrix x_{mod N} you can make use of the Matlab command toeplitz.
Exercise 10
Write an m le [z_{m}]=matched filter[z_{e},x], which performs the above matrixvector multiplication in order to compute the output of the matched lter. Make use of the Matlab command toeplitz to implement above matrix structure. To test the implementation use z_{e} = [x_{0}(0); :::; x_{0}(N 1)]^{T} and compute z_{m} = xz_{e}. Make a plot of the result. Is this in line with your sketch from Exercise 7?
Exercise 11
We will now use this implementation of the matched lter to determine the pulse location of the pulsar signal. First generate epoch folded signals using the function epoch folding on the noisy signal y for the Kvalues that you found in Exercise 6. Use each of these epoch folded signals as input z_{e} in the function matched filter. For each Kvalue, make a plot of the epochfolded signal and the matched lter output. Which value of K do you need such that you can detect the pulse using epoch folding? Which value of K do you need to detect the pulse if you combine epochfolding and the matched lter? Give the relative position of the detected pulse in terms of seconds compared to the template.
9