// AMPLITUDE MODULATION

// Author: Federico Miyara
// Date: 2020-03-06

clear
close

// Sample rate
Fs = 44100
// Duration
T = 5
// Time vector
t = 0:1/Fs:T;

// FFT window length
N = 65536;

caso = 3
select caso
case 1
    // Arbitrary frequencies
    fp = 400
    fm = 40
case 2
    // Frequencies that can be exactly represented by 
    // FFT lines
    fp = floor(400/(Fs/N))*Fs/N
    fm = floor(50/(Fs/N))*Fs/N
case 3
    // Arbitrary frequencies
    fp = 400
    fm = 5
case 4
    fp = 400
    fm = 20   
end

// Modulation index
m = 0.5

// Carrier
xp = sin(2*%pi*fp*t);

// Amplitude envelope
A = (1 + m*sin(2*%pi*fm*t));

// Amplitude-modulated signal 
x = A.*xp;

// FFT spetrum of signal
X = fft(x(1:N));
X = abs([X(1)/N, X(2:N/2)*2/N]);
// Frequency vector  
f = (0:N/2-1)*Fs/N;

playsnd(0.1*x(1:2*Fs), Fs)

// Amplitude-modulated signal and its spetrum plots
scf(1);
clf(1);

subplot(2,1,1)
plot2d(t', [x',A'])
bounds = [0,-1.6;0.2,1.6];
gca().data_bounds = bounds;
gca().x_label.text = 't  [s]';
gca().y_label.text = "x(t)";

subplot(2,1,2)
plot2d(f, X)
gca().x_label.text = 'f  [Hz]';
gca().y_label.text = "X(f)";
bounds = [0,0;1000,1];
gca().data_bounds = bounds;

// Select values of modulation index
m = [0:0.1:0.9, 0.95, 0.99, 1]';
// Corresponding modulation depth values in dB
MD = 20*log10((1 + m)./(1 - m));

// Modulation depth as a function of modulation index
disp([m, MD])

// Select values of modulation depth 
MD = [0,3,6,10,20,30,40]';
// Corresponding modulation indices
m = (10^(MD/20) - 1)./(10^(MD/20) + 1);
disp([m, MD])

// Calculation for plot 
m = [0:0.01:0.99];
MD = 20*log10((1 + m)./(1 - m))

// modulation depth as a function of modulación index
figure(2)
gcf().background = 8;
plot(m, MD)
gca().grid = [1,1];