Continuous time signals are basically analog signals, which occur all around us in nature. Even though most of the signal processing is done in digital domain these days, we still need to convert the output waveform back to continuous time analog signals for the human beings to perceive them.
In this post, I am going to describe how we can use MATLAB to approximate a continuous time periodic signal using exponential Fourier series. Although, bear in mind that, for a Fourier series to exists, the signals must satisfy certain criteria, which might be discussed in upcoming posts.
The general form of an exponential Fourier series is,
First of all, lets determine the expression for computing the co-efficients of the FS. We know that the period of the square wave is T = 2 sec. So,
% Approximating a periodic square wave using exponential fourier series clear all; clc; clf; N = 8000; % # of sample pts T = 2; % Period of x(t) Wo = 2*pi/T; % Fundamental freq. Ts = 0.001; % Sampling period t = 0:Ts:(N-1)*Ts; % Time vector x_t = zeros(1,N); % Initializing K = 5; % # of terms in FS for k = -K:1:K; W = k*Wo; if k == 0 X_k = 0; else X_k = [1/(i*2*pi*k)]*[exp(-i*2*pi*k) - 2*exp(-i*pi*k) + 1]; end x_t = x_t + X_k*exp(i*W*t); end x_t_real = real(x_t); plot(t, x_t_real); grid on; xlabel('Time') ylabel('Amplitude') title('Approximation of a periodic waveform using exponential FS')
The waveform generated is,
Obviously, this doesn’t appear to be a perfect square wave. There are squiggles at places where we expect a constant amplitude of 1 or -1. Even the transition at the point of discontinuity does not appear to be a perfect vertical line, rather it has some slope, implying a delay such as that of a ramp function. Bear in mind, that an approximation can be as good as our model, and that the chance of approximating a perfect square wave in this case is next to zero. So, the obvious question which might arise is that what can we do to have a better approximation?
The answer to this question is a simple one. All we need to do is consider more terms while approximating the signal. As we increase the number of terms (change K in the script), the approximation appears to get better and better, and this is shown in fig 3.
As we can observe from the above plots, when we increase the number of terms in our FS, the approximation appears to get much better. Although, we still notice that the squiggles at the jump discontinuity still persists. This is something known as Gibbs Phenomenon, which will be discussed in one of the upcoming posts.
Now that we have looked at some of the aspects of the exponential FS, let’s see what’s going on in the MATLAB script shown above.
clear all; clc; clf; N = 8000; % # of sample pts T = 2; % Period of x(t) Wo = 2*pi/T; % Fundamental freq. Ts = 0.001; % Sampling period t = 0:Ts:(N-1)*Ts; % Time vector
This part of the code defines all the required parameters to compute the FS co-efficients and approximate x(t). This part of the code be compared to variable declaration and initialization blocks found in Java or C++. For our approximation, we will be taking N samples of x(t) at every Ts second for the entire length. Hence, we will have T/Ts samples for each period.
x_t = zeros(1,N); % Initializing K = 5; % # of terms in FS for k = -K:1:K; W = k*Wo; if k == 0 X_k = 0; else X_k = [1/(i*2*pi*k)]*[exp(-i*2*pi*k) - 2*exp(-i*pi*k) + 1]; end x_t = x_t + X_k*exp(i*W*t); end x_t_real = real(x_t); plot(t, x_t_real); grid on; xlabel('Time') ylabel('Amplitude') title('Approximation of a periodic waveform using exponential FS')
The declaration of x(t) creates a row vector of N elements, with each element initialized to zero. Let’s take a look at the crazy looking FOR loop. The very first statement in the FOR loop, k = -K:1:K, declares the condition for which the loop is valid. For K = 5, the first iteration of the loop starts at k = -5. At the end of each iteration, the value of k is incremented by 1. At the end of the last iteration, k equals the number of terms in approximation.
One might think, that the approximation is 11 terms and not 5 terms as indicated by the value of K. So what the heck is going on? Just remember that an exponential FS is a complex series, and like all complex quantities, each term in an exponential FS has a complex conjugate. The terms on the negative side are just conjugate components to their counterparts on the positive side. Moreover, the term related to k = 0, can be though of as the DC Offset, with respect to the time axis in this case. If X_k were to be a different than 0 for k = 0, then the waveform would have be simply shifted up for X_k > 0, and shifted down for X_k < 0. For this particular example, I have used 0 DC offset, so that the output is centered around the x-axis.
Obviously, this is not the only we in which we can approximate the periodic waveform. There are plenty of algorithms which one might explore. Rather than deriving the expression for the co-efficient by hand, we can just use the built in MATLAB integration function and embed that in the scripting file.