Frontmatter

If you are publishing this notebook on the web, you can set the parameters below to provide HTML metadata. This is useful for search engines and social media.

using Plots, WAV, FFTW, PlutoUI
2.8 s

Simple programs, scripts, debugging

md"""## Simple programs, scripts, debugging"""
176 μs

First we look at ...

  • Accessing MATLAB

    • Locally

    • Via Citrix

    • Via Terminal/Command Line

  • Getting used to the different windows in MATLAB

  • Simple commands and the command window

  • Scripts vs. functions vs. function files

md"""First we look at ...

- Accessing MATLAB
* Locally
* Via Citrix
* Via Terminal/Command Line
- Getting used to the different windows in MATLAB
- Simple commands and the command window
- Scripts vs. functions vs. function files
"""
3.0 ms

Then we demonstrate/practice some examples. Recall that the Fibonacci recurrences is given as:

$$\begin{align} F_0 &= 0\\ F_1 &= 1 \\ F_n &= F_{n-1} + F_{n-2} \end{align}$$

Can we write a function that returns the n-th term in this sequence?

md""" Then we demonstrate/practice some examples. Recall that the Fibonacci recurrences is given as:
```math
\begin{align}
F_0 &= 0\\
F_1 &= 1 \\
F_n &= F_{n-1} + F_{n-2}
\end{align}
```

Can we write a function that returns the n-th term in this sequence?
"""
252 μs

Simple Fibonacci

md"""#### Simple Fibonacci"""
173 μs
fib_simple (generic function with 1 method)
fib_simple(n) = n < 2 ? n : fib_simple(n-1) + fib_simple(n-2)
491 μs

Implement the above function in MATLAB. Then click on the arrow in the below cell to expand it. Compare.

md"""Implement the above function in MATLAB. Then click on the arrow in the below cell to expand it. Compare."""
180 μs
Simple Fibonacci - MATLAB

function fib = my_fib(n)
% Generates the n-th term in the Fibonacci sequence

if n < 2
    fib = n ;
else
    fib = my_fib(n-1) + my_fib(n-2);
end
end

Foldable("Simple Fibonacci - MATLAB", md"""
```octave
function fib = my_fib(n)
% Generates the n-th term in the Fibonacci sequence

if n < 2
fib = n ;
else
fib = my_fib(n-1) + my_fib(n-2);
end
end
```
""")
3.1 ms

There is a problem with the above function ... it is slow.

In the version above, the function calls itself twice to calculate the two previous terms, which are then added to find the current term. This process continues until the n-th term is calculated. However, this version has exponential time complexity and can be quite slow for larger values of n, especially since it recalculates the same terms multiple times.

Better Fibonacci

A better approach would be an iterative one.

md""" There is a problem with the above function ... it is slow.

In the version above, the function calls itself twice to calculate the two previous terms, which are then added to find the current term. This process continues until the n-th term is calculated. However, this version has exponential time complexity and can be quite slow for larger values of n, especially since it recalculates the same terms multiple times.

#### Better Fibonacci

A better approach would be an iterative one.

"""
321 μs
fibonacci (generic function with 1 method)
function fibonacci(n)
n_minus2, n_minus1 = 0,1
for i=1:n
n_minus2, n_minus1 = n_minus1, n_minus1 + n_minus2
end
return n_minus2
end
972 μs

Can you implement the above in MATLAB?

Note: MATLAB doesn't do swap assignments like above, so will need to make use of a temporary variable

md"""
Can you implement the above in MATLAB?


**Note:** MATLAB doesn't do swap assignments like above, so will need to make use of a temporary variable"""
261 μs
Better Fibonacci - MATLAB

function fib = my_fib(n)
% Generates the n-th term in the Fibonacci sequence

if n < 2
    fib = n;
else
    n_minus1 = 1;
    n_minus2 = 0;
    for i = 2:n
        fib = n_minus1 + n_minus2;
        n_minus2 = n_minus1;
        n_minus1 = fib;
    end
end
end

Foldable("Better Fibonacci - MATLAB", md"""
```matlab
function fib = my_fib(n)
% Generates the n-th term in the Fibonacci sequence

if n < 2
fib = n;
else
n_minus1 = 1;
n_minus2 = 0;
for i = 2:n
fib = n_minus1 + n_minus2;
n_minus2 = n_minus1;
n_minus1 = fib;
end
end
end
```
""")
183 μs

The above function uses a loop to iteratively calculate each term in the sequence. The first two terms, 0 and 1, are handled as special cases. For each subsequent term, the function adds the two previous terms to find the current term.

Advanced topic

If you still want to implement recursive Fibonacci, you can make use of memoization. this is the process of storing already calculated values in memory so we can reuse them instead of performing the computation again.

Memoized Fibonacci

md"""

The above function uses a loop to iteratively calculate each term in the sequence. The first two terms, 0 and 1, are handled as special cases. For each subsequent term, the function adds the two previous terms to find the current term.

_Advanced topic_

If you still want to implement recursive Fibonacci, you can make use of memoization. this is the process of storing already calculated values in memory so we can reuse them instead of performing the computation again.


#### Memoized Fibonacci"""
359 μs
Memoized Fibonacci - MATLAB

function fib = fibonacci(n, memo)
% Generates the n-th term in the Fibonacci sequence using memoization

if n <= 0
    fib = 0;
elseif n == 1
    fib = 1;
elseif n <= length(memo) && memo(n) ~= 0
    fib = memo(n);
else
    memo(n) = fibonacci(n-1, memo) + fibonacci(n-2, memo);
    fib = memo(n);
end

Foldable("Memoized Fibonacci - MATLAB",md"""
```matlab
function fib = fibonacci(n, memo)
% Generates the n-th term in the Fibonacci sequence using memoization

if n <= 0
fib = 0;
elseif n == 1
fib = 1;
elseif n <= length(memo) && memo(n) ~= 0
fib = memo(n);
else
memo(n) = fibonacci(n-1, memo) + fibonacci(n-2, memo);
fib = memo(n);
end
```
""")
160 μs

In this version, an additional argument memo is added to the function. The memo argument is a vector that stores the results of previous calculations. Before calculating a term, the function checks if the result is already stored in memo. If so, it returns the stored result. If not, it calculates the result and stores it in memo for future use.

To use this function, you would call fibonacci(n, zeros(1, n)), where n is the desired term and zeros(1, n) is the initial memo vector of zeros. The size of memo is set to n because the maximum term that needs to be calculated is n.

md"""
In this version, an additional argument memo is added to the function. The memo argument is a vector that stores the results of previous calculations. Before calculating a term, the function checks if the result is already stored in memo. If so, it returns the stored result. If not, it calculates the result and stores it in memo for future use.

To use this function, you would call `fibonacci(n, zeros(1, n))`, where `n` is the desired term and `zeros(1, n)` is the initial memo vector of zeros. The size of memo is set to `n` because the maximum term that needs to be calculated is `n`.
"""
299 μs

Plots & FFT

md"""## Plots & FFT"""
179 μs

In this section we transition to talking about the fft function. For this we will use audio files as our signals.

Load up an audio file and plot its waveform. Download the audios file from Ed. They are called sound_t.wav and sound_b.wav.

md"""
In this section we transition to talking about the `fft` function. For this we will use audio files as our signals.

Load up an audio file and plot its waveform. Download the audios file from Ed. They are called `sound_t.wav` and `sound_b.wav`.
"""
252 μs
y, fs, nbits, opt = wavread("files/sound_b.wav");
383 ms
plot(y, label=false, title="Plot of the waveform")
76.0 ms

Here

  • y - Waveform

  • fs - Sampling frequency

  • nbits and opts - Encoding options.

As you can see from examining the fs= 16000 Hz

md""" Here

* `y` - Waveform
* `fs` - Sampling frequency
* `nbits` and `opts` - Encoding options.

As you can see from examining the `fs=` $(Int32(fs)) Hz
"""
22.0 ms

Take the Fourier Transform and plot the magnitude and phase components

md""" Take the Fourier Transform and plot the magnitude and phase components"""
175 μs
begin
YF = fft(y)
freqs = fftfreq(length(y), fs)
mags = abs.(YF)
phase = angle.(YF);
end;
106 ms

In the above fftfreq is a function in Julia that will generate the frequency vector for us. MATLAB doesn't have an equivalent so you must generate the frequency vector yourself and in addition use the fftshift function!!

md""" In the above `fftfreq` is a function in Julia that will generate the frequency vector for us. MATLAB doesn't have an equivalent so you must generate the frequency vector yourself and in addition use the `fftshift` function!!
"""
198 μs
begin
p_mags = plot(freqs, mags, title = "Magnitude spectrum", label=false)
p_phase = plot(freqs, phase, title = "Phase spectrum", label=false)
plot(p_mags, p_phase, layout=(1, 2), size=(800,300))
end
191 ms

As we can see the phase plot isn't all that informative. But what is this negative frequency?

Check: Lecture note 8 for details.

We are interested in the portion 0 to 8 kHz because we only have reliable information at up to half the sampling frequency (here 16 kHz). So let us redo the plot.

md""" As we can see the phase plot isn't all that informative. But what is this negative frequency?

**Check:** [Lecture note 8](https://courses.grainger.illinois.edu/bioe205/sp2023/lectures/lec08/#computing_the_fourier_transform) for details.

We are interested in the portion 0 to 8 kHz because we only have reliable information at up to half the sampling frequency (here 16 kHz). So let us redo the plot.
"""
362 μs
plot(freqs, mags, title="Magnitude spectrum", xlim=(0, 8000), label=false)
17.7 ms
MATLAB Version

% READ the file sound_t.wave after downloading it from EdStem. 
% Below code intentionally leaves out generating the frequency vector. 
% Exercise!!

[y, Fs] = audioread("sound_b.wav");

FY = fft(y);
figure
subplot(2, 2, [1 2])
plot(y)
title("Waveform")
subplot(2, 2, 3)
plot(fftshift(abs(FY)))
title("Magnitude spectrum")
subplot(2, 2, 4)
plot(fftshift(angle(FY)))
title("Phase spectrum")

Foldable("MATLAB Version", md"""
```octave
% READ the file sound_t.wave after downloading it from EdStem.
% Below code intentionally leaves out generating the frequency vector.
% Exercise!!

[y, Fs] = audioread("sound_b.wav");

FY = fft(y);
figure
subplot(2, 2, [1 2])
plot(y)
title("Waveform")
subplot(2, 2, 3)
plot(fftshift(abs(FY)))
title("Magnitude spectrum")
subplot(2, 2, 4)
plot(fftshift(angle(FY)))
title("Phase spectrum")
```
""")
157 μs

Audio engineering?

In the next bit we are going to zero out the complex coefficients corresponding to a few frequencies and then reconstruct the wave form use the Inverse Fourier Transform.

Specifically we are zero-ing out the coefficients above 800 Hz.

md""" ### Audio engineering?

In the next bit we are going to zero out the complex coefficients corresponding to a few frequencies and then reconstruct the wave form use the Inverse Fourier Transform.

Specifically we are zero-ing out the coefficients above 800 Hz.
"""
2.4 ms
begin
pos_zero_indices = (freqs .> 800)
neg_zero_indices = (freqs .< -800)
Fcopy = deepcopy(YF)
end;
149 ms

Now we reconstruct the signal

md""" Now we reconstruct the signal"""
182 μs
begin
modded_y = real(ifft(Fcopy))
filename = "files/zeroed_5th.wav"
wavwrite(modded_y, filename, Fs=fs)
end
83.9 ms
DownloadButton(read(filename),basename(filename))
181 μs
Exercise: can you do the same in MATLAB?
md"""##### Exercise: can you do the same in MATLAB?"""
2.3 ms

Folding code

begin
struct Foldable{C}
title::String
content::C
end
function Base.show(io, mime::MIME"text/html", fld::Foldable)
write(io,"<details><summary>$(fld.title)</summary><p>")
show(io, mime, fld.content)
write(io,"</p></details>")
end
md"""`Folding code`"""
end
4.9 ms