fseries.m version 01Jul23
(Fourier series)


Author: Paul Mennen
Email:  paul@mennen.org


Overview

In 1807 the French mathematician Joseph Fourier showed that any periodic function can be written as the sum of infinitely many harmonically related sine and cosine waves. This implies that the periodic function can be approximated by a finite sum and that this approximation improves as more terms are added. Later the German mathematician Peter Dirichlet generalized and formalized these results and proved that the infinite series converges. This is a fundamental result that has been extended in many ways that are essential to most areas of signal processing. There are several ways of writing the Fourier series, but one simple way is:


where f(x) is any (piecewise continuous) function that is periodic between -π and π.

The coefficients an and bn may be computed using these integrals:


Your calculus skills must be sharp to use these formulas to compute the an and bn coefficients, but at least for common periodic functions such as the ones used with fseries.m many people have already done this and we can easily look up the results. When you do that, notice that for most simple functions either the an or bn will all be zero which means only one of the sums in the expression for f(x) above is needed.

This program demonstrates the principles of the Fourier series by approximating these seven periodic functions:


A square wave can be approximated by using only the odd harmonics. This trace is the result of summing together 8 terms of the series (i.e. n = 1,3,5,7,9,11,13,15). Notice that the segments that should be flat for a perfect square wave have ripples. These ripples get smaller as we add more terms to the series. Also, note that the rising and falling edges have overshoot (known as the Gibbs phenomenon). The duration of the overshoot (but not the amplitude) gets smaller as we add more terms to the series.


The triangle wave is also approximated with only the odd harmonics. This is also the sum of the first 8 non-zero terms of the series (as are the remaining periodic functions plotted below). This approximation doesn't exhibit an overshoot, but the corners are rounded. These corners become sharper as we add more terms to the series.


The sawtooth waveform requires all the harmonics. This shows the sum of the first 8 terms (n = 1,2,3,4,5,6,7,8).


The pulse (or non-symmetrical square wave) also requires all the harmonics. Here I have chosen the positive portion of the pulse to be 1/6 of the total period.


The quadratic function also requires all the harmonics. Once we make it periodic it doesn't look much like the quadratic function we normally think about, but the usual quadratic (without the repetition) is not periodic and so it can't be represented by a Fourier series.


The absolute value of the sine function is sometimes referred to as a full wave rectified sine wave since this is the signal you get when you pass an AC signal through a full wave rectifier. This function only requires the even harmonics (i.e. n = 2,4,6,8,10,12,14,16).


This function consists of only the positive portion of a sine wave (the maximum value of sin(x) and zero). This is the signal we get when you pass an AC signal through a half-wave rectifier. Unlike all the previous functions, this one is actually the sum of 9 harmonics (n = 1,2,4,6,8,10,12,14,16). The first one (n=1) is a sine function and the remaining eight (even harmonics) are cosine functions.




Here is what the fseries figure window looks like after startup (i.e. with its default settings). The first of seven waveforms is selected (square) showing two cycles of the periodic function with 100 points per cycle (i.e. 200 points total).
The first trace (in green) is the fundamental frequency with trace ID "1st".
The second trace (in purple) has trace ID "3rd" and is the sum of the fundamental frequency plus the 3rd harmonic.
The third trace (in cyan) has trace ID "5th and is the sum of the fundamental and the 3rd and 5th harmonics.
...
The last trace (in white) has trace ID "15th" and is the sum of all the odd harmonics between 1 and 15. It is rendered with a thicker trace so that the best approximation of the square wave is more obvious.

Right click on the Waveform popup to cycle thru all seven periodic functions - or left click on the popup to choose the function you want from a menu. You can also adjust the number of cycles displayed as well as the number of points plotted per cycle. The number of points plotted for each trace is the product of the number of cycles and the number of points per cycle so it ranges from a minimum of 25 points (1 cycle) to a maximum of 25600 points (16 cycles with 1600 points/cycle).

When you click on one of the eight traces, the x and y values of the clicked position appear below the plot in the x and y edit boxes. However, sometimes it is useful to be able to see the y values for all the traces at the same time. The usual way of doing this is to enable the multi-cursor. You can do this by right clicking on the y edit box, and selecting the multiCursor option. You might want to try that here but you will find that the multiCursor is not very useful for this particular application. That's because the 8 traces are all quite similar which means that plt will be placing 8 cursors very close to each other. The cursor values will overlap making them quite difficult to read. But fortunately, plt has a couple of other ways of accomplishing this. The one used here is called the traceid cursor. Note that 'tidcuR' is included in the plt Options list. This enables the traceid cursor which causes the y values of all 8 traces to be displayed in the TraceID box. Try clicking and dragging the cursor left or right and notice that all 8 trace IDs change as you do that. (Also try moving the cursor by left or right clicking on the X axis label ("Cycles") and holding the mouse button down.) Finally, there is one more way of viewing the data for all the traces at once. Click on the "Data" tag in the menu box near the lower left corner of the figure. That will create a new figure that shows 24 successive values for al 8 traces at once starting at the current cursor position. As you move the cursor by either of the methods just mentioned, the new data window will be automatically updated so it always shows the data starting at the cursor location. Also if you want to see more than 24 successive values, you can make the data window taller.


For this figure, I have selected the Sawtooth function and disabled all the traces except the 1st, 2nd, and 8th. Now let's click on the "start" button and observe that the plot starts moving. The amplitude of the selected waveform will vary periodically at a speed that is related to the value entered for the "Speed" parameter (which has been set to 32 in this example). The speed parameter doesn't change how often the display is updated, but rather how much the amplitude changes between every plot update. Below the plot, we see the yellow text "156.50 updates/sec" which tells us how often the plot is updating. There is no pause between display updates, so this number gives us a measure of how fast the computer can update the display. You can use the number to compare the performance of different computers, operating systems, and Matlab versions. The display will continue to update forever (modulating the waveform amplitude in the process) ... until you click on the start button again (which actually will have been renamed to "stop"). This feature of the application (moving plots) complicate the code noticeably, so you may be wondering why this feature was added to the application (especially since the moving plots don't allow us to learn any more about the Fourier series this application is meant to study). The two reasons for the moving plots are one, that it makes the application more visually appealing, and two (and more importantly), it demonstrates some of the techniques involved in creating real-time plots, which are important for many signal processing applications.

Notice that the start/stop button doesn't look like the typical uicontrol button. That's because we have used a pseudo popup control for this, using a mode of the pseudo popup referred to as a "super button". We could have used a standard uicontrol button, but the programming for the pseudo popup control is simpler.

The second new control is the "Speed" popup which is right below the Waveform popup. The speed popup doesn't change the display update rate (although it may seem to). Instead, it changes how much the amplitude is changed at each display update.

The number of display updates per second is calculated and displayed below the plot. This is done every 200 display updates. (It would be a waste of time doing this every update and would also be changing so fast it would be difficult to read.) Every time the updates per second is calculated, the color of the characters alternates between yellow and green to make it obvious that a new value was calculated. When experimenting with the coding of the display loop, you can use this readout to decide if a change improves or degrades the efficiency of the loop. You can also use it to judge the speed differences between various computers, operating systems, and Matlab versions.

Notice that the position coordinates for all the controls have been grouped in one place (lines 26-29). This makes it easy to make adjustments to the GUI with the plt move function. How to do this is described in the GUI building with plt section.

It's easy to make a plot move. Simply create a loop that modifies the 'xData' and/or 'yData' properties of a line (or set of lines). Complications arise when trying to make this idea practical in a typical application. One needs a way of starting and stopping the plot motion. Also one usually wants to change what is happening in the plot when various controls are modified - even while the display loop is in progress. Also, it would be nice for the user to be able to close the figure window while the plot is moving without generating any error messages (Errors will occur if the display loop tries to modify data of a line that no longer exists because the figure was closed.) A final requirement is that the display update should be as quick as possible (i.e. the display loop should be as concise and as efficient as possible). These requirements require some care in how we construct the display loop. This application shows one way of managing those requirements.

The callback function (gCB) responds to changes in any of the controls and also includes the display loop, so, of course, this function is the one that requires the most thought. The callback function behaves quite differently depending on whether the display is static (i.e. in the stopped condition) or dynamic (in the running condition with moving plots). If the display is static, the main display loop (consisting of the 31 lines between line 56 and 88) is run exactly once. On the other hand, if you make any change to one of the controls while the display is dynamic, then the callback function doesn't do anything. (It is called, but it returns right away on line 49.) Instead, we are depending on the main loop which is still running. Eventually, the main loop will recognize that one of the parameters has changed and take appropriate action (by executing the code between lines 61 and 82).

As I mentioned, the display loop consists of 31 lines - however only the last 3 lines (85-87) are executed every time thru the loop. This is done to make the update rate as high as possible. Those 3 lines are sufficient to update the trace values and change the amplitude. The rest of the loop is run only once every 200 display updates and is responsible for detecting that the user has changed a parameter and making the appropriate adjustments based on that change. This is why when you change a parameter while the display is running you may notice a slight delay before the display is affected by that change.

Sometimes it is useful to have a program enter its dynamic (i.e. moving image) state immediately after startup without the need to hit a start button. To provide this option for this application, we add an input argument so that fseries(0) starts the program in its static (non-running) state and fseries(1) starts the program in its dynamic state. Zero is the default option, so calling fseries with no arguments starts in the static state. Also, any character string is treated as nonzero, so something like fseries('go') starts in its dynamic state. One remaining conundrum is how to start the plots moving when the input argument is nonzero. The first idea is as follows:

plt('pop',S.go,'index',go+1);

So when the input argument (go) is zero, the index parameter for the start/stop superbutton is set to 1 (the index for the static state) and when the input argument is 1, the index parameter is set to 2 (the index for the running state). However, this doesn't actually enter the running state because setting the index in this way does not call the pseudo popup's callback function. However, if we use the negative of the index we want, then the callback is called. This leads us to what works (which is also shown in line 42 of the program:

plt('pop',S.go,'index',-go-1); % start running if go is 1

So you may be wondering why this line is commented out in the code. Although that line would work, it does have one drawback. For this command to complete, the callback has to terminate. However since the callback is essentially an infinite loop, it will never terminate (unless we hit the stop button). This means after typing fseries(1) in the command window and pressing the Enter key the >> that we normally see in the command window telling us that we can enter another Matlab command does not appear. This is somewhat inconvenient, since we may want to use the command window for other calculations or to start other applications while fseries.m is still running. The solution to this conundrum is to use the following line, which appears in line 43 of our program:

funcStart({@plt 'pop' S.go 'index' -go-1});

This has the same effect as line 42 except that now the command prompt returns immediately after typing fseries(1) into the command window. More details about the funcStart auxiliary function can be found in the Auxiliary functions section.

Copyright © 2023
Paul Mennen