// C.f. la série 9 sur les sons sinusoïdaux, les battements et les fréquences variables
// Le but est de générer et traiter des sons avec Scilab
// Références :
// https://help.scilab.org/docs/5.3.3/fr_FR/wavread.html
// https://help.scilab.org/doc/5.5.2/en_US/section_d11cd0f9e362390f953e7199c5bb4b3a.html
// https://atoms.scilab.org/

funcprot(0); // pour ne pas afficher le message : "redefining function".

function ex1()
//============
// Création d'un son pur, un sinus à 100 [Hz], qui dure 4 secondes.
// Avec une fréquence d'échantillonnage de Fs = 22'050 [Hz]
Fs = 22050; // [Hz]
tMax = 4; // [s]
nu = 100; // [Hz]
/* Amélioration ci-dessous
tt = 0; // [s]
jj = 0;
while (tt < tMax) do
    tt = jj / Fs;
    jj = jj+1
    at(jj) = tt;
    yy(jj) = 0.8*sin(2 * %pi * nu * tt);
end
*/
at = linspace(0, tMax, tMax*Fs);
yy = 0.8 * sin(2*%pi*nu*at); // Un sinus à fréquence nu

// Graphique du signal en fonction du temps.
// La période est bien de 0.01 [s]  = 1 / nu, comme attendu.
scf(1);
clf();
plot(at, yy);

playsnd(yy, Fs); // Joue le son

// Si on désire sauver le son dans un fichier.
//savewave("c:\Travail\Pur_100Hz.wav", yy, Fs, 16);
endfunction
//ex1();

function ex2()
//============
// Création d'un son pur, un sinus à 400 [Hz], qui dure 2 secondes.
// Avec une fréquence d'échantillonage de Fs = 22'050 [Hz]
Fs = 22050; // [Hz]
tMax = 2; // [s]
nu = 400; // [Hz]
at = linspace(0, tMax, tMax*Fs);
yy = 0.8 * sin(2*%pi*nu*at); // Un sinus à fréquence nu

// Graphique du signal en fonction du temps.
// La période est bien de 0.0025 [s]  = 1 / nu, comme attendu.
scf(1);
clf();
plot(at, yy);

playsnd(yy, Fs); // Joue le son
endfunction
//ex2();

function ex3()
//============
// Création d'un son pur, un sinus à 400 [Hz], qui dure 2 secondes.
// Avec une fréquence d'échantillonage de Fs = 22'050 [Hz]
// Avec un fondu en entrée de 0,5 [s] et un fondu en sortie de 0,5 [s].
Fs = 22050; // [Hz]
tMax = 2; // [s]
nu = 400; // [Hz]
at = linspace(0, tMax, tMax*Fs);
yy = 0.8 * sin(2*%pi*nu*at); // Un sinus à fréquence nu

/* Amélioration ci-dessous
// Fondu en ouverture de 0,5 seconde
nLen = 0.5 * Fs;
for jj=1:nLen do
    yy(jj) = yy(jj) * jj / nLen;
end

// Fondu en fermeture de 0,5 seconde
nLenTot = size(yy, 2);
nLen = 0.5 * Fs;
for jj = (nLenTot - nLen):nLenTot do
    yy(jj) = yy(jj) * (nLenTot - jj) / nLen;
end
*/

// Fondu en entrée de 0,5 seconde
yy(1:0.5*Fs) = yy(1:0.5*Fs) .* linspace(0, 1, 0.5*Fs);

// Fondu en sortie de 0,5 seconde
yy($-0.5*Fs:$) = yy($-0.5*Fs:$) .* linspace(1, 0, 0.5*Fs+1);

scf(1);
clf();
plot(at, yy, Fs);

playsnd(yy)
endfunction
ex3();

function ex4()
//============
// Création de deux sons purs, un sinus à 440 [Hz], 
// un de 440,5 [Hz] qui durent 6 secondes.
// En les additionnant, on crée un battement.
// L'intervalle de temps entre deux minimums de battement
// est égale à l'inverse de la différence des fréquences, soit 2 [s].

Fs = 22050; // [Hz]
tMax = 6; // [s]
nu1 = 400.0; // [Hz]
nu2 = 400.5; // [Hz]
at = linspace(0, tMax, tMax*Fs);
y1 = 0.4 * sin(2*%pi*nu1*at); // Un sinus à fréquence nu1
y2 = 0.4 * sin(2*%pi*nu2*at); // Un sinus à fréquence nu2

// Somme des deux signaux, pour avoir un battement
yy = y1 - y2;

scf(1);
clf();
plot(at, yy, 'm');

playsnd(yy, Fs);
endfunction
//ex4();

function ex5()
//============
// Création d'un son sinusoïdal de fréquence variable, qui dure 10 secondes.
// Avec une fréquence d'échantillonnage de Fs = 22'050 [Hz]
Fs = 22050; // [Hz]
tMax = 10; // [s]
nu0 = 100; // [Hz]
/* Amélioration ci-dessous
tt = 0; // [s]
jj = 0;
while (tt < tMax) do
    tt = jj / Fs;
    jj = jj+1
    at(jj) = tt;
    alpha = nu0 + 15*tt
    yy(jj) = 0.4*sin(2 * %pi * alpha * tt);
end
*/

at = linspace(0, tMax, tMax*Fs);
alpha = linspace(nu0, 15*tMax, tMax*Fs);
yy = 0.4 * sin(2*%pi .* alpha .* at); // Un sinus à fréquence variable

// Graphique du signal en fonction du temps.
// La période initiale est de 0.01 [s], donc fréqu. = 100 [Hz].
// La période finale est de 0.0025 [s], donc fréqu. = 400 [Hz].
scf(1);
clf();
plot(at, yy);

playsnd(yy, Fs); // Joue le son
endfunction
//ex5();
// On remarque que la fréquence commence à 100 [Hz] pour
// terminer à 400 [Hz] et non à 250 [Hz].
// alpha(t) ne correspond pas à la "fréquence instantanée".
// c'est la dérivée (d/dt) (alpha(t)*t) qui correspond à la "fréquence instantanée".

function ex6()
//============
// Création d'un son sinusoïdal de fréquence variable, qui dure 10 secondes.
// Avec une fréquence d'échantillonnage de Fs = 22'050 [Hz]
// La fréquence initiale vaut 100 [Hz]
// La fréquence finale vaut 200 [Hz]
Fs = 22050; // [Hz]
tMax = 10; // [s]
nu0 = 100; // [Hz]
delta_nu = 100; // [Hz]
/* Amélioration ci-dessous
tt = 0; // [s]
jj = 0;
while (tt < tMax) do
    tt = jj / Fs;
    jj = jj+1
    at(jj) = tt;
    phi_t = nu0*tt + delta_nu/2 * tt^2 / tMax;
    yy(jj) = 0.4*sin(2 * %pi * phi_t);
end
// phi'(t) = nu0 + delta_nu * t / tMax
*/

at = linspace(0, tMax, tMax*Fs);
phi_t = nu0 * at + delta_nu/(2 * tMax) .* at.^2;  // phi'(t) = nu0 + delta_nu * t / tMax
yy = 0.4 * sin(2*%pi .* phi_t); // Un sinus à fréquence variable

// Graphique du signal en fonction du temps.
// La période initiale est de 0.01 [s], donc fréqu. = 100 [Hz].
// La période finale est de 0.005 [s], donc fréqu. = 200 [Hz].
scf(1);
clf();
plot(at, yy);

playsnd(yy, Fs); // Joue le son
endfunction
//ex6();

function ex7()
//============
// Création d'un son sinusoïdal de fréquence variable, qui dure 17 secondes.
// Avec une fréquence d'échantillonnage de Fs = 2*44'100 [Hz]
// La fréquence initiale vaut 5'000 [Hz]
// La fréquence finale vaut 22'000 [Hz]
Fs = 4*22050; // [Hz]
tMax = 17; // [s]
nu0 = 5000; // [Hz]
delta_nu = 17000; // [Hz]
/* Amélioration ci-dessous
tt = 0; // [s]
jj = 0;
while (tt < tMax) do
    tt = jj / Fs;
    jj = jj+1
    at(jj) = tt;
    phi_t = nu0*tt + delta_nu/2 * tt^2 / tMax;
    yy(jj) = 0.4*sin(2 * %pi * phi_t)
end
// phi'(t) = nu0 + delta_nu * t / tMax
*/

at = linspace(0, tMax, tMax*Fs);
phi_t = nu0 * at + delta_nu/(2 * tMax) .* at.^2;  // phi'(t) = nu0 + delta_nu * t / tMax
yy = 0.4 * sin(2*%pi .* phi_t); // Un sinus à fréquence variable

// Graphique du signal en fonction du temps.
scf(1);
clf()
plot(at, yy);

//savewave("c:\Travail\Tacs2019_SinWaves_5_22_kHz_17s.wav", yy, Fs, 16);

disp("Départ");
playsnd(yy, Fs); // Joue le son
disp("Arrêt");
endfunction
//ex7();

function ex8()
//============
// Création d'un son sinusoïdal de fréquence variable, qui dure 10 secondes.
// Avec une fréquence d'échantillonnage de Fs = 22'050 [Hz]
// La fréquence initiale vaut 0 [Hz]
// La fréquence finale vaut  100 [Hz]
Fs = 22050; // [Hz]
tMax = 10; // [s]
nu0 = 0; // [Hz]
delta_nu = 10*tMax; // [Hz]

at = linspace(0, tMax, tMax*Fs);
phi_t = nu0 * at + delta_nu/(2 * tMax) .* at.^2;  // phi'(t) = nu0 + delta_nu * t / tMax
yy = 0.4 * sin(2*%pi .* phi_t); // Un sinus à fréquence variable

// Graphique du signal en fonction du temps.
scf(1);
clf()
plot(at, yy);

//savewave("c:\Travail\Tacs2019_SinWaves_0_100_Hz_10s.wav", yy, Fs, 16);

disp("Départ");
playsnd(yy, Fs); // Joue le son
disp("Arrêt");
endfunction
//ex8();

function ex8b()
//============
// Création d'un son en crénaux de fréquence variable, qui dure .. secondes.
// Avec une fréquence d'échantillonnage de Fs = 22'050 [Hz]
// La fréquence initiale vaut .. [Hz]
// La fréquence finale vaut    0 [Hz]
Fs = 22050; // [Hz]
tMax = 40; // [s]
nu0 = 40; // [Hz]
delta_nu = -1*tMax; // [Hz]
/* Amélioration ci-dessous
tt = 0; // [s]
jj = 0;
while (tt < tMax) do
    tt = jj / Fs;
    jj = jj+1
    at(jj) = tt;
    phi_t = nu0*tt + delta_nu/2 * tt^2 / tMax;
    yy(jj) = 0.4*sign(sin(2 * %pi * phi_t));
end
// phi'(t) = nu0 + delta_nu * t / tMax
*/

at = linspace(0, tMax, tMax*Fs);
phi_t = nu0 * at + delta_nu/(2 * tMax) .* at.^2;  // phi'(t) = nu0 + delta_nu * t / tMax
yy = 0.4 * sign(sin(2*%pi .* phi_t)); // Un sinus à fréquence variable

// Graphique du signal en fonction du temps.
scf(1);
clf()
plot(at, yy);

//savewave("c:\Travail\Tacs2019_CarreWaves_40_0_Hz_40s.wav", yy, Fs, 16);

disp("Départ");
playsnd(yy, Fs); // Joue le son
disp("Arrêt");
endfunction
//ex8b();

function ex9()
//============
// Création d'un son en crénaux de fréquence variable, qui dure 45 secondes.
// Avec une fréquence d'échantillonnage de Fs = 4*22'050 [Hz]
// Multiplication par 10 de la fréquence toutes les 10 secondes.
// La fréquence initiale vaut 1 [Hz]
// La fréquence finale vaut environ 10^4,5 = 31'600 [Hz]
Fs = 4*22050; // [Hz]
tMax = 45; // [s]
at = linspace(0, tMax, tMax*Fs);
phi_t = (10 / log(10)) * exp(0.1 * log(10) .* at); // phi'(t) = exp(0.1 * log(10) * tt) = 10^(t/10)
yy = 0.4 * sign(sin(2*%pi .* phi_t)); // Un sinus à fréquence variable

// Graphique du signal en fonction du temps.
scf(1);
clf()
plot(at, yy);

disp("Départ");
playsnd(yy, Fs); // Joue le son
disp("Arrêt");
endfunction
//ex9();

function ex9b()
//============
// Série 9, exercice 9.4*
// Création d'un son en crénaux de fréquence variable, qui dure 10 secondes.
// Avec une fréquence d'échantillonnage de Fs = 22'050 [Hz]
// La fréquence initiale vaut 150 [Hz]
// La fréquence finale vaut environ 1 [Hz]
Fs = 22050; // [Hz]
tMax = 22; // [s]
at = linspace(0, tMax, tMax*Fs);
phi_t = 150*(10 / log(10)) * exp(-0.1 * log(10) .* at); // phi'(t) = -150*exp(0.1 * log(10) * tt) = -150*10^(-t/10)
yy = 0.4 * sign(sin(2*%pi .* phi_t)); // Un sinus à fréquence variable

// Graphique du signal en fonction du temps.
scf(1);
clf();
plot(at, yy);

disp("Départ");
playsnd(yy, Fs); // Joue le son
disp("Arrêt");
endfunction
//ex9b();


function ex10()
//============
// Création d'un signal avec modulation d'amplitude
// Le signal est de haute fréquence, 
// son amplitude est modulée par le signal de beaucoup plus basse fréquence.
// Une fois transmis, on récupère l'amplitude du signal,
// Qui contient le son à jouer.

Fs = 16*22050; // [Hz]
tMax = 3; // [s]
nu0 = 22050; // [Hz] Fréquence de la porteuse
at = linspace(0, tMax, tMax*Fs);

nu_signal = 440; // [Hz] Fréquence du signal

// Signal = porteuse dont l'amplitude est modulée par le signal
yy = (0.5 + 0.4 * sin(2*%pi * nu_signal * at)) .* sin(2*%pi * nu0 * at)

// On extrait l'enveloppe supérieure du signal reçu,
// Pour récupérer le signal.
// Il faudrait le lisser pour l'éméliorer.
kmin = 8;
kmax = size(yy,2) - kmin;
for kk=kmin+1:kmax do
    ymax(kk) = max(yy(kk-kmin:kk+kmin))
end
ymax(1:kmin) = ymax(1);
ymax(kmax:size(at,2)) = ymax($);

// Enveloppe du pas du signal, pas utile ici.
for kk=kmin+1:kmax do
    ymin(kk) = min(yy(kk-kmin:kk+kmin))
end
ymin(1:kmin) = ymin(1);
ymin(kmax:size(at,2)) = ymin($);

// Graphique du signal en fonction du temps.
scf(1);
clf()
plot(at, yy);

// Graphique de l'enveloppe du signal en fonction du temps.
plot(at, ymax, 'r');
plot(at, ymin, 'm');

disp("Départ");
playsnd(ymax, Fs); // Joue le son
disp("Arrêt");
endfunction
//ex10();

function ex11()
//============
// Création d'une modulation de fréquence ou de phase.
// La difficulté est de garder une phase correcte pour démoduler le signal reçu,
// pour en extraire le son.
// Avec une fréquence d'échantillonnage de Fs

nu0 = 20; // [Hz] Fréquence de la porteuse
Fs = 16*nu0; // [Hz]
tMax = 3; // [s]
at = linspace(0, tMax, tMax*Fs);

nu1 = 0.5; // [Hz] Fréquence du signal
nu2 = 0.7; // [Hz] Fréquence du signal

// Son à transmettre
as = 0.6 + 0.2*sin(2*%pi * nu1 * at) + 0.24*sin(2*%pi * nu2 * at);
as = 0.2*sin(2*%pi * nu1 * at) + 0.24*sin(2*%pi * nu2 * at);

// Signal transmis = fréquence de la porteuse dont la fréquence est modulée par le son à transmettre
//phi = 2 * asin(as/2);
phi = as; // Si l'amplitude est faible, l'arcsinus(x) ~= x
yy = sin(2*%pi * nu0 * at + phi);

// Porteuse, pour décoder le signal
decal = 0.1;
y_porteuse = sin(2*%pi * nu0 * at + decal);

// Son, obtenu en combinant le signal transmis et la fréquence de la porteuse.
y_signal = yy - y_porteuse;

// Graphique du signal en fonction du temps.
scf(1);
clf()
plot(at, y_signal);
plot(at, as - decal, 'r');

// Graphique de l'onde transmise en fonction du temps
// L'amplitude est constante, la fréquence varie très légèrement.
// 
scf(2);
clf()
plot(at, yy);
//plot(at, y_porteuse, 'r');

disp("Départ");
//playsnd(ymax, Fs); // Joue le son
disp("Arrêt");
endfunction
//ex11();















