Skip to content
Snippets Groups Projects
Commit 05aa9ad8 authored by Robert Lanzafame's avatar Robert Lanzafame
Browse files

WS 2.3 ready

parent b441dfa2
No related branches found
No related tags found
1 merge request!25Main
Pipeline #259471 passed
Stage: deploy
Expand all Hide whitespace changes
Inline Side-by-side
content/Week_2_3/WS_2_3_DFT_you_try_meow.ipynb
+ 10
− 10
View file @ 05aa9ad8
%% Cell type:markdown id:4af6b45f-f7aa-46a8-9d79-c1579ae0a799 tags:
%% Cell type:markdown id:81475e62 tags:
# WS11: Discrete Fourier Transform (DFT): You Try Meow (Miauw)
# WS 2.3: Discrete Fourier Transform (DFT): You Try Meow (Miauw)
<h1 style="position: absolute; display: flex; flex-grow: 0; flex-shrink: 0; flex-direction: row-reverse; top: 60px;right: 30px; margin: 0; border: 0">
<style>
.markdown {width:100%; position: relative}
article { position: relative }
</style>
<img src="https://gitlab.tudelft.nl/mude/public/-/raw/main/tu-logo/TU_P1_full-color.png" style="width:100px" />
<img src="https://gitlab.tudelft.nl/mude/public/-/raw/main/mude-logo/MUDE_Logo-small.png" style="width:100px" />
</h1>
<h2 style="height: 10px">
</h2>
*[CEGM1000 MUDE](http://mude.citg.tudelft.nl/): Week 2.3, Signal Processing. For: November 29, 2023*
*[CEGM1000 MUDE](http://mude.citg.tudelft.nl/): Week 2.3, Signal Processing. For: November 27, 2024*
%% Cell type:markdown id:eff9791a tags:
The goal of this workshop to work with the _Discrete Fourier Transform_ (DFT), implemented in Python as the _Fast Fourier Transform_ (FFT) through `np.fft.fft`, and to understand and interpret its output.
The notebook consists of two parts:
- The first part (Task 0) is a demonstration of the use of the DFT (_you read and execute the code cells_),
- The second part is a simple exercise with the DFT (_you write the code_).
To start off, let's do a quick quiz question: _what is the primary purpose of the DFT?_
_Run the code cell below to find the answer._
%% Cell type:code id:36a05bbd-69a1-4805-a4ad-e84de89de560 tags:
``` python
%%html
<iframe src="https://tudelft.h5p.com/content/1292126914399042257/embed" aria-label="Meow" width="1088" height="637" frameborder="0" allowfullscreen="allowfullscreen" allow="autoplay *; geolocation *; microphone *; camera *; midi *; encrypted-media *"></iframe><script src="https://tudelft.h5p.com/js/h5p-resizer.js" charset="UTF-8"></script>
```
%% Cell type:markdown id:03f3b3c6-d3f6-4121-9f39-c30c8bf70ea5 tags:
That's right! We convert our signal into the frequency domain.
And if you would like an additional explanation of the key frequencies, you can find it [here](https://medium.com/@kovalenko.alx/fun-with-fourier-591662576a77).
%% Cell type:markdown id:bfca25a5-75ce-4837-9387-01f95be10bd0 tags:
<div style="background-color:#facb8e; color: black; vertical-align: middle; padding:15px; margin: 10px; border-radius: 10px; width: 95%"> <p>Note the use of <code>zip</code>, <code>stem</code>, <code>annotate</code> and the modulo operator
<code>%</code>. Refer to PA11 if you do not understand these tools. Furthermore, note that the term _modulus_ is also used here (and in the textbook), which is another term for _absolute value._</p></div>
<code>%</code>. Refer to PA 2.3 if you do not understand these tools. Furthermore, note that the term _modulus_ is also used here (and in the textbook), which is another term for _absolute value._</p></div>
%% Cell type:code id:da13fbf3 tags:
``` python
import numpy as np
from matplotlib import pyplot as plt
```
%% Cell type:markdown id:13234855 tags:
## Task 0: Demonstration of DFT using pulse function
In the first part of this notebook, we use $x(t)=\Pi(\frac{t}{4})$, and its Fourier transform $X(f)=4 \,\textrm{sinc}(4f)$, as an example (see the first worked example in Chapter 3 on the Fourier transform). The pulse lasts for 4 seconds in the time domain; for convenience, below it is not centered at $t$ = 0, but shifted (delayed) to the right.
The signal $x(t)$ clearly is non-periodic and is an energy signal; apart from a short time span of 'activity' it is zero elsewhere.
We create a pulse function $x(t)$ in discrete time $x_n$ by numpy. The total signal duration is $T$ = 20 seconds (observation or record length). The sampling interval is $\Delta t$ = 1 second. There are $N$ = 20 samples, and each sample represents 1 second, hence $N \Delta t = T$.
Note that the time array starts at $t$ = 0 s, and hence, the last sample is at $t$ = 19 s (and not at $t$ = 20 s, as then we would have 21 samples).
### Task 0.1: Visualize the Signal
%% Cell type:code id:b33e1be2 tags:
``` python
t = np.arange(0,20,1)
xt = np.concatenate((np.zeros(8), np.ones(4), np.zeros(8)))
plt.plot(t, xt,'o')
plt.stem(t[8:12], xt[8:12])
plt.xticks(ticks=np.arange(0,21,5), labels=np.arange(0,21,5))
plt.xlabel('time [s]')
plt.ylabel('xn');
```
%% Cell type:markdown id:a338eede tags:
### Task 0.2: Evaluate (and visualize) the DFT
We use `numpy.fft.fft` to compute the one-dimensional Discrete Fourier Transform (DFT) of $x_n$, which takes a signal as argument and returns an array of coefficients (refer to the [documentation](https://numpy.org/doc/stable/reference/generated/numpy.fft.fft.html) as needed).
The DFT converts $N$ samples of the time domain signal $x_n$, into $N$ samples of the frequency domain. In this case, it produces $X_k$ with $k = 0, 1, 2, \ldots, N-1$, with $N$ = 20, which are complex numbers.
%% Cell type:markdown id:26d0f558-caa7-41f2-9477-60d74940fc10 tags:
<div style="background-color:#AABAB2; color: black; vertical-align: middle; padding:15px; margin: 10px; border-radius: 10px; width: 95%">
<p>
<b>Task 0.2:</b>
Read the code cell below before executing it and identify the following:
<ol>
<li>Where is the DFT computed, and what is the output?</li>
<li>Why is the modulus (absolute value) used on the DFT output?</li>
<li>What are the values are used for the x and y axes of the plot?</li>
<li>How is frequency information found and added to the plot (mathematically)?</li>
</ol>
Once you understand the figure, continue reading to understand <em>what do these 20 complex numbers mean, and how should we interpret them?</em>
</p>
</div>
%% Cell type:code id:3acb465b tags:
``` python
abs_fft = np.abs(np.fft.fft(xt))
index_fft = np.arange(0,20,1)
plt.plot(index_fft, abs_fft, 'o')
freq = np.arange(0, 1, 0.05)
for x,y in zip(index_fft, abs_fft):
if x%5 == 0 or x==19:
label = f"f={freq[x]:.2f} Hz"
plt.annotate(label,
(x,y),
textcoords="offset points",
xytext=(0,10),
fontsize=10,
ha='center')
plt.ylim(0,5)
plt.xlim(-2,21)
plt.xlabel('fft-index')
plt.ylabel('$|X_k|$')
plt.stem(index_fft, abs_fft);
```
%% Cell type:markdown id:09c72cdd tags:
The frequency resolution $\Delta f$ equals one-over-the-measurement-duration, hence $\Delta f = 1/T$. With that knowledge, we can reconstruct the frequencies expressed in Hertz.
The spectrum of the sampled signal is periodic in the sampling frequency $f_s$ which equals 1 Hz ($\Delta t$ = 1 s). Therefore, it is computed just for one period $[0,f_s)$. The last value, with index 19, represents the component with frequency $f$ = 0.95 Hz. A spectrum is commonly presented and interpreted as double-sided, so in the above graph we can interpret the spectral components with indices 10 to 19 corresponding to negative frequencies.
### Task 0.3: Identify negative frequencies
%% Cell type:code id:df1f483b tags:
``` python
abs_fft = np.abs(np.fft.fft(xt))
plt.stem(index_fft, abs_fft)
plt.plot(index_fft, abs_fft, 'o')
freq = np.concatenate((np.arange(0, 0.5, 0.05), np.arange(-0.5, 0, 0.05)))
for x,y in zip(index_fft, abs_fft):
if x%5 == 0 or x==19:
label = f"f={freq[x]:.2f} Hz"
plt.annotate(label,
(x,y),
textcoords="offset points",
xytext=(0,10),
fontsize=10,
ha='center')
plt.ylim(0,5)
plt.xlim(-2,21)
plt.xlabel('fft-index')
plt.ylabel('$|X_k|$');
```
%% Cell type:markdown id:2d87a9a6 tags:
Now we can interpret the DFT of the $x(t)=\Pi(\frac{t}{4})$. The sampling interval is $\Delta t$ = 1 second, and the obsesrvation length is T=20 seconds, and we have N=20 samples.
- We can recognize a bit of a sinc function.
- The DFT is computed from 0 Hz, for positive frequencies up to $\frac{fs}{2} = 0.5$ Hz, after which the negative frequencies follow from -0.5 to -0.05 Hz.
- $X(f)=4 \textrm{sinc}(4f)$ has its first null at $f=0.25$Hz.
### Task 0.4: Create symmetric plot
For convenient visualization, we may want to explicitly shift the negative frequencies to the left-hand side to create a symmetric plot. We can use `numpy.fft.fftshift` to do that. In other words, the zero-frequency component appears in the center of the spectrum. Now, it nicely shows a symmetric sprectrum. It well resembles the sinc-function (taking into account that we plot the modulus/absolute value). The output of the DFT still consists of $N$ = 20 elements. To enable this, we set up a new frequency array (in Hz) on the interval $[-fs/2,fs/2)$.
%% Cell type:markdown id:e1eef6b6-1e90-43db-b059-6521b2a840d3 tags:
<div style="background-color:#AABAB2; color: black; vertical-align: middle; padding:15px; margin: 10px; border-radius: 10px; width: 95%">
<p>
<b>Task 0.4:</b>
Read the code cell below before executing it and identify how the plot is modified based on the (new) specification of frequency. Note that it is more than just the <code>freq</code> variable!
</p>
</div>
%% Cell type:code id:a8aac894 tags:
``` python
abs_fft_shift = np.abs(np.fft.fftshift(np.fft.fft(xt)))
freq = np.arange(-0.5, 0.5, 0.05)
plt.stem(freq, abs_fft_shift)
plt.plot(freq, abs_fft_shift, 'o')
plt.ylabel('|Xk|')
plt.xlabel('frequency [Hz]');
```
%% Cell type:markdown id:569ca2d3 tags:
### Task 0.5: Showing spectrum only for positive frequencies
In practice, because of the symmetry, one typically plots only the right-hand side of the (double-sided) spectrum, hence only the part for positive frequencies $f \geq 0$. This is simply a matter of preference, and a way to save some space.
%% Cell type:markdown id:a26e361d-c5fa-47a2-b32c-416cc3d10659 tags:
<div style="background-color:#AABAB2; color: black; vertical-align: middle; padding:15px; margin: 10px; border-radius: 10px; width: 95%">
<p>
<b>Task 0.5:</b>
Can you identify what has changed (in the code and visually), compared to the previous plot?
</p>
</div>
%% Cell type:code id:b9608638 tags:
``` python
N=len(xt)
abs_fft = np.abs(np.fft.fft(xt))
freq = np.arange(0.0, 1.0, 0.05)
plt.plot(freq[:int(N/2)], abs_fft[:int(N/2)], 'o')
plt.stem(freq[:int(N/2)], abs_fft[:int(N/2)])
plt.ylabel('$|X_k|$')
plt.xlabel('frequency [Hz]');
```
%% Cell type:markdown id:83c4ca67-9234-4741-a56f-294a12bcd98a tags:
<div style="background-color:#AABAB2; color: black; vertical-align: middle; padding:15px; margin: 10px; border-radius: 10px; width: 95%">
<p>
<b>Task 0.6:</b>
Confirm that you understand how we have arrived at the plot above, which illustrates the magnitude (amplitude) spectrum for frequencies $f \in [0,f_s/2)$, rather than $[0,f_s)$.
</p>
</div>
%% Cell type:markdown id:9f5558b4 tags:
## Task 1: Application of DFT using simple cosine
%% Cell type:markdown id:46dd5c0c tags:
It is always a good idea, in spectral analysis, to run a test with a very simple, basic signal. In this way you can test and verify your coding and interpretation of the results.
Our basic signal is just a plain cosine. We take the amplitude equal to one, and zero initial phase, so the signal reads $x(t) = \cos(2 \pi f_c t)$, with $f_c$ = 3 Hz in this exercise. With such a simple signal, we know in advance how the spectrum should look like. Namely just a spike at $f$ = 3 Hz, and also one at $f$ = -3 Hz, as we're, for mathematical convenience, working with double sided spectra. The spectrum should be zero at all other frequencies.
As a side note: the cosine is strictly a periodic function, not a-periodic (as above); still, the Fourier transform of the cosine is defined as two Dirac delta pulses or peaks (at 3 Hz and -3 Hz). You may want to check out the second worked example in Chapter 3 on the Fourier transform: Fourier transform in the limit.
%% Cell type:markdown id:87de3cc5 tags:
<div style="background-color:#AABAB2; color: black; vertical-align: middle; padding:15px; margin: 10px; border-radius: 10px; width: 95%">
<p>
<b>Task 1:</b>
Create a sampled (discrete time) cosine signal by sampling at $f_s$ = 10 Hz, for a duration of $T$ = 2 seconds (make sure you use exactly $N$ = 20 samples). Plot the sampled signal, compute its DFT and plot its magnitude spectrum $|X_k|$ with proper labeling of the axes (just like we did in the first part of this notebook). Include a plot of the spectrum of the sampled cosine signal using only the positive frequencies (up to $f_s/2$, as in the last plot of the previous task).
<em>Note: you are expected to produce three separate plots.</em>
</p>
</div>
%% Cell type:code id:28ad034e tags:
``` python
YOUR_CODE_HERE_PLOT_1
```
%% Cell type:code id:af87d120 tags:
``` python
YOUR_CODE_HERE_PLOT_2
```
%% Cell type:code id:6fd9014c tags:
``` python
YOUR_CODE_HERE_PLOT_3
```
%% Cell type:markdown id:6a4a306d-a372-460e-b312-e8fa138a9797 tags:
%% Cell type:markdown id:0b3584c9 tags:
**End of notebook.**
<h2 style="height: 60px">
</h2>
<h3 style="position: absolute; display: flex; flex-grow: 0; flex-shrink: 0; flex-direction: row-reverse; bottom: 60px; right: 50px; margin: 0; border: 0">
<style>
.markdown {width:100%; position: relative}
article { position: relative }
</style>
<a rel="license" href="http://creativecommons.org/licenses/by-nc-sa/4.0/">
<img alt="Creative Commons License" style="border-width:; width:88px; height:auto; padding-top:10px" src="https://i.creativecommons.org/l/by-nc-sa/4.0/88x31.png" />
<a rel="license" href="http://creativecommons.org/licenses/by/4.0/">
<img alt="Creative Commons License" style="border-width:; width:88px; height:auto; padding-top:10px" src="https://i.creativecommons.org/l/by/4.0/88x31.png" />
</a>
<a rel="TU Delft" href="https://www.tudelft.nl/en/ceg">
<img alt="TU Delft" style="border-width:0; width:100px; height:auto; padding-bottom:0px" src="https://gitlab.tudelft.nl/mude/public/-/raw/main/tu-logo/TU_P1_full-color.png"/>
<img alt="TU Delft" style="border-width:0; width:100px; height:auto; padding-bottom:0px" src="https://gitlab.tudelft.nl/mude/public/-/raw/main/tu-logo/TU_P1_full-color.png" />
</a>
<a rel="MUDE" href="http://mude.citg.tudelft.nl/">
<img alt="MUDE" style="border-width:0; width:100px; height:auto; padding-bottom:0px" src="https://gitlab.tudelft.nl/mude/public/-/raw/main/mude-logo/MUDE_Logo-small.png"/>
<img alt="MUDE" style="border-width:0; width:100px; height:auto; padding-bottom:0px" src="https://gitlab.tudelft.nl/mude/public/-/raw/main/mude-logo/MUDE_Logo-small.png" />
</a>
</h3>
<span style="font-size: 75%">
&copy; Copyright 2023 <a rel="MUDE Team" href="https://studiegids.tudelft.nl/a101_displayCourse.do?course_id=65595">MUDE Teaching Team</a> TU Delft. This work is licensed under a <a rel="license" href="http://creativecommons.org/licenses/by-nc-sa/4.0/">Creative Commons Attribution-NonCommercial-ShareAlike 4.0 International License</a>.
&copy; Copyright 2024 <a rel="MUDE" href="http://mude.citg.tudelft.nl/">MUDE</a> TU Delft. This work is licensed under a <a rel="license" href="http://creativecommons.org/licenses/by/4.0/">CC BY 4.0 License</a>.
......
content/Week_2_3/WS_2_3_solution.ipynb
+ 8
− 8
View file @ 05aa9ad8
%% Cell type:markdown id:4af6b45f-f7aa-46a8-9d79-c1579ae0a799 tags:
# WS11: Discrete Fourier Transform (DFT): You Try Meow (Miauw)
# WS 2.3: Discrete Fourier Transform (DFT): You Try Meow (Miauw)
<h1 style="position: absolute; display: flex; flex-grow: 0; flex-shrink: 0; flex-direction: row-reverse; top: 60px;right: 30px; margin: 0; border: 0">
<style>
.markdown {width:100%; position: relative}
article { position: relative }
</style>
<img src="https://gitlab.tudelft.nl/mude/public/-/raw/main/tu-logo/TU_P1_full-color.png" style="width:100px" />
<img src="https://gitlab.tudelft.nl/mude/public/-/raw/main/mude-logo/MUDE_Logo-small.png" style="width:100px" />
</h1>
<h2 style="height: 10px">
</h2>
*[CEGM1000 MUDE](http://mude.citg.tudelft.nl/): Week 2.3, Signal Processing. For: November 29, 2023*
*[CEGM1000 MUDE](http://mude.citg.tudelft.nl/): Week 2.3, Signal Processing. For: November 27, 2024*
%% Cell type:markdown id:eff9791a tags:
The goal of this workshop to work with the _Discrete Fourier Transform_ (DFT), implemented in Python as the _Fast Fourier Transform_ (FFT) through `np.fft.fft`, and to understand and interpret its output.
The notebook consists of two parts:
- The first part (Task 0) is a demonstration of the use of the DFT (_you read and execute the code cells_),
- The second part is a simple exercise with the DFT (_you write the code_).
To start off, let's do a quick quiz question: _what is the primary purpose of the DFT?_
_Find the answer [here](https://tudelft.h5p.com/content/1292126914399042257/embed)._
%% Cell type:markdown id:03f3b3c6-d3f6-4121-9f39-c30c8bf70ea5 tags:
That's right! We convert our signal into the frequency domain. And if you would like an additional explanation of the key frequencies, you can find it [here](https://medium.com/@kovalenko.alx/fun-with-fourier-591662576a77).
%% Cell type:markdown id:bfca25a5-75ce-4837-9387-01f95be10bd0 tags:
<div style="background-color:#facb8e; color: black; vertical-align: middle; padding:15px; margin: 10px; border-radius: 10px; width: 95%"> <p>Note the use of <code>zip</code>, <code>stem</code>, <code>annotate</code> and the modulo operator
<code>%</code>. Refer to PA11 if you do not understand these tools. Furthermore, note that the term _modulus_ is also used here (and in the textbook), which is another term for _absolute value._</p></div>
%% Cell type:code id:da13fbf3 tags:
``` python
import numpy as np
from matplotlib import pyplot as plt
```
%% Cell type:markdown id:13234855 tags:
## Task 0: Demonstration of DFT using pulse function
In the first part of this notebook, we use $x(t)=\Pi(\frac{t}{4})$, and its Fourier transform $X(f)=4 \,\textrm{sinc}(4f)$, as an example (see the first worked example in Chapter 3 on the Fourier transform). The pulse lasts for 4 seconds in the time domain; for convenience, below it is not centered at $t$ = 0, but shifted (delayed) to the right.
The signal $x(t)$ clearly is non-periodic and is an energy signal; apart from a short time span of 'activity' it is zero elsewhere.
We create a pulse function $x(t)$ in discrete time $x_n$ by numpy. The total signal duration is $T$ = 20 seconds (observation or record length). The sampling interval is $\Delta t$ = 1 second. There are $N$ = 20 samples, and each sample represents 1 second, hence $N \Delta t = T$.
Note that the time array starts at $t$ = 0 s, and hence, the last sample is at $t$ = 19 s (and not at $t$ = 20 s, as then we would have 21 samples).
### Task 0.1: Visualize the Signal
%% Cell type:code id:b33e1be2 tags:
``` python
t = np.arange(0,20,1)
xt = np.concatenate((np.zeros(8), np.ones(4), np.zeros(8)))
plt.plot(t, xt,'o')
plt.stem(t[8:12], xt[8:12])
plt.xticks(ticks=np.arange(0,21,5), labels=np.arange(0,21,5))
plt.xlabel('time [s]')
plt.ylabel('xn');
```
%% Output
%% Cell type:markdown id:a338eede tags:
### Task 0.2: Evaluate (and visualize) the DFT
We use `numpy.fft.fft` to compute the one-dimensional Discrete Fourier Transform (DFT) of $x_n$, which takes a signal as argument and returns an array of coefficients (refer to the [documentation](https://numpy.org/doc/stable/reference/generated/numpy.fft.fft.html) as needed).
The DFT converts $N$ samples of the time domain signal $x_n$, into $N$ samples of the frequency domain. In this case, it produces $X_k$ with $k = 0, 1, 2, \ldots, N-1$, with $N$ = 20, which are complex numbers.
%% Cell type:markdown id:26d0f558-caa7-41f2-9477-60d74940fc10 tags:
<div style="background-color:#AABAB2; color: black; vertical-align: middle; padding:15px; margin: 10px; border-radius: 10px; width: 95%">
<p>
<b>Task 0.2:</b>
Read the code cell below before executing it and identify the following:
<ol>
<li>Where is the DFT computed, and what is the output?</li>
<li>Why is the modulus (absolute value) used on the DFT output?</li>
<li>What are the values are used for the x and y axes of the plot?</li>
<li>How is frequency information found and added to the plot (mathematically)?</li>
</ol>
Once you understand the figure, continue reading to understand <em>what do these 20 complex numbers mean, and how should we interpret them?</em>
</p>
</div>
%% Cell type:code id:3acb465b tags:
``` python
abs_fft = np.abs(np.fft.fft(xt))
index_fft = np.arange(0,20,1)
plt.plot(index_fft, abs_fft, 'o')
freq = np.arange(0, 1, 0.05)
for x,y in zip(index_fft, abs_fft):
if x%5 == 0 or x==19:
label = f"f={freq[x]:.2f} Hz"
plt.annotate(label,
(x,y),
textcoords="offset points",
xytext=(0,10),
fontsize=10,
ha='center')
plt.ylim(0,5)
plt.xlim(-2,21)
plt.xlabel('fft-index')
plt.ylabel('$|X_k|$')
plt.stem(index_fft, abs_fft);
```
%% Output
%% Cell type:markdown id:09c72cdd tags:
The frequency resolution $\Delta f$ equals one-over-the-measurement-duration, hence $\Delta f = 1/T$. With that knowledge, we can reconstruct the frequencies expressed in Hertz.
The spectrum of the sampled signal is periodic in the sampling frequency $f_s$ which equals 1 Hz ($\Delta t$ = 1 s). Therefore, it is computed just for one period $[0,f_s)$. The last value, with index 19, represents the component with frequency $f$ = 0.95 Hz. A spectrum is commonly presented and interpreted as double-sided, so in the above graph we can interpret the spectral components with indices 10 to 19 corresponding to negative frequencies.
### Task 0.3: Identify negative frequencies
%% Cell type:code id:df1f483b tags:
``` python
abs_fft = np.abs(np.fft.fft(xt))
plt.stem(index_fft, abs_fft)
plt.plot(index_fft, abs_fft, 'o')
freq = np.concatenate((np.arange(0, 0.5, 0.05), np.arange(-0.5, 0, 0.05)))
for x,y in zip(index_fft, abs_fft):
if x%5 == 0 or x==19:
label = f"f={freq[x]:.2f} Hz"
plt.annotate(label,
(x,y),
textcoords="offset points",
xytext=(0,10),
fontsize=10,
ha='center')
plt.ylim(0,5)
plt.xlim(-2,21)
plt.xlabel('fft-index')
plt.ylabel('$|X_k|$');
```
%% Output
%% Cell type:markdown id:2d87a9a6 tags:
Now we can interpret the DFT of the $x(t)=\Pi(\frac{t}{4})$. The sampling interval is $\Delta t$ = 1 second, and the obsesrvation length is T=20 seconds, and we have N=20 samples.
- We can recognize a bit of a sinc function.
- The DFT is computed from 0 Hz, for positive frequencies up to $\frac{fs}{2} = 0.5$ Hz, after which the negative frequencies follow from -0.5 to -0.05 Hz.
- $X(f)=4 \textrm{sinc}(4f)$ has its first null at $f=0.25$Hz.
### Task 0.4: Create symmetric plot
For convenient visualization, we may want to explicitly shift the negative frequencies to the left-hand side to create a symmetric plot. We can use `numpy.fft.fftshift` to do that. In other words, the zero-frequency component appears in the center of the spectrum. Now, it nicely shows a symmetric sprectrum. It well resembles the sinc-function (taking into account that we plot the modulus/absolute value). The output of the DFT still consists of $N$ = 20 elements. To enable this, we set up a new frequency array (in Hz) on the interval $[-fs/2,fs/2)$.
%% Cell type:markdown id:e1eef6b6-1e90-43db-b059-6521b2a840d3 tags:
<div style="background-color:#AABAB2; color: black; vertical-align: middle; padding:15px; margin: 10px; border-radius: 10px; width: 95%">
<p>
<b>Task 0.4:</b>
Read the code cell below before executing it and identify how the plot is modified based on the (new) specification of frequency. Note that it is more than just the <code>freq</code> variable!
</p>
</div>
%% Cell type:code id:a8aac894 tags:
``` python
abs_fft_shift = np.abs(np.fft.fftshift(np.fft.fft(xt)))
freq = np.arange(-0.5, 0.5, 0.05)
plt.stem(freq, abs_fft_shift)
plt.plot(freq, abs_fft_shift, 'o')
plt.ylabel('|Xk|')
plt.xlabel('frequency [Hz]');
```
%% Output
%% Cell type:markdown id:569ca2d3 tags:
### Task 0.5: Showing spectrum only for positive frequencies
In practice, because of the symmetry, one typically plots only the right-hand side of the (double-sided) spectrum, hence only the part for positive frequencies $f \geq 0$. This is simply a matter of preference, and a way to save some space.
%% Cell type:markdown id:a26e361d-c5fa-47a2-b32c-416cc3d10659 tags:
<div style="background-color:#AABAB2; color: black; vertical-align: middle; padding:15px; margin: 10px; border-radius: 10px; width: 95%">
<p>
<b>Task 0.5:</b>
Can you identify what has changed (in the code and visually), compared to the previous plot?
</p>
</div>
%% Cell type:code id:b9608638 tags:
``` python
N=len(xt)
abs_fft = np.abs(np.fft.fft(xt))
freq = np.arange(0.0, 1.0, 0.05)
plt.plot(freq[:int(N/2)], abs_fft[:int(N/2)], 'o')
plt.stem(freq[:int(N/2)], abs_fft[:int(N/2)])
plt.ylabel('$|X_k|$')
plt.xlabel('frequency [Hz]');
```
%% Output
%% Cell type:markdown id:83c4ca67-9234-4741-a56f-294a12bcd98a tags:
<div style="background-color:#AABAB2; color: black; vertical-align: middle; padding:15px; margin: 10px; border-radius: 10px; width: 95%">
<p>
<b>Task 0.6:</b>
Confirm that you understand how we have arrived at the plot above, which illustrates the magnitude (amplitude) spectrum for frequencies $f \in [0,f_s/2)$, rather than $[0,f_s)$.
</p>
</div>
%% Cell type:markdown id:9f5558b4 tags:
## Task 1: Application of DFT using simple cosine
%% Cell type:markdown id:46dd5c0c tags:
It is always a good idea, in spectral analysis, to run a test with a very simple, basic signal. In this way you can test and verify your coding and interpretation of the results.
Our basic signal is just a plain cosine. We take the amplitude equal to one, and zero initial phase, so the signal reads $x(t) = \cos(2 \pi f_c t)$, with $f_c$ = 3 Hz in this exercise. With such a simple signal, we know in advance how the spectrum should look like. Namely just a spike at $f$ = 3 Hz, and also one at $f$ = -3 Hz, as we're, for mathematical convenience, working with double sided spectra. The spectrum should be zero at all other frequencies.
As a side note: the cosine is strictly a periodic function, not a-periodic (as above); still, the Fourier transform of the cosine is defined as two Dirac delta pulses or peaks (at 3 Hz and -3 Hz). You may want to check out the second worked example in Chapter 3 on the Fourier transform: Fourier transform in the limit.
%% Cell type:markdown id:87de3cc5 tags:
<div style="background-color:#AABAB2; color: black; vertical-align: middle; padding:15px; margin: 10px; border-radius: 10px; width: 95%">
<p>
<b>Task 1:</b>
Create a sampled (discrete time) cosine signal by sampling at $f_s$ = 10 Hz, for a duration of $T$ = 2 seconds (make sure you use exactly $N$ = 20 samples). Plot the sampled signal, compute its DFT and plot its magnitude spectrum $|X_k|$ with proper labeling of the axes (just like we did in the first part of this notebook). Include a plot of the spectrum of the sampled cosine signal using only the positive frequencies (up to $f_s/2$, as in the last plot of the previous task).
<em>Note: you are expected to produce three separate plots.</em>
</p>
</div>
%% Cell type:code id:28ad034e tags:
``` python
### SOLUTION
fc=3
fs=10
dt=1/fs
T=2
N=T*fs
t=np.arange(0,T,dt)
xt=np.cos(2*np.pi*fc*t)
plt.plot(t,xt)
plt.xlabel('time [s]')
plt.ylabel('xn');
```
%% Output
%% Cell type:code id:af87d120 tags:
``` python
### SOLUTION
abs_fft = np.abs(np.fft.fft(xt))
freq=np.arange(0,fs,1/T)
plt.stem(freq, abs_fft)
plt.plot(freq, abs_fft, 'o')
plt.ylabel('|Xk|')
plt.xlabel('frequency [Hz]');
```
%% Output
%% Cell type:code id:6fd9014c tags:
``` python
### SOLUTION
abs_fft = np.abs(np.fft.fft(xt))
freq=np.arange(0,fs,1/T)
plt.stem(freq[:int(N/2)], abs_fft[:int(N/2)])
plt.plot(freq[:int(N/2)], abs_fft[:int(N/2)], 'o')
plt.ylabel('|Xk|')
plt.xlabel('frequency [Hz]');
```
%% Output
%% Cell type:markdown id:6a4a306d-a372-460e-b312-e8fa138a9797 tags:
%% Cell type:markdown id:35961e4d tags:
**End of notebook.**
<h2 style="height: 60px">
</h2>
<h3 style="position: absolute; display: flex; flex-grow: 0; flex-shrink: 0; flex-direction: row-reverse; bottom: 60px; right: 50px; margin: 0; border: 0">
<style>
.markdown {width:100%; position: relative}
article { position: relative }
</style>
<a rel="license" href="http://creativecommons.org/licenses/by-nc-sa/4.0/">
<img alt="Creative Commons License" style="border-width:; width:88px; height:auto; padding-top:10px" src="https://i.creativecommons.org/l/by-nc-sa/4.0/88x31.png" />
<a rel="license" href="http://creativecommons.org/licenses/by/4.0/">
<img alt="Creative Commons License" style="border-width:; width:88px; height:auto; padding-top:10px" src="https://i.creativecommons.org/l/by/4.0/88x31.png" />
</a>
<a rel="TU Delft" href="https://www.tudelft.nl/en/ceg">
<img alt="TU Delft" style="border-width:0; width:100px; height:auto; padding-bottom:0px" src="https://gitlab.tudelft.nl/mude/public/-/raw/main/tu-logo/TU_P1_full-color.png"/>
<img alt="TU Delft" style="border-width:0; width:100px; height:auto; padding-bottom:0px" src="https://gitlab.tudelft.nl/mude/public/-/raw/main/tu-logo/TU_P1_full-color.png" />
</a>
<a rel="MUDE" href="http://mude.citg.tudelft.nl/">
<img alt="MUDE" style="border-width:0; width:100px; height:auto; padding-bottom:0px" src="https://gitlab.tudelft.nl/mude/public/-/raw/main/mude-logo/MUDE_Logo-small.png"/>
<img alt="MUDE" style="border-width:0; width:100px; height:auto; padding-bottom:0px" src="https://gitlab.tudelft.nl/mude/public/-/raw/main/mude-logo/MUDE_Logo-small.png" />
</a>
</h3>
<span style="font-size: 75%">
&copy; Copyright 2023 <a rel="MUDE Team" href="https://studiegids.tudelft.nl/a101_displayCourse.do?course_id=65595">MUDE Teaching Team</a> TU Delft. This work is licensed under a <a rel="license" href="http://creativecommons.org/licenses/by-nc-sa/4.0/">Creative Commons Attribution-NonCommercial-ShareAlike 4.0 International License</a>.
&copy; Copyright 2024 <a rel="MUDE" href="http://mude.citg.tudelft.nl/">MUDE</a> TU Delft. This work is licensed under a <a rel="license" href="http://creativecommons.org/licenses/by/4.0/">CC BY 4.0 License</a>.
......
src/teachers/Week_2_3/WS_2_3_DFT_you_try_meow.html 0 → 100644
+ 7980
− 0
View file @ 05aa9ad8
This diff is collapsed.
src/teachers/Week_2_3/WS_2_3_solution.html 0 → 100644
+ 8086
− 0
View file @ 05aa9ad8
This diff is collapsed.
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment