% ZZS - 3.cvičení: Aditivní syntéza, R.Čmejla, 19.října 2015
% -----------------------------------------------------------
% PŘÍKLAD 1: Syntéza základních periodických průběhů a vliv fáze
%    i)  Syntéza z N harmonických podle Fourierova rozvoje
%    ii) Syntéza z N harmonických podle Fourierova rozvoje s náhodnou fází
%        Výsledné průběhy porovnejme v časové oblasti i poslechem.
%
%        Jedná se o ilustraci větší citlivosti lidského sluchu na obsah
%        harmonických složek a menší citlivosti na vliv fáze ve spektru.
% 

% a) Sled impulsů
clear, close all
f=440; fs=16000; doba=.5; t=0:1/fs:doba;

zvuk_1a(1,:)=cos(2*pi*f*t); 
zvuk_1b(1,:)=cos(2*pi*f*t+2*pi*rand);
for k=2:10
    zvuk_1a(k,:)=cos(k*2*pi*f*t);            
    zvuk_1b(k,:)=cos(k*2*pi*f*t+2*pi*rand);
    subplot(321), plot(t(1:200),sum(zvuk_1a(:,1:200)));
    subplot(322), plot(t(1:200),sum(zvuk_1b(:,1:200)));
    soundsc(sum(zvuk_1a),fs); pause(1.5*doba)
    soundsc(sum(zvuk_1b),fs); pause(1.5*doba)
end;
pause

% b) Obdélníkový průběh
zvuk_2a(1,:)=sin(2*pi*f*t);
zvuk_2b(1,:)=sin(2*pi*f*t+2*pi*rand);
for k=3:2:18
    zvuk_2a(k,:)=(1/k)*sin(k*2*pi*f*t);     
    zvuk_2b(k,:)=(1/k)*sin(k*2*pi*f*t+2*pi*rand);
    subplot(323), plot(t(1:200),sum(zvuk_2a(:,1:200)));
    subplot(324), plot(t(1:200),sum(zvuk_2b(:,1:200)));
    soundsc(sum(zvuk_2a),fs); pause(1.2*doba)
    soundsc(sum(zvuk_2b),fs); pause(1.2*doba)
end;
pause

% c) Pilový průběh
zvuk_3a(1,:)=sin(2*pi*f*t); 
zvuk_3b(1,:)=sin(2*pi*f*t+2*pi*rand);
for k=2:18
    zvuk_3a(k,:)=(1/k)*sin(k*2*pi*f*t);  
    zvuk_3b(k,:)=(1/k)*sin(k*2*pi*f*t+2*pi*rand);
    subplot(325), plot(t(1:200),sum(zvuk_3a(:,1:200)))
    subplot(326), plot(t(1:200),sum(zvuk_3b(:,1:200))), 
    soundsc(sum(zvuk_3a),fs),pause(1.2*doba)
    soundsc(sum(zvuk_3b),fs),,pause(1.2*doba)
end;

% -----------------------------------------------------------
% PŘÍKLAD 2: Modelování zvukové barvy 
% Syntéza hudebních nástrojů výchází ze spektrální analýzy, 
% kdy z frekvenčního spektra signálu odečteme příslušné amplitudy.

% Pro píšťaly, flétny a dřevěné hudební nástroje platí, že v jejich 
% frekvenčním spektru jsou především liché harmonické, 
% kdežto např. pro žestě, že mají zvýrazněné sudé harmonické.
%%
figure(2)
c1	= 261.63; d1	= 293.66; e1	= 329.63; 
f1	= 349.23; g1	= 392;	  a1	= 440;  
h1	= 493.88; c2	= 523.25;

noty	=	[c1 d1 e1 f1 g1 a1 h1  c2]; % skladba
doba 	= 	[.5 .5 .5 .5 .5 .5 .5 .75]; % trvani not v sekundach

signal=[];                             % inicializace promennych
zeste1=[];zeste2=[];zeste3=[];zeste4=[];
zeste5=[];zeste6=[];zeste7=[];zeste=[];
fletna1=[];fletna2=[];fletna3=[];fletna4=[];
fletna5=[];fletna6=[];fletna7=[];fletny=[];

fs=11025;                              % vzorkovaci frekvence

for i=1:length(noty)
    zeste1= [zeste1  0.2*sin(1*2*pi*noty(i)*(0:1/fs:doba(i)))]; % 1.osc.
    zeste2= [zeste2  0.6*sin(2*2*pi*noty(i)*(0:1/fs:doba(i)))]; % 2.osc.
    zeste3= [zeste3  0.4*sin(3*2*pi*noty(i)*(0:1/fs:doba(i)))]; % 3.osc.
    zeste4= [zeste4  0.6*sin(4*2*pi*noty(i)*(0:1/fs:doba(i)))]; % 4.osc.
    zeste5= [zeste5  0.4*sin(5*2*pi*noty(i)*(0:1/fs:doba(i)))]; % 5.osc.
    fletna1=[fletna1 0.8*sin(1*2*pi*noty(i)*(0:1/fs:doba(i)))]; % 1.osc.
    fletna2=[fletna2 0.0*sin(2*2*pi*noty(i)*(0:1/fs:doba(i)))]; % 2.osc.
    fletna3=[fletna3 0.4*sin(3*2*pi*noty(i)*(0:1/fs:doba(i)))]; % 3.osc.
    fletna4=[fletna4 0.0*sin(4*2*pi*noty(i)*(0:1/fs:doba(i)))]; % 4.osc.
    fletna5=[fletna5 0.2*sin(5*2*pi*noty(i)*(0:1/fs:doba(i)))]; % 5.osc.
end;

zeste=zeste1+zeste2+zeste3+zeste4+zeste5;
fletny=fletna1+fletna2+fletna3+fletna4+fletna5;
sound([zeste fletny],fs)

subplot(211), plot(fletny(1:499))
subplot(212), plot(zeste(1:499))

% -----------------------------------------------------------
% PŘÍKLAD 3: Harmonická analýza v MATLABu
>> analyza('banjo')
>> [X,Y]=ginput(10)

function analyza(soubor)
% funkce analyza(soubor) vykresli amplitudove
% spektrum *.wav souboru.
[signal,fs] = wavread(soubor);    

N = length(signal);               
c = fft(signal)/N;                    
A = 2*abs(c(2:floor(N/2)));         
f = (1:floor(N/2)-1)*fs/N;         
plot(f,A,'r')

% -----------------------------------------------------------
% PŘÍKLAD 4: Aditivní syntéza banja


% -----------------------------------------------------------
% PŘÍKLAD 5: Aditivní syntéza zvonku
 
fs =44100;          % vzorkovaci frekvence
T1 =0.06;           % doba jednotlivých úderů  
T2 =0.480;          % delka posledniho uderu v sérii
f01=180;            % zakladni frekvence 1.zvonku
f02=181;            % zakladni frekvence 2.zvonku
A=[.01 .1  1  1];   % amplitudy ctyr oscilatoru

K=[  5 10 20 40];   % nasobky zakl.frekvence        
                    % jednotlivych oscilatoru 
M=2;                % pocet serii zvoneni
N=4;                % celkovy pocet uderu = 2*N+1


% -----------------------------------------------------------
% PŘÍKLAD 6: Další příklady z aditivní syntézy:
% Syntetizujte následující zvuky 
% a zobrazte je v časové i frekvenční oblasti

    frekvenční složka = k * f0

    obálka            = A * t^n * exp(-t/tau)


    KLARINET                 f0 = 262 Hz
    ------------------------------------
    k       | A       | n       | tau
    ------------------------------------
    1.0000  | 1.0000  | 0.7500  | 0.7000
    2.0000  |      0  | 0.7500  | 0.3000
    3.0000  | 0.5000  | 0.7500  | 0.7000
    4.0000  |      0  | 0.7500  | 0.7000
    5.0000  | 0.2000  | 0.7500  | 0.7000


    DRNKNUTÍ STRUNY          f0 = 262 Hz
    ------------------------------------
    k       | A       | n       | tau
    ------------------------------------
    1.0000  | 1.0000  | 0.0100  | 0.7000
    2.0000  | 1.0000  | 0.0100  | 0.3000
    3.0000  | 0.5000  | 0.0100  | 0.1000
    4.0000  | 0.3000  | 0.0100  | 0.1000
    5.0000  | 0.1000  | 0.0100  | 0.1000


    BICÍ                     f0 = 262 Hz
    ------------------------------------
    k       | A       | n       | tau
    ------------------------------------
    1.0000  | 1.0000  | 0.0100  | 0.7000
    1.5800  | 1.0000  | 0.0100  | 0.3000
    3.0000  | 1.0000  | 0.0100  | 0.1000
    2.2400  | 0.3000  | 0.0100  | 0.4000
    2.5500  | 0.3000  | 0.0100  | 0.1000


    ZVON                     f0 = 262 Hz
    ------------------------------------
    k       | A       | n       | tau
    ------------------------------------
    1.0000  | 1.0000  | 0.0100  | 0.7000
    1.5800  | 1.0000  | 0.0100  | 0.3000
    3.0000  | 1.0000  | 0.0100  | 0.1000
    2.2400  | 0.3000  | 0.0100  | 0.2000
    2.5500  | 0.3000  | 0.0100  | 0.1000

% -----------------------------------------------------------
% PŘÍKLAD 7: Aditivní syntéza neharmonických signálů (parciály)
% Syntetizujte tympány a zvon s parametry:

    TYMPÁNY  f0 = 55 Hz         T = 2 s
    frekvenční složky = k * f0
    obálky            = A * exp(-2.8*t) * interp1(X,Y,t)
                        X = [0 .002 T*.99 T]
                        Y = [0 1  .9   0]                          
    ------------------
    k       | A      |
    ------------------
    1.0000  | 1.0000 |
    1.5958  | 0.4518 |
    2.1417  | 0.2740 |
    2.3000  | 0.2427 |
    2.6583  | 0.1897 |
    2.9250  | 0.1613 |
    3.1625  | 0.1412 |
    3.5083  | 0.1184 |
    3.6042  | 0.1131 |
    3.6542  | 0.1105 |
    4.0667  | 0.0921 |
    4.2375  | 0.0859 |
    4.6083  | 0.0745 |
    4.8417  | 0.0685 |
    4.9125  | 0.0668 |
    5.1417  | 0.0618 |
    5.4250  | 0.0564 |
    5.5500  | 0.0543 |
    5.9875  | 0.0477 |
    6.1667  | 0.0454 |
    6.2208  | 0.0447 |
    6.5417  | 0.0411 |
    6.7583  | 0.0388 |
    6.8625  | 0.0378 |
    7.3417  | 0.0337 |
    7.4833  | 0.0327 |
    7.9083  | 0.0297 |
    8.0875  | 0.0286 |
    8.6792  | 0.0254 |
    9.2583  | 0.0227 |

============================================
    ZVON  f0 = 116 Hz         T = 7 s
    frekvenční složky = k * f0
    obálka            = A * exp(-0.8*t) 
    ------------------
    k       | A      |
    ------------------
    1.0000  | 1.0000 |
    2.0000  | 2.7143 |
    2.3760  | 1.4286 |
    3.0600  | 0.4286 |
    4.0000  | 2.0000 |
    4.9400  | 0.2857 |
    5.2140  | 0.7143 |
    5.3000  | 1.0571 |
    5.9820  | 2.8571 |
    6.7340  | 0.5143 |
    8.2740  | 0.8571 |
    8.9740  | 0.2857 |
    9.6580  | 0.4286 |
   10.7700  | 0.8571 |
   11.7260  | 0.2857 |
   13.4180  | 0.2857 |
   16.1540  | 0.1429 |
   17.0940  | 0.0571 |
   18.0340  | 0.0286 |
   19.0600  | 0.1000 |
   22.0520  | 0.0143 |
   24.7860  | 0.0429 |

% -----------------------------------------------------------
% PŘÍKLAD 8: Další příklady z aditivní syntézy s obálkami typu ADSR

    FLÉTNY  f0 = 440 Hz         T = 1 s
            f1 = f0*2^(3/12) Hz T = 1 s
            f2 = f0*2^(7/12) Hz T = 1 s
    frekvenční složky = k * f0
    obálky        = A * interp1(X,Y,t)
                    X = [0 .2 .9 1]
                    Y = [0 1  .9 0]
    ----------------------------
    k       | A       | A      |
    ----------------------------
    1.0000  | 1.0000  | 1.0000 |
    2.0000  | 0.2500  | 0.0500 |
    3.0000  | 0.0625  | 0.0250 |
    4.0000  | 0.0156  | 0.0005 |
    5.0000  | 0.0625  | 0.0010 |
    6.0000  | 0.0156  | 0.0005 |

    DECHOVÉ NÁSTROJE   f0 = 440 Hz         T = 2 s
    frekvenční složky = k * f0
    obálky            = A * interp1(X,Y,t)
                    X = [0 .1 T*.9 T]
                    Y = [0 1  .9   0]
    ----------------------------
    k       | A"plný" | A"dutý"|
    ----------------------------
    1.0000  | 1.0000  | 1.0000 |
    2.0000  | 0.2500  | 0.0000 |
    3.0000  | 0.1111  | 0.1111 |
    4.0000  | 0.0625  | 0.0000 |
    5.0000  | 0.0400  | 0.0400 |
    6.0000  | 0.0278  | 0.0000 |
    7.0000  | 0.0204  | 0.0204 |
    8.0000  | 0.0156  | 0.0000 |
    9.0000  | 0.0123  | 0.0123 |
   10.0000  | 0.0100  | 0.0000 |
   11.0000  | 0.0083  | 0.0083 |
   12.0000  | 0.0069  | 0.0000 |
   13.0000  | 0.0059  | 0.0059 |
   14.0000  | 0.0051  | 0.0000 |
   15.0000  | 0.0044  | 0.0044 |

% PŘÍKLAD 9: Příklad z aditivní syntézy s obálkami typu ADSR

    TRUBKA  f0 = 440 Hz         T = 3 s
    frekvenční složky = k * f0
    obálky            = A * interp1(X,Y,t)
                            obálka Y  = [0 1 .8 0]
    ------------------------------------------------
    k       | A       |     obálka X               |
    ------------------------------------------------
    1.0000  | 1.0000  | 0    0.0071    0.9875    1 |
    2.0000  | 1.8222  | 0    0.0286    0.9750    1 |
    3.0000  | 1.4444  | 0    0.0643    0.9625    1 |
    4.0000  | 0.5556  | 0    0.1143    0.9500    1 |
    5.0000  | 0.2222  | 0    0.1786    0.9375    1 |
    6.0000  | 0.0889  | 0    0.2571    0.9250    1 |
    7.0000  | 0.0667  | 0    0.3500    0.9125    1 |
    8.0000  | 0.0444  | 0    0.4571    0.9000    1 |
    9.0000  | 0.0222  | 0    0.5786    0.8875    1 |