Correlation

Step-by-step illustration of correlation between two signals with all intermediate plots

Correlation is a foundation over which the whole structure of digital communications is built. In fact, correlation is the heart of a digital communication system, not only for data detection but for parameter estimation of various kinds as well. Throughout, we will find recurring reminders of this fact.

As a start, consider from the article on Discrete Fourier Transform that each DFT output S[k] is just a sum of term-by-term products between an input signal and a cosine/sine wave, which is actually a computation of correlation. Later, we will learn that to detect the transmitted bits at the receiver, correlation is utilized to select the most likely candidate. Moreover, estimates of timing, frequency and phase of a received signal are extracted through judicious application of correlation for synchronization, as well as channel estimates for equalization. Don’t worry if you did not understand the last sentence, as we will have plenty of opportunity on this website to learn about these topics.

By definition, correlation is a measure of similarity between two signals. In our everyday life, we recognize something by running in our heads its correlation with what we know. Correlation plays such a vital and deep role in diverse areas of our life, be it science, sports, economics, business, marketing, criminology or psychology, that a complete book can be devoted to this topic.


"The world is full of obvious things which nobody by any chance ever observes."

Holmes to Watson – The Hound of the Baskervilles

For all of Sherlock Holmes’ inferences, his next step after observation was always correlation. For example, he accurately described Dr James Mortimer’s dog through correlating some observations with templates in his mind:


"\cdots and the marks of his teeth are very plainly visible. The dog’s jaw, as shown in the space between these marks, is too broad in my opinion for a terrier and not broad enough for a mastiff. It may have been — yes, by Jove, it is a curly-haired spaniel."

Sherlock Holmes – The Hound of the Baskervilles

As in the case of convolution, we start with real signals and the case of complex signals will be discussed later.

Correlation of Real Signals


The objective of correlation between two signals is to measure the degree to which those two signals are similar to each other. Mathematically, correlation between two signals s[n] and h[n] is defined as

(1)   \begin{equation*}         r_{sh}[n] = \sum \limits _{m = -\infty} ^{\infty} s[m] h[m-n]     \end{equation*}

where the above sum is computed for each n from -\infty to +\infty and the subscript s and h represent the order in which the signals appear on right hand side. We denote correlation operation by "\text{corr}" as

(2)   \begin{equation*}       r_{sh}[n] = s[n] ~\text{corr}~ h[n]     \end{equation*}

Note that unlike convolution,

(3)   \begin{equation*}         s[n] ~\text{corr}~ h[n] \neq  h[n] ~\text{corr}~ s[n]     \end{equation*}

This can be verified by plugging p = m-n in Eq (1) which yields m = n+p and hence

    \begin{align*}       \sum \limits _{p = -\infty} ^{\infty} s[n+p] h[p]  &= \sum \limits _{p = -\infty} ^{\infty}  h[p] s[p+n] \\                                                          & \neq \sum \limits _{m = -\infty} ^{\infty} h[m] s[m-n]     \end{align*}

Nevertheless, it can be deduced that Eq (1) is equivalent to

(4)   \begin{equation*}              r_{sh}[n] = \sum \limits _{m = -\infty} ^{\infty} s[m+n] h[m]     \end{equation*}

Now we can say that

    \begin{equation*} r_{sh}[n] = r_{hs}[-n] \end{equation*}

In terms of conveying information, there is not much difference and one is just a flipped version of the other.

Correlation computation


Comparing Eq (1) with convolution Eq, it is evident that

(5)   \begin{equation*}         r_{sh}[n] = s[n] * h[-n]     \end{equation*}

Therefore, from a viewpoint of conventional method, computing correlation between two signals is very similar to their convolution, except that there is no flipping of one signal. This is because h[-n] flips the signal once and convolution flips it again, hence bringing the original signal back.

From the viewpoint of intuitive method, it is clear that a negative sign with the NOW, -n, turns the future into past, and the past into future. Consequently, the last sample of the original signal arrives first, since it has become the farthest past.

Except this difference, correlation of real signals is very similar to convolution and the discussion on convolution accordingly applies here as well. For complex signals, there is another remarkable difference between the two, which we discuss next.

An example of correlation between the same two signals as in convolution example, s[n] = [2\hspace{1mm}-\hspace{-1mm}1\hspace{2mm} 1] and h[n] = [-1\hspace{2mm} 1\hspace{2mm} 2], is shown in Figure below, where the result r[n] is shown for each n.

Step-by-step illustration of correlation between two signals with all intermediate plots

Correlation of Complex Signals


Correlation between two complex signals s[n] and h[n] can be understood through writing Eq (1) in IQ form. However, another difference from convolution is that one signal is conjugated as

(6)   \begin{align*}         r_{sh}[n] &= s[n] ~\text{corr}~ h^*[n] \nonumber \\                   &= \sum \limits _{m = -\infty} ^{\infty} s[m] h^*[m-n]       \end{align*}

where conjugate of a signal was defined in the article on complex numbers. The above equation can be decomposed as in this Eq,

(7)   \begin{align*}          (r_{sh}[n])_I\: &= s_I[n] ~\text{corr}~ h_I[n] + s_Q[n] ~\text{corr}~ h_Q[n] \\         (r_{sh}[n])_Q &= s_Q[n] ~\text{corr}~ h_I[n] - s_I[n] ~\text{corr}~ h_Q[n]       \end{align*}

The actual computations can be written as

(8)   \begin{align*}          (r_{sh}[n])_I\: &= \sum \limits _{m = -\infty} ^{\infty} s_I[m] h_I[m-n] + \sum \limits _{m = -\infty} ^{\infty} s_Q[m] h_Q[m-n] \\         (r_{sh}[n])_Q &= \sum \limits _{m = -\infty} ^{\infty} s_Q[m] h_I[m-n] - \sum \limits _{m = -\infty} ^{\infty} s_I[m] h_Q[m-n]       \end{align*}

Due to the identity \cos A \cos B + \sin A \sin B = \cos (A-B), a positive sign in I term indicates that phases of the two aligned-axes terms are actually getting subtracted. Obviously, the identity applies in above equations only if magnitude can be extracted as common term, but the concept of phase-alignment still holds. Similarly, the identity \sin A \cos B - \cos A \sin B = \sin (A-B) implies that phases of the two cross-axes terms are also getting subtracted in Q expression. Hence, a complex correlation can be described as a process that

  • computes 4 real correlations: I ~\text{corr}~ I, Q ~\text{corr}~ Q, Q ~\text{corr}~ I and I ~\text{corr}~ Q
  • subtracts by phase anti-aligning the 2 aligned-axes correlations (I \diamond I + Q \diamond Q) to obtain the I component
  • subtracts the 2 cross-axes correlations (Q ~\text{corr}~ II ~\text{corr}~ Q) to obtain the Q component.

Now it can be inferred why a conjugate was required in the definition of complex correlation but not complex convolution. The purpose of correlation is to extract the degree of similarity between two signals, and whenever A is close to B,

    \begin{align*}       \cos(A-B) &\approx 1 \\       \sin(A-B) &\approx 0     \end{align*}

thus maximizing the correlation output.

Correlation and Frequency Domain


Just like convolution, there is an interesting interpretation of correlation in frequency domain. As before, DFT works with circular shifts only due to the way both time and frequency domain sequences are defined within a range 0 \le n, k \le N-1.

As always, we utilize the definition of DFT by applying Eq (1) to DFT definition. The derivation is similar to convolution.

Circular correlation between two signals in time domain is equivalent to multiplication of the first signal with conjugate of the second signal in frequency domain because

    \begin{equation*}         s[n] ~\boxed{\text{corr}}~ h[n] = s[n] \circledast h^*[-n]     \end{equation*}

which leads to

(9)   \begin{equation*}         s[n] ~\boxed{\text{corr}}~ h[n]  ~\xrightarrow{\text{\large{F}}} ~ S[k] \cdot H^*[k]      \end{equation*}

and the relation between h^*[-n] and H^*[k] can be established through the DFT definition.

Cross and Auto-Correlation


The correlation discussed above between two different signals is naturally called cross-correlation. When a signal is correlated with itself, it is called auto-correlation. It is defined by setting h[n] = s[n] in Eq (6) as

(10)   \begin{align*}         r_{ss}[n] &= s[n] ~\text{corr}~ s^*[n] \nonumber \\                   &= \sum \limits _{m = -\infty} ^{\infty} s[m] s^*[m-n]      \end{align*}

An interesting fact to note is that

(11)   \begin{align*}         r_{ss}[0] &= \sum \limits _{m = -\infty} ^{\infty} s[m] s^*[m] \nonumber \\                   &= \sum \limits _{m = -\infty} ^{\infty} |s[m]|^2 \nonumber \\                   &= E_s      \end{align*}

which is the energy of the signal s[n].

Remember that another signal can have a large amount of energy such that the result of its cross-correlation with s[n] can be greater than the auto-correlation of s[n], which intuitively should not happen. Normalized correlation \overline{r}_{sh}[n] is defined in Eq (12) in a way that the maximum value of 1 can only occur for correlation of a signal with itself.

(12)   \begin{align*}         \overline{r}_{sh}[n] &= \frac{\sum \limits _{m = -\infty} ^{\infty} s[m] h^*[m-n]}{\sqrt{\sum \limits _{m = -\infty} ^{\infty} |s[m]|^2} \cdot \sqrt{\sum \limits _{m = -\infty} ^{\infty} |h[m]|^2}} \nonumber \\         &= \frac{1}{\sqrt{E_s} \cdot \sqrt{E_h}} ~~ \sum \limits _{m = -\infty} ^{\infty} s[m] h^*[m-n]      \end{align*}

In this case, regardless of the energy in the other signal, its normalized cross-correlation with another signal cannot be greater than the normalized auto-correlation of a signal due to both energies appearing in the denominator.

Spectral Density


Taking the DFT of auto-correlation of a signal and utilizing Eq (9), we get

(13)   \begin{equation*}         r_{ss}[n] ~\xrightarrow{\text{\large{F}}} ~ S[k] \cdot S^*[k] = |S[k]|^2     \end{equation*}

The expression |S[k]|^2 is called Spectral Density, because from Parseval relation in the article on DFT Examples that relates the signal energy in time domain to that in frequency domain,

    \begin{equation*}         E_s = \sum _{n=0} ^{N-1} |s[n]|^2 = \frac{1}{N} \sum _{n=0} ^{N-1} |S[k]|^2     \end{equation*}

Thus, energy of a signal can be obtained by summing the energy |S[k]|^2 in each frequency bin (up to a normalizing constant 1/N). Accordingly, |S[k]|^2 can be termed as energy per spectral bin, or spectral density.

From the above discussion, there are two ways to find the spectral density of a signal:

  1. Take the magnitude squared of the DFT of a signal.
  2. Take the DFT of the signal auto-correlation.

If you found this article useful, you might want to subscribe to my email list below to receive new articles.


Linear Systems

A linear system with scaled input and output

A linear system implies that if two inputs are scaled and summed together to form a new input, the new output of the system is also a scaled sum of their individual outputs.

[Scaling] For scaling to hold, if

    \begin{align*}               \textmd{input} \quad s_1[n] ~& \xrightarrow{\hspace*{2cm}}~ r_1[n] \quad \textmd{output}             \end{align*}

then

    \begin{align*}                 \textmd{input}\quad \alpha s_1[n] ~& \xrightarrow{\hspace*{2cm}}~ \alpha r_1[n] \quad \textmd{output}             \end{align*}

where \alpha is any scalar.

[Addition] When two such inputs are added together, the output should be the sum of their individual outputs as

    \begin{align*}               \textmd{input}\quad s_1[n] + s_2[n] ~& \xrightarrow{\hspace*{2cm}}~ r_1[n] + r_2[n] \quad \textmd{output}             \end{align*}

A linear system combines the above two properties as

(1)   \begin{align*}       \textmd{input}\quad \alpha_1 s_1[n] + \alpha_2 s_2[n] ~& \xrightarrow{\hspace*{2cm}}~ \alpha_1 r_1[n] + \alpha_2 r_2[n] \quad \textmd{output}     \end{align*}

Below, we discuss examples of a linear and a non-linear system.

Example


Consider a system

    \begin{equation*}         r[n] = \frac{3}{7} \cdot s[n]     \end{equation*}

The output of this system, as a response to an input s_1[n] = \sin(2\pi 0.1n), is

    \begin{equation*}       r_1[n] = \frac{3}{7}\sin(2\pi 0.1n)     \end{equation*}

Similarly, response to a different signal s_2[n] = \sin(2\pi 0.3n) is

    \begin{equation*}       r_2[n] = \frac{3}{7}\sin(2\pi 0.3n)     \end{equation*}

When this system is given the input s_1[n] + s_2[n], the output is

    \begin{align*}         r[n] &= \frac{3}{7}\left\{\sin(2\pi 0.1n) + \sin(2\pi 0.3n)\right\} \\              &= \frac{3}{7}\sin(2\pi 0.1n) + \frac{3}{7}\sin(2\pi 0.3n) \\              &= r_1[n] + r_2[n]     \end{align*}

Hence, it is a linear system.

On the other hand, when the same input s_1[n] + s_2[n] is given to another system

    \begin{equation*}         r[n] = s^2[n]     \end{equation*}

and using the identity \sin(A)\sin(B) = \frac{1}{2}\{\cos(A-B) - \cos(A+B)\}, the output is

    \begin{align*}         r[n] &= \left\{\sin(2\pi 0.1n) + \sin(2\pi 0.3n) \right\}^2 \\              &= \sin^2(2\pi 0.1n) + \sin^2(2\pi 0.3n) + 2\sin(2\pi 0.1 n) \sin(2\pi 0.3 n) \\              &= \sin^2(2\pi 0.1n) + \sin^2(2\pi 0.3n) + \cos(2\pi 0.2n) - \cos(2\pi 0.4n) \\              & \neq r_1[n] + r_2[n] \\              &= \sin^2(2\pi 0.1n) + \sin^2(2\pi 0.3n)     \end{align*}

Clearly, it is a non-linear system.

From above example, it is also clear that input sinusoids do not interact with each other in linear systems, and hence output frequencies were the same as the input frequencies. In a non-linear system, however, input sinusoids interacted with each other to produce frequencies that were not present in either of the input signal, as shown in Figure below for r[n] = s^2[n]. Note from |S[k]| that s(t) actually consists of only two frequencies at bins 1 and 2, but the output [r[n] of the system is composed of other frequencies as well shown in |R[k]|.

A non-linear system shown with output containing extra frequency components as compared to input

The Discrete Fourier Transform, DFT, is a linear operation as it is evident from DFT definition that any scaling and addition of two or more input signals will result in a DFT output that is a scaled and summed version of their individual DFT outputs.

The Concept of Phase

Changing phase with frequency indices

We have explained what Discrete Fourier Transform (DFT) is. Also, we have covered that the concept of frequency is related to rotational speed of a complex sinusoid. Subsequently, a frequency axis was defined first in continuous domain and then in discrete domain. The tools to view the spectrum of a signal in frequency domain are IQ and magnitude-phase plots. We defined the magnitude and phase of a complex number before. Similar definitions hold for complex signals as

(1)   \begin{equation*}     |S[k]| = \sqrt{S_I[k]^2 + S_Q[k]^2}      \end{equation*}

and the phase \measuredangle S[k] is defined as four-quadrant inverse tangent as in this Eq with S_I[k] and S_Q[k] in place of V_I and V_Q, respectively.

Magnitude-phase plots are usually used more than IQ plots because a magnitude plot shows the strength of a complex sinusoid at each frequency in the whole spectrum. In focusing on the magnitude plot, sometimes it is easy to miss a great deal of information provided by the phase plot. Although phase is relatively unimportant for some particular areas such as most audio applications due to relative insensitivity of human ear to phase, it plays a significant role in wireless communications as we will see in later chapters.

As an example of what happens when phase information is neglected, this is how my daughter writes some English letters:

  • L \rightarrow \Gamma
  • V \rightarrow \Lambda
  • Z \rightarrow a reversed \Gamma

All her symbols above are almost correct. Nevertheless, their phase is distorted leading to incorrect results. In addition, further confusion can develop if the relationship of phase with time domain is not clearly understood.

Phase in frequency domain has a special relationship with initial sample of a signal in time domain. Intuitively, a waveform is about how a signal changes in time IQ-plane while its first samples on I and Q axes indicate the starting point of that waveform. On the other hand, frequency is about rotational speed of complex sinusoids constructing that signal (the higher the frequency, the farther it is on frequency axis but not rotating) and phase indicates their orientation on frequency IQ-plane. Naturally, this orientation of each such complex sinusoid depends on where its initial sample is.

Example


As an example, consider this Figure. In time domain, when the starting sample of a cosine is changed by 1/4 of a period, it becomes a sine wave. Correspondingly in frequency domain, the cosine changes its phase by 1/4 of 360^{\circ} =90^{\circ} (to be exact, -90^ {\circ} and +90^ {\circ} on positive and negative frequency axis, respectively). Mathematically, remember that \cos(A-90^{\circ}) = \sin(A).

Consider the inphase component of a complex sinusoid shifted by n_0 samples:

    \begin{align*}         \cos 2\pi \frac{k}{N} \left(n-n_0 \right) &= \cos \left( 2\pi \frac{k}{N} n - 2\pi \frac{k}{N} n_0 \right)  \\                                                                      &= \cos \left( 2\pi \frac{k}{N} n + \Delta \theta (k) \right)     \end{align*}

Since n_0 is a constant above, \Delta \theta (k) = - 2 \pi (k/N)n_0 can be seen as the phase shift incurred by a delay of n_0 samples. This result is known as the shifting property of the DFT which holds true for circular shifts in time. It states that a (circular) time shift of an input signal s[n] results in a corresponding phase shift at each frequency of its DFT S[k].

Effect of time shift on DFT – Magnitude and Phase In the light of discussion above, the DFT of s[(n-n_0) \:\textmd{mod}\: N] has its magnitude unchanged. However its phase is rotated by \Delta \theta (k) = -2\pi (k/N) n_0. We denote the rotated DFT by \widetilde{S}[k].

(2)   \begin{equation*}           \begin{aligned}             |\widetilde{S}[k]| &= |S[k]| \\             \measuredangle \widetilde{S}[k] &= \measuredangle S[k] - 2\pi \frac{k}{N} n_0           \end{aligned}         \end{equation*}

Effect of time shift on DFT – I and Q It is straightforward to prove through DFT definition that the DFT of s[(n-n_0) \:\textmd{mod}\: N] is given by

(3)   \begin{align*}         \widetilde{S}_I[k]\: &= S_I[k] \cos 2\pi \frac{k}{N} n_0 + S_Q[k] \sin 2\pi \frac{k}{N} n_0 \\         \widetilde{S}_Q[k] &= S_Q[k] \cos 2\pi \frac{k}{N} n_0 - S_I[k] \sin 2\pi \frac{k}{N} n_0       \end{align*}

As a verification, comparing with this Eq, \widetilde{S}[k] is nothing but rotations of complex numbers S[k] by angles \Delta \theta (k) = -2\pi (k/N) n_0 for each k = -N/2, \cdots,-1,0,1, \cdots,N/2-1.

The converse of the above argument is also true. A phase shift at a discrete frequency bin of the DFT informs us about (circular) time shift of that sinusoid. The conclusions from above are summarized in the note below.

Time shift ~\xrightarrow{\text{{F}}}~ Phase shift


If a signal s[n] is circularly right shifted in time by n_0 samples (i.e., samples are moved n_0 places to the right, with elements that fall off at one end of the sequence appearing at the other end), then the magnitude of its DFT |S[k]| remains unchanged. However, the phase of its DFT \measuredangle S[k] gets rotated by \Delta \theta (k) = -2\pi (k/N) n_0 for each k = -N/2, \cdots,-1,0,1, \cdots,N/2-1. Similarly, for a circular left shift in time by n_0 samples incurs a phase rotation of \Delta \theta (k) = +2\pi (k/N) n_0 for all k.

(4)   \begin{align*}               \textmd{Time shift} \quad s[(n-n_0)  \:\textmd{mod}\: N] ~& \xrightarrow{\hspace*{1cm}}~ -2\pi \frac{k}{N} n_0 \quad \textmd{Phase shift} \\               \textmd{Time shift} \quad s[(n+n_0)  \:\textmd{mod}\: N] ~& \xrightarrow{\hspace*{1cm}}~ +2\pi \frac{k}{N} n_0 \quad \textmd{Phase shift} \\         \end{align*}

Therefore, a time delay (going back in time from NOW, or waveform shift to the right) rotates the original spectral phase in negative direction (clockwise). On the other hand, a time advance (future travel from NOW, or waveform shift to the left) rotates the original spectral phase in positive direction (anticlockwise). In terms of intuitive method of time shifting, it makes perfect sense: Traveling in the past should decrease the DFT phase, and vice versa.

In examples discussed above, the amount of phase shift has a linear relation with the frequency index k as evident from the term \Delta \theta (k) = -2\pi (k/N) n_0. This is the concept of linear phase where the phase of each complex sinusoid is directly proportional to the frequency of that sinusoid. Intuitively, if sinusoids with different frequencies get delayed by the same number of samples, then they naturally end up with different phases at the end of that common sample duration.

The magnitude and direction of rotation for these frequencies is symbolically shown for s[n+n_0] in Figure below. Remember that this is a frequency IQ-plane, unlike time IQ-plane.

Changing phase with frequency indices

Example


The Figure below shows a unit impulse signal s[n] and its DFT S[k] along with its circularly time shifted version s[(n-1) \:\textmd{mod}\: 5] and its DFT \widetilde{S}[k]. This DFT is computed shortly in later article.

Note that phase shift for each frequency bin k is different for each k according to Eq (4). To be exact, \Delta \theta (k) = 2\pi (k/N) \cdot (-1) for k = -2 to k = 2 and N = 5 turns out to be

    \begin{align*}     k &= ~~~0 \quad \rightarrow \quad \Delta \theta (0) ~~\:= 2\pi \frac{0}{5} \cdot (-1) \times \frac{180^\circ}{\pi} = 0^\circ \\     k &= \pm 1 \quad \rightarrow \quad \Delta \theta (\pm 1) = 2\pi \frac{\pm 1}{5} \cdot (-1) \times \frac{180^\circ}{\pi} = \mp 72^\circ \\     k &= \pm 2 \quad \rightarrow \quad \Delta \theta (\pm 2) = 2\pi \frac{\pm 2}{5} \cdot (-1) \times \frac{180^\circ}{\pi} = \mp 144^\circ     \end{align*}

The phase rotations of 72^\circ and 144^\circ are illustrated in the figure. Also observe that for a right shift, the phase rotations are clockwise for positive k and anticlockwise for negative k.

A time shift in a sequence generates a phase shift in frequency domain to all frequency indices

The understanding of phase above cannot be overemphasized. It is relatively straightforward to see magnitude plots and diagnose the behavior of signals and systems. However, true insights can only be developed through grasping the implications of phase rotations.

If you found this article useful, you might want to subscribe to my email list below to receive new articles.


A Digital Signal

A digital signal and its underlying continuous waveform

We have talked about obtaining a discrete-time signal through sampling the time-axis and obtaining a discrete frequency set through sampling the frequency axis. The same concept can be applied to the amplitude-axis, where the signal amplitude can be sampled to take only a finite set of discrete values. This discrete-time discrete-valued signal is called a digital signal, as opposed to an analog signal that is continuous in time and continuous in amplitude.

The above Figure shows how a digital signal having amplitudes over a fixed set of values can be obtained through slicing the underlying continuous amplitudes. For example, an amplitude of 2.2 can be rounded to 2, 1.4 to 1 and so on depending on the desired resolution.

Computers can only work with digital signals because discrete-time signals — though defined only for finite values on time-axis — can have infinite values on the amplitudes-axis. Just as computer memory is finite and can store only a known number of time values, its width is also finite (e.g., 8 bits) and can store only a fixed number of amplitudes (e.g., for an 8-bit wide memory, we can have 256 values for amplitudes).

Pulse Amplitude Modulation (PAM)

Block diagram of a pulse amplitude modulator and demodulator

In the article on modulation – from numbers to signals, we said that the Pulse Amplitude Modulation (PAM) is an amplitude scaling of the pulse p(nT_S) according to the symbol value. What happens when this process of scaling the pulse amplitude by symbols is repeated for every symbol during each interval T_M? Clearly, a series of bits \underline{b} (1010 in our initial example) can be transmitted by choosing a rectangular pulse and scaling it with appropriate symbols.

    \begin{align*}         m = 0 &\rightarrow b = 1 &\rightarrow a[0] = +A \\         m = 1 &\rightarrow b = 0 &\rightarrow a[1] = -A \\         m = 2 &\rightarrow b = 1 &\rightarrow a[2] = +A \\         m = 3 &\rightarrow b = 0 &\rightarrow a[3] = -A     \end{align*}

Our next step is forming a cumulative waveform from these individual symbol-scaled pulses. Remember from the article on transforming a signal that mathematical expression for a signal p(nT_S) delayed by an amount T_M or L samples is given as p(nT_S-LT_S) = p(nT_S-T_M), where T_M is the symbol duration and L is samples/symbol defined as T_M/T_S. Since the same pulse is scaled by the symbol value during each T_M,

  • At time instant 0T_M, the output is a[0] p(nT_S-0T_M).
  • At time instant 1T_M, the output is a[1] p(nT_S-1T_M).
  • At time instant 2T_M, the output is a[2] p(nT_S-2T_M).
  • At time instant 3T_M, the output is a[3] p(nT_S-3T_M).

And so on. Finally, their addition gives the expression for a general PAM waveform:

(1)   \begin{align*}         s(nT_S) &= a[0]p(nT_S-0T_M) + a[1]p(nT_S-1T_M) + \nonumber \\              & \hspace{1in} a[2] p(nT_S-2T_M) + a[3]p(nT_S-3T_M) + \cdots \nonumber \\              &= \sum _{m} a[m] p(nT_S-mT_M)      \end{align*}

After digital to analog conversion (DAC), the continuous-time signal s(t) can be expressed as

(2)   \begin{align*}         s(t) &= \sum _{m} a[m] p(t-mT_M)      \end{align*}

As an example, a 2-PAM waveform is illustrated in Figure below with red dashed curve being the underlying continuous-time signal.

A Pulse Amplitude Modulation waveform is constructed by summation of individually modulated bits

In a similar manner to 2-PAM, a 4-PAM waveform based on 4 symbols can be constructed by scaling the pulse amplitude by different symbol values during each T_M, as illustrated in Figure below.

A 4-PAM digital signal with underlying continuous waveform

Constellation Diagram


Just like a constellation of stars, a constellation diagram shows the actual symbol values representing a set of \log_2 M bits. We have already encountered constellation diagrams before (e.g., in the article on a simple communication system). A general constellation diagram for M-PAM is shown in Figure below.

Constellation diagram for general Pulse Amplitude Modulation

Average Symbol Energy


The average symbol energy in a constellation is given by the average of all individual symbol energies. For M=2,

    \begin{align*}         E_{M-\text{PAM}} &= \frac{A^2}{2}\left\{(+1)^2 + (-1)^2 \right\} = A^2     \end{align*}

And for M=4,

    \begin{align*}         E_{M-\text{PAM}} &= \frac{A^2}{4}\left\{(-3)^2 + (-1)^2 + (+1)^2 + (+3)^2\right\} = 5A^2     \end{align*}

For a general M,

    \begin{align*}         E_{M-\text{PAM}} &= \frac{A^2}{M}\left\{(-M+1)^2 + \cdots (-3)^2 + (-1)^2 + \right.\\                          &~ \left. \hspace{2in} (+1)^2 + (+3)^2 + \cdots + (M-1)^2\right\} \\                          &= 2\frac{A^2}{M}\left\{1^2 + 3^2 + \cdots + (M-1)^2\right\}     \end{align*}

The term in the brackets is the sum of squares of first (M-1+1)/2 = M/2 odd integers. Using the formula k(2k-1)(2k+1)/3 for the first k odd integers squared, we get

(3)   \begin{align*}         E_{M-\text{PAM}} &= 2\frac{A^2}{M} \cdot \frac{M/2 \cdot (M-1) \cdot (M+1)}{3} \nonumber \\                          &= \frac{M^2-1}{3} A^2      \end{align*}

The main purpose of this text is to find the answer to the following question: which operations are necessary to perform at the Tx and Rx sides such that we can detect that symbol with the highest probability? To explore the answer, we start with a simple PAM modulator and detector.

PAM Modulator


At this stage, we are ready to build a conceptual PAM modulator. The block diagram is drawn in Figure below in which Tx signal is generated in the following way.

  • Every T_b seconds, a new bit arrives at the input forming a serial bit stream.
  • A serial-to-parallel (S/P) converter collects \log_2 M such bits every T_M = \log_2 M \times T_b seconds that are used as an address to access a Look-Up Table (LUT) that stores M symbol values specified by the constellation.
  • To produce a PAM waveform, the symbol sequence a[m] is converted to a discrete-time impulse train through upsampling by L, where L is samples/symbol defined as ratio of symbol time to sample time T_M/T_S, or equivalently sample rate to symbol rate F_S/R_M.
  • As explained in the post on sample rate conversion, upsampling inserts L-1 zeros between each symbol after which the interpolated intermediate samples can be raised from dead with the help of a lowpass filter that suppresses all the spectral replicas except the primary one. We will see in the next section that a proper pulse shaping filter is also a lowpass filter and hence another extra lowpass filter is not actually required.
  • The generated discrete-time signal s(nT_S) is converted to a continuous-time signal s(t) by a DAC.

The mathematical derivation for the PAM modulator was shown in Eq (1) and Eq (2).

Block diagram of a pulse amplitude modulator and demodulator

PAM Detector


The received signal r(t) is the same as the transmitted signal s(t) but with the addition of additive white Gaussian noise (AWGN) w(t). The symbols are detected through the following steps illustrated in Figure above.

  • Through an analog to digital converter (ADC), r(t) is sampled at a rate of F_S samples/s to produce a sequence of T_S-spaced samples r(nT_S).
  • Next, r(nT_S) is processed through a matched filter at the Rx side to generate z(nT_S). As discussed earlier, the output of the matched filter z(nT_S) is a continuous correlation of the symbol-scaled pulse shape with an unscaled and time-reversed pulse shape.
  • This output is downsampled by L at optimal sampling instants

    (4)   \begin{equation*}                     \tcbhighmath{n = mL = m \frac{T_M}{T_S}}                 \end{equation*}

    to produce T_M-spaced numbers z(mT_M) back from the signal.

  • The minimum distance decision rule is employed to find the symbol estimates \hat a[m].

Take a special note of Eq (4). It will be employed over and over again.

Notice that a symbol is the basic building block of a digital communication system. Consequently, symbol time T_M is the basic unit of measurement along the time axis of such a system. When Figure on correlator outputs depicts sampling the output just once at optimum instant 0 out of -(L-1) \le n \le 0, it is true for all integer multiples of symbol time as well, i.e., the output is sampled just once for every symbol duration at optimum instants T_M, 2T_M, \cdots.

Key samples


Carefully examine the key samples at T_M, 2T_M, \cdots. These are the samples we are looking for in the waveform for the detection purpose. Even when the waveform has suffered from all the distortions the real world has to offer, locating these samples and mapping them back to the constellation is a beautiful process, actual details of which we will encounter throughout this text.

Let us discuss the mathematical details of this process. The Tx signal in Eq (1) is expressed as (the reason for using a different variable i instead of m will shortly become clear)

    \begin{equation*}         s(nT_S) = \sum _{i} a[i] p(nT_S-iT_M)     \end{equation*}

In a noiseless case, this signal is input to matched filter h(nT_S) = p(-nT_S) and the output is written as

(5)   \begin{align*}         z(nT_S) &= s(nT_S) * p(-nT_S) \nonumber\\                 &= \sum \limits _k \Big(\sum \limits _i a[i] p(kT_S - iT_M)\Big) p\Big(-(nT_S - kT_S)\Big)  \nonumber\\                 &= \sum \limits _i a[i] \cdot \sum \limits _k  p(kT_S - iT_M) p\Big(kT_S -iT_M - (nT_S - iT_M)\Big)  \nonumber\\                 &= \sum \limits _i a[i] r_p(nT_S - iT_M)      \end{align*}

where r_p(nT_S) comes into play from the definition of auto-correlation function. To generate symbol decisions, T_M-spaced samples of the matched filter output are required at n = mL = mT_M/T_S. Downsampling the matched filter output generates

(6)   \begin{align*}         z(mT_M) &= z(nT_S) \bigg| _{n = mL = mT_M/T_S} \nonumber \\                 &= \sum \limits _i a[i] r_p(mT_M - iT_M) = \sum \limits _i a[i] r_p\{(m-i)T_M\} \nonumber \\                 &= a[m] \nonumber     \end{align*}

This is because for a rectangular pulse shape, the matched filter output is triangular with maximum occurring at i=m and zero at the next symbol location (since we wanted to denote our current symbol with m, we opted for variable i at the start of this derivation), see Figure on correlation of a rectangular pulse.

Observe that the system shown in PAM system block diagram is a multirate system. In the PAM detector, for example, the ADC and the matched filter operate at the sample rate F_S. After the output of the matched filter is downsampled by L, the symbol decisions are made at the symbol rate R_M. Furthermore, there are some hidden assumptions in the PAM detector:

[Resampling] The ADC in general does not produce an integer number of samples per symbol, i.e., T_M/T_S is not an integer. As we will see later, a resampling system is required in the Rx chain that changes the sample rate from the ADC rate to a rate that is an integer multiple of the symbol rate.

[Symbol Timing Synchronization] The peak sample at the end of each symbol duration is not known in advance at the Rx and in fact does not necessarily coincide with a generated sample as well. This is because ADC just samples the incoming continuous waveform without any information about the symbol boundaries. This is a symbol timing synchronization problem which we will learn later.

If you found this article useful, you might want to subscribe to my email list below to receive new articles.