User Tools
matlab
Differences
This shows you the differences between two versions of the page.
|
matlab [2015/05/15 09:18] schultz [Using the fft command] |
matlab [2023/01/18 16:12] (current) scott Figure scale |
||
|---|---|---|---|
| Line 1: | Line 1: | ||
| - | ====== Matlab Basic Functions ====== | + | ====== MATLAB Basic Functions ====== |
| * len | * len | ||
| * plot | * plot | ||
| Line 11: | Line 11: | ||
| You can type <code> help 'name'</code> to get a description of the particular function. | You can type <code> help 'name'</code> to get a description of the particular function. | ||
| - | Let's do an example to use these basic functions. We are going to plot 4 different signals on 2 separate plots. | + | Let's do an example to use these basic functions. We are going to plot 4 different signals on 2 separate plots. |
| * Plot 1 will contain a square wave and sine waves with amplitude of 2 and frequency of f=1 kHz. | * Plot 1 will contain a square wave and sine waves with amplitude of 2 and frequency of f=1 kHz. | ||
| * Plot 2 will contain a square wave and sine waves with amplitude of 2 and frequency of f=10 kHz. | * Plot 2 will contain a square wave and sine waves with amplitude of 2 and frequency of f=10 kHz. | ||
| - | Here is the Matlab code for the first plot. | + | Here is the MATLAB code for the first plot. |
| <code> | <code> | ||
| f1=1e3; | f1=1e3; | ||
| Line 34: | Line 34: | ||
| </code> | </code> | ||
| - | Here is a description of the Matlab code. | + | Here is a description of the MATLAB code. |
| - We use a variable for the signal frequency (//f1//) and the sampling frequency. | - We use a variable for the signal frequency (//f1//) and the sampling frequency. | ||
| - | - We create a vector for the time (//t//). It has an increment of the sampling time (1/fs) and a total length of 10 periods (//10/f1//). | + | - We create a vector for the time (//t//). It has an increment of the sampling time (1/fs) and a total length of 10 periods (//10/f1//). |
| - We create the vectors for the sine and square waves. | - We create the vectors for the sine and square waves. | ||
| - We use the //subplot// function to plot the signals in the upper half of the plot window. | - We use the //subplot// function to plot the signals in the upper half of the plot window. | ||
| - | - The //hold on// makes it so that the next plot command overlays onto the same plot. //hold off// does the opposite. | + | - The //hold on// makes it so that the next plot command overlays onto the same plot. //hold off// does the opposite. |
| - Use //xlabel// and //ylabel// to add in axes. | - Use //xlabel// and //ylabel// to add in axes. | ||
| - | The second part of the code is similar and just changes the frequency, the time, and the plot location. Here is the resulting Matlab code. | + | The second part of the code is similar and just changes the frequency, the time, and the plot location. Here is the resulting MATLAB code. |
| <code> | <code> | ||
| f2=10e3; | f2=10e3; | ||
| Line 61: | Line 61: | ||
| The resulting plot should look like the following. | The resulting plot should look like the following. | ||
| - | {{ :example1.jpg?100 |}} | + | {{ :example1.jpg?600 |}} |
| - | ====== Matlab Filtering Functions ====== | + | ====== MATLAB Filtering Functions ====== |
| * butter: Butterworth filter design | * butter: Butterworth filter design | ||
| * filter: Time domain filtering of a signal | * filter: Time domain filtering of a signal | ||
| Line 79: | Line 79: | ||
| </code> | </code> | ||
| - | on the Matlab command line. Read through the description. Let’s look at parts of the description and what it means. | + | on the MATLAB command line. Read through the description. Let’s look at parts of the description and what it means. |
| **Description**: “butter Butterworth digital and analog filter design.” | **Description**: “butter Butterworth digital and analog filter design.” | ||
| Line 103: | Line 103: | ||
| </code> | </code> | ||
| - | on the Matlab command line. Read through the description. Let’s look at parts of the description and what it means. | + | on the MATLAB command line. Read through the description. Let’s look at parts of the description and what it means. |
| **Description**: “One-dimensional digital filter.” | **Description**: “One-dimensional digital filter.” | ||
| Line 132: | Line 132: | ||
| You should get a plot that looks like | You should get a plot that looks like | ||
| - | {{ :matlab_help_fig1.jpg?100 |}} | + | {{ :matlab_help_fig1.jpg?600 |}} |
| - | If we create a low pass filter with a cutoff frequency of //fcut=1.5kHz// then the signal should pass through the filter. Let’s make some mistakes and see what happens and see if we can figure out what we did wrong. We are doing this to learn about debugging. | + | If we create a low pass filter with a cutoff frequency of //fcut=1.5kHz// then the signal should pass through the filter. Let’s make some mistakes and see what happens and see if we can figure out what we did wrong. We are doing this to learn about debugging. |
| **Mistake**: Pass the cutoff frequency as the actual frequency. | **Mistake**: Pass the cutoff frequency as the actual frequency. | ||
| Line 162: | Line 162: | ||
| </code> | </code> | ||
| The resulting plot should look something like. | The resulting plot should look something like. | ||
| - | {{ :matlab_help_fig2.jpg?100 |}} | + | {{ :matlab_help_fig2.jpg?600 |}} |
| - | The signal y has a lower amplitude. This mean that it is being filtered. We must have made an error in the cutoff frequency. If we make the following change | + | The signal y has a lower amplitude. This mean that it is being filtered. We must have made an error in the cutoff frequency. If we make the following change |
| <code> | <code> | ||
| Wn=Wn/(2*Fs); | Wn=Wn/(2*Fs); | ||
| </code> | </code> | ||
| The following plot shows that it is now worse because we were supposed to divide by Fs/2 not 2*Fs. | The following plot shows that it is now worse because we were supposed to divide by Fs/2 not 2*Fs. | ||
| - | {{ :matlab_help_fig3.jpg?100 |}} | + | {{ :matlab_help_fig3.jpg?600 |}} |
| The correct value is | The correct value is | ||
| Line 177: | Line 177: | ||
| Resulting in the correct plot. | Resulting in the correct plot. | ||
| - | {{ :matlab_help_fig4.jpg?100 |}} | + | {{ :matlab_help_fig4.jpg?600 |}} |
| - | To finish up let’s make sure that we can use the butter command to perform a bandpass filter. Create a signal that is x=sin(2*pi*1e3*t) + sin(2*pi*2e3*t) + sin(2*pi*3e3*t) and has the same sampling frequency of Fs=10kHz. Use a bandpass filter to pass the f=2kHz signal through and attenuate the other two. The resulting plots should look like the following. | + | To finish up let’s make sure that we can use the butter command to perform a bandpass filter. Create a signal that is x=sin(2*pi*1e3*t) + sin(2*pi*2e3*t) + sin(2*pi*3e3*t) and has the same sampling frequency of Fs=10kHz. Use a bandpass filter to pass the f=2kHz signal through and attenuate the other two. The resulting plots should look like the following. |
| - | {{ :matlab_help_fig5.jpg?100 |}} | + | {{ :matlab_help_fig5.jpg?600 |}} |
| ===== Using the fft command ===== | ===== Using the fft command ===== | ||
| - | MATLAB has a built-in function for performing numerical Fourier transforms called the fast Fourier transform (FFT). The FFT is actually a discrete Fourier transform (DFT). We can relate the DFT to the continuous time domain. This [[https://www.youtube.com/watch?v=z7X6jgFnB6Y| youtube video]] explains the conversion between the FFT and frequency. MATLAB also has some help on using the FFT. In the command window type | + | MATLAB has a built-in function for performing numerical Fourier transforms called the fast Fourier transform (FFT). The FFT is actually a discrete Fourier transform (DFT). We can relate the DFT to the continuous time domain. This [[https://www.youtube.com/watch?v=z7X6jgFnB6Y| youtube video]] explains the conversion between the FFT and frequency. MATLAB also has some help on using the FFT. In the command window type |
| <code> help fft</code> | <code> help fft</code> | ||
| Then click on the | Then click on the | ||
| Line 194: | Line 194: | ||
| Let's go through another example besides the one provided in the MATLAB help to make sure that you understand how to use the FFT. | Let's go through another example besides the one provided in the MATLAB help to make sure that you understand how to use the FFT. | ||
| - | You need to create a frequency domain plot of a sinusoidal signal with frequency of f1=3 kHz and amplitude of 1. We are using a sinusoidal signal because we know that the Fourier transform of a sine wave is a delta function. So the FFT should produce a spike at a frequency of 3 kHz and have an amplitude of 1. | + | You need to create a frequency domain plot of a sinusoidal signal with frequency of f1=3 kHz and amplitude of 1. We are using a sinusoidal signal because we know that the Fourier transform of a sine wave is a delta function. So the FFT should produce a spike at a frequency of 3 kHz and have an amplitude of 1. |
| Here is the MATLAB code to create the time signal | Here is the MATLAB code to create the time signal | ||
| Line 206: | Line 206: | ||
| Here is the plot of our signal. | Here is the plot of our signal. | ||
| - | {{ :fft_example1.jpg?100 |}} | + | {{ :fft_example1.jpg?600 |}} |
| Now you need to create the Fourier transform of the signal using the FFT command. | Now you need to create the Fourier transform of the signal using the FFT command. | ||
| Line 213: | Line 213: | ||
| plot(abs(Y)) | plot(abs(Y)) | ||
| </code> | </code> | ||
| - | Notice that we plotted the absolute value of Y since it is complex. | + | Notice that we plotted the absolute value of Y since it is complex. |
| - | {{ :fft_example2.jpg?100 |}} | + | {{ :fft_example2.jpg?600 |}} |
| Notice that there are two spikes and that the amplitude is not equal to our amplitude of 1. | Notice that there are two spikes and that the amplitude is not equal to our amplitude of 1. | ||
| - | Let's start with the frequency. We only want the positive frequencies so we are only going to use the first half of the vector. We know that maximum frequency of a sampled signal is Fs/2. So we create an x-axis vector that goes from 0 to Fs/2 and has a length of half of Y. Here is the resulting MATLAB code. | + | Let's start with the frequency. We only want the positive frequencies so we are only going to use the first half of the vector. We know that maximum frequency of a sampled signal is Fs/2. So we create an x-axis vector that goes from 0 to Fs/2 and has a length of half of Y. Here is the resulting MATLAB code. |
| <code> | <code> | ||
| L=length(Y); %the total length of our fft | L=length(Y); %the total length of our fft | ||
| Line 226: | Line 226: | ||
| plot(f,abs(Y1)) | plot(f,abs(Y1)) | ||
| </code> | </code> | ||
| - | The linspace command creates a vector with equally spaced elements between a starting and an ending value and with a specified number of elements. Here is the resulting plot. | + | The linspace command creates a vector with equally spaced elements between a starting and an ending value and with a specified number of elements. Here is the resulting plot. |
| - | {{ :fft_example3.jpg?100 |}} | + | {{ :fft_example3.jpg?600 |}} |
| You can zoom in on the spike and verify that it has a frequency close to 3kHz. | You can zoom in on the spike and verify that it has a frequency close to 3kHz. | ||
| - | Our signal does not have an amplitude greater than 100 so we obviously need to scale the amplitude. We simply need to divide the signal by the number of elements. | + | Our signal does not have an amplitude greater than 100 so we obviously need to scale the amplitude. We simply need to divide the signal by the number of elements. |
| <code> | <code> | ||
| Y2=Y1/(L/2); | Y2=Y1/(L/2); | ||
| Line 236: | Line 236: | ||
| </code> | </code> | ||
| Here is the resulting plot. | Here is the resulting plot. | ||
| - | {{ :fft_example4.jpg?100 |}} | + | {{ :fft_example4.jpg?600 |}} |
| Zoom in on the spike to make sure that it is centered close to 3kHz. | Zoom in on the spike to make sure that it is centered close to 3kHz. | ||
| - | //NOTE: As you learned by watching the [[https://www.youtube.com/watch?v=z7X6jgFnB6Y| youtube video]], the FFT creates an array of frequency domain spectrum samples. The number of samples is the same as the length of your data (N) and the frequency range is based on the sampling frequency fmax=Fs/2. This means that the frequency bin is Fs/(2*N).// | + | //NOTE: As you learned by watching the [[https://www.youtube.com/watch?v=z7X6jgFnB6Y| youtube video]], the FFT creates an array of frequency domain spectrum samples. The number of samples is the same as the length of your data (N) and the frequency range is based on the sampling frequency fmax=Fs/2. This means that the frequency bin is Fs/(2*N).// |
| - | In our example the number of samples is L/2=167 and the frequency band is Fs/2=50kHz. This means that our frequency resolution is 299.4 Hz. If we increase the number of points in our time domain signal without increasing the total time of the signal then we are actual changing the sampling frequency, which would increase our total frequency band but not change our frequency resolution. In order to change the resolution we need to increase the total time width of our signal without changing the sampling frequency. | + | In our example the number of samples is L/2=167 and the frequency band is Fs/2=50kHz. This means that our frequency resolution is 299.4 Hz. If we increase the number of points in our time domain signal without increasing the total time of the signal then we are actual changing the sampling frequency, which would increase our total frequency band but not change our frequency resolution. In order to change the resolution we need to increase the total time width of our signal without changing the sampling frequency. |
| This MATLAB code produces a sine wave with 10 periods and a corresponding frequency resolution of 299 Hz. | This MATLAB code produces a sine wave with 10 periods and a corresponding frequency resolution of 299 Hz. | ||
| Line 262: | Line 262: | ||
| </code> | </code> | ||
| Here are the resulting time domain and frequency domain plots. | Here are the resulting time domain and frequency domain plots. | ||
| - | {{ :fft_example5.jpg?100 |}} | + | {{ :fft_example5.jpg?600 |}} |
| This MATLAB code produces a sine wave with 100 periods and a corresponding frequency resolution of 29.9 Hz. | This MATLAB code produces a sine wave with 100 periods and a corresponding frequency resolution of 29.9 Hz. | ||
| Line 282: | Line 282: | ||
| </code> | </code> | ||
| Here are the resulting plots. | Here are the resulting plots. | ||
| - | {{ :fft_example6.jpg?100 |}} | + | {{ :fft_example6.jpg?600 |}} |
| Line 292: | Line 292: | ||
| **Description**: “freqz Frequency response of digital filter.” | **Description**: “freqz Frequency response of digital filter.” | ||
| - | **Why we care**: We want the frequency response of a digital filter. So it looks like we found the correct MATLAB command. | + | **Why we care**: We want the frequency response of a digital filter. So it looks like we found the correct MATLAB command. |
| **Description**: “[H,W] = freqz(B,A,N) returns the N-point complex frequency response vector H and the N-point frequency vector W in radians/sample of the filter:” | **Description**: “[H,W] = freqz(B,A,N) returns the N-point complex frequency response vector H and the N-point frequency vector W in radians/sample of the filter:” | ||
| - | **Why we care**: H is the complex frequency response. N is the number of points that are calculated. | + | **Why we care**: H is the complex frequency response. N is the number of points that are calculated. |
| - | We don't really want the plot in terms of radians/sample. We would prefer to create the plot in terms of frequency. So we are going to scan down the description until we find the form of freqz that is in terms of frequency. | + | We don't really want the plot in terms of radians/sample. We would prefer to create the plot in terms of frequency. So we are going to scan down the description until we find the form of freqz that is in terms of frequency. |
| - | **Description**:"H = freqz(...,F,Fs) returns the complex frequency response at the frequencies designated in vector F (in Hz), where Fs is the sampling frequency (in Hz)." | + | **Description**: "H = freqz(...,F,Fs) returns the complex frequency response at the frequencies designated in vector F (in Hz), where Fs is the sampling frequency (in Hz)." |
| - | **Why we care**: This tells us how to create the plot in terms of frequency rather than samples. | + | **Why we care**: This tells us how to create the plot in terms of frequency rather than samples. |
| Here is an example to show how to use the freqz command. | Here is an example to show how to use the freqz command. | ||
| Line 332: | Line 332: | ||
| Here is the resulting plot. | Here is the resulting plot. | ||
| - | {{ :freqz_1.jpg?200 |}} | + | {{ :freqz_1.jpg?600 |}} |
| - | Notice how the filter is flat over the band between the two frequencies and then slopes down outside of this band. This is the shape of the frequency response of a bandpass filter. | + | Notice how the filter is flat over the band between the two frequencies and then slopes down outside of this band. This is the shape of the frequency response of a bandpass filter. |
| - | Often the data is scaled to be in decibels so that it is easier to distinguish small changes. When you do this you often need to scale the y-axis of the plot to your desired range of interest. | + | Often the data is scaled to be in decibels so that it is easier to distinguish small changes. When you do this you often need to scale the y-axis of the plot to your desired range of interest. |
| <code> | <code> | ||
| H_db=20*log10(abs(H)); | H_db=20*log10(abs(H)); | ||
| Line 342: | Line 342: | ||
| axis([0 5000 -50 3]) | axis([0 5000 -50 3]) | ||
| </code> | </code> | ||
| - | Here is the plot in decibel scale. You can zoom in on the 1kHz and 2kHz frequencies and see that these are the 3dB points. | + | Here is the plot in decibel scale. You can zoom in on the 1kHz and 2kHz frequencies and see that these are the 3dB points. |
| - | {{ :freqz_3.jpg?200 |}} | + | {{ :freqz_3.jpg?600 |}} |
| You can change the order of the filter and see that the filter becomes steeper and change the frequency corners and see the change in the width of the filter. | You can change the order of the filter and see that the filter becomes steeper and change the frequency corners and see the change in the width of the filter. | ||
matlab.1431703137.txt.gz · Last modified: 2015/05/15 09:18 by schultz
Page Tools
Except where otherwise noted, content on this wiki is licensed under the following license: CC Attribution-Share Alike 3.0 Unported
