%--------------------------------------------------------------------------
%           Výuková podpora přemětu Modulované signály
%--------------------------------------------------------------------------
%                       Program QAM MODULACE 
%--------------------------------------------------------------------------
% 
% Výpočet a grafické vykreslení modulovaného vysokofrekvenčního signálu 
% pomocí digitální kvadraturní amplitudové QAM modulace (Quadrature Amplitude Modulation)
%
% Volitelný modulační (informační) signál a vysokofrekvenční nosná
% 
% Počet možných stavů zvolené modulace QAM m(t) je 16 
%
% VŠB - Technická univerzita Ostrava
% Fakulta elektrotechniky a informatiky
% Katedra měřicí a řídicí techniky
% 17.listopadu 15
% Ostrava - Poruba
% 708 33
%
% Vypracoval Zdeněk Macháček 2010
%--------------------------------------------------------------------------

%--------------------------------------------------------------------------
% Vyčištění všech proměnných paměti a uzavření všech grafických oken
clear all;
close all;

%--------------------------------------------------------------------------
% Specifikace časového průběhu-čas t[s] a nosné frekvence fc[Hz]RF signálu
f_c=1000; % nosná frekvence vysokofrekvenčního RF signálu
cas = 1; % čas průběhu simulace 
pocet_bodu=10; %počet bodů v čase během jedné periody fc
deltat=1/(f_c*pocet_bodu); % zanedbatelný časový okamžik
t=0:deltat:cas-deltat; % časový průběh
T_m = 0.05; % perioda informačního signálu m(t)

%--------------------------------------------------------------------------
% Definování a výpočet digitálního modulačního signálu m(t)
% v daném časovém okamžiku periody čtyřbitové číslo reprezentující m(t) 

t1=0:deltat:T_m-deltat;
m1(1:length(t1))=0;

t2=T_m:deltat:2*T_m-deltat;
m2(1:length(t2))=1;

t3=2*T_m:deltat:3*T_m-deltat;
m3(1:length(t3))=2;

t4=3*T_m:deltat:4*T_m-deltat;
m4(1:length(t4))=3;

t5=4*T_m:deltat:5*T_m-deltat;
m5(1:length(t5))=4;

t6=5*T_m:deltat:6*T_m-deltat;
m6(1:length(t6))=5;

t7=6*T_m:deltat:7*T_m-deltat;
m7(1:length(t7))=6;

t8=7*T_m:deltat:8*T_m-deltat;
m8(1:length(t8))=7;

t9=8*T_m:deltat:9*T_m-deltat;
m9(1:length(t9))=8;

t10=9*T_m:deltat:10*T_m-deltat;
m10(1:length(t10))=9;

t11=10*T_m:deltat:11*T_m-deltat;
m11(1:length(t11))=10;

t12=11*T_m:deltat:12*T_m-deltat;
m12(1:length(t12))=11;

t13=12*T_m:deltat:13*T_m-deltat;
m13(1:length(t13))=12;

t14=13*T_m:deltat:14*T_m-deltat;
m14(1:length(t14))=13;

t15=14*T_m:deltat:15*T_m-deltat;
m15(1:length(t15))=14;

t16=15*T_m:deltat:16*T_m-deltat;
m16(1:length(t16))=15;

t17=16*T_m:deltat:17*T_m-deltat;
m17(1:length(t17))=0;

t18=17*T_m:deltat:18*T_m-deltat;
m18(1:length(t18))=1;

t19=18*T_m:deltat:19*T_m-deltat;
m19(1:length(t19))=2;

t20=19*T_m:deltat:cas-deltat;
m20(1:length(t20))=3;


m_t = [m1 m2 m3 m4 m5 m6 m7 m8 m9 m10 m11 m12 m13 m14 m15 m16 m17 m18 m19 m20]; % modulační signál m(t)

%--------------------------------------------------------------------------
% QAM modulace - definování a výpočet pásmového signálu v(t)
% amplitudova modulacni slozka R_t = zmena stavu 
% fazova modulacni slozka Theta_t = změna stavu 
% soufázová modulační složka x_t = změna stavu
% kvadraturní modulační složka y_t = změna stavu


C = 0.2; % konstantní hodnota změny soufázových a kvadraturních modulačních složek

for i=1:length(m_t)
    if (m_t(i)==2) || (m_t(i)==6) || (m_t(i)==10) || (m_t(i)==14)
        x_t(i) = C;       % soufázová modulační složka 
    elseif (m_t(i)==3) || (m_t(i)==7) || (m_t(i)==11) || (m_t(i)==15)
        x_t(i) = 3*C;       % soufázová modulační složka 
    elseif (m_t(i)==1) || (m_t(i)==5) || (m_t(i)==9) || (m_t(i)==13)
        x_t(i) = -C;       % soufázová modulační složka 
    elseif (m_t(i)==0) || (m_t(i)==4) || (m_t(i)==8) || (m_t(i)==12)
        x_t(i) = -3*C;       % soufázová modulační složka 
    end;

    if (m_t(i)==4) || (m_t(i)==5) || (m_t(i)==6) || (m_t(i)==7)
        y_t(i) = C;       % kvadraturní modulační složka 
    elseif (m_t(i)==0) || (m_t(i)==1) || (m_t(i)==2) || (m_t(i)==3)
        y_t(i) = 3*C;       % kvadraturní modulační složka 
    elseif (m_t(i)==8) || (m_t(i)==9) || (m_t(i)==10) || (m_t(i)==11)
        y_t(i) = -C;       % kvadraturní modulační složka 
    elseif (m_t(i)==12) || (m_t(i)==13) || (m_t(i)==14) || (m_t(i)==15)
        y_t(i) = -3*C;       % kvadraturní modulační složka 
    end;
    
end;

% Pásmový signál v(t)
v_t=x_t.*cos(2*pi*f_c*t)-y_t.*sin(2*pi*f_c*t);
 
R_t=sqrt(x_t.^2+y_t.^2); % Amplitudova modulacni slozka 
Theta_t=atan2(y_t,x_t);  % Fázová modulacni slozka

%--------------------------------------------------------------------------
%Výpočet frekvenčního spektra pásmového signálu pomocí funkce algoritmu FFT
N=length(v_t); % počet hodnot pásmového signálu v(t) 
v_f = (fftshift(fft(v_t)))./N; % komplexní vektor frekvenčního spektra v(t) 
k=-N/2:N/2-1; % pomocný výpočet symetrického pole osy x
f = k.*pocet_bodu*f_c./N; % výpočet x-ové osy - frekvence f[Hz]

V_amp = abs(v_f);           %amplitudové frekvenční spektrum v(t)
V_vykon = V_amp.^2;         %výkonové frekvenční spektrum v(t)

% Je-li nulová hodnota amplitudy na dané frekvenci je také fáze = 0
for q=1:N
    if ((abs(real(v_f(q))) < 1e-1) && (abs(imag(v_f(q))) < 1e-1)) 
        v_f(q)=0;   
    end;
end;
V_faze = angle(v_f);        %fázové frekvenční spektrum v(t)

%--------------------------------------------------------------------------
% GRAFY - ZOBRAZENÍ SIGNÁLŮ V ČASOVÉ OBLASTI
%--------------------------------------------------------------------------
% Vykreslení grafického časového průběhu modulačního signálu
figure;
subplot(4,1,1);
plot(t,m_t);
title('Modulacní digitální signál m(t)');
ylabel('m(t)');
xlabel('t[s]');
 
%--------------------------------------------------------------------------
% Vykreslení grafického časového průběhu pásmového (modulovaného) signálu
subplot(4,1,2);
plot(t,v_t);
title('QAM modulace - Pásmový signál v(t)');
ylabel('v(t)');
xlabel('t[s]');

%--------------------------------------------------------------------------
% Vykreslení grafického časového průběhu amplitudové modulacní slozky R(t)
subplot(4,1,3);
plot(t,R_t)            
xlabel('t[s]')          % Popis x-ove osy
ylabel('R(t)')          % Popis y-ove osy
title('QAM modulace - Amplitudová modulacní slozka R(t)')    % Popis grafu

%--------------------------------------------------------------------------
% Vykreslení grafického časového průběhu fázové modulacní slozky Theta(t)
subplot(4,1,4);
plot(t,Theta_t)            
xlabel('t[s]')          % Popis x-ove osy
ylabel('\Theta(t)')          % Popis y-ove osy
title('QAM modulace - Fázová modulacní slozka \Theta(t)')    % Popis grafu

%--------------------------------------------------------------------------
% GRAFY - ZOBRAZENÍ SIGNÁLŮ VE FREKVENČNÍ OBLASTI
%--------------------------------------------------------------------------
%vykresleni amplitudového frekvenčního spektra 
figure;
subplot(3,1,1);
stem(f,V_amp)
xlabel('f[Hz]')
ylabel('^F^R|W_m|')
grid on;
title('QAM modulace - Amplitudové frekvencni spektrum pasmoveho (modulovaneho)signalu v(t)')
%--------------------------------------------------------------------------
%vykresleni fázového frekvenčního spektra
subplot(3,1,2);
stem(f,V_faze)
xlabel('f[Hz]')
ylabel('\Theta_m')
grid on
title('QAM modulace - Fazove frekvencni spektrum pasmoveho (modulovaneho)signalu v(t)')
%--------------------------------------------------------------------------
%vykresleni výkonového frekvenčního spektra
subplot(3,1,3);
stem(f,V_vykon)
xlabel('f[Hz]')
ylabel('^F^R|P_m|') 
grid on
title('QAM modulace - Vykonove frekvencni spektrum pasmoveho (modulovaneho)signalu v(t)')

%--------------------------------------------------------------------------
% GRAFY - ZOBRAZENÍ KONSTELACE MODULACE - IQ DIAGRAM
%--------------------------------------------------------------------------

%Výpočet a vykreslení konstelace modulace pásmového signálu
figure;
m_t_last = -1;
for i=1:length(m_t)
    if (m_t(i) ~= m_t_last)
        
        x = R_t.*cos(Theta_t);       	% soufázová modulační složka   
        y = R_t.*sin(Theta_t);       	% kvadraturní modulační složka   
        
        if ((y(i) < 1e-1) && ((y(i) > -1e-1))) % nepřesnost výsledku z důvodu zaokrouhlování
            y(i)=0;
        end;
        
        osax=[0,x(i)];
        osay=[0,y(i)];
        plot(osax,osay)
        hold on
        stem(x(i),y(i),'-wo','LineWidth',0.1,'MarkerEdgeColor','b','MarkerFaceColor','w','MarkerSize',8)
        m_t_last = m_t(i);
        if m_t(i)==0
            text(x(i),y(i)-0.05,'  data m(t)=0')
        elseif m_t(i)==1 
            text(x(i),y(i)-0.05,'  data m(t)=1')
        elseif m_t(i)==2 
            text(x(i),y(i)-0.05,'  data m(t)=2')
        elseif m_t(i)==3 
            text(x(i),y(i)-0.05,'  data m(t)=3')
        elseif m_t(i)==4 
            text(x(i),y(i)-0.05,'  data m(t)=4')
        elseif m_t(i)==5 
            text(x(i),y(i)-0.05,'  data m(t)=5')
        elseif m_t(i)==6 
            text(x(i),y(i)-0.05,'  data m(t)=6')
        elseif m_t(i)==7 
            text(x(i),y(i)-0.05,'  data m(t)=7')
        elseif m_t(i)==8 
            text(x(i),y(i)-0.05,'  data m(t)=8')
        elseif m_t(i)==9 
            text(x(i),y(i)-0.05,'  data m(t)=9')
        elseif m_t(i)==10 
            text(x(i),y(i)-0.05,'  data m(t)=10')
        elseif m_t(i)==11 
            text(x(i),y(i)-0.05,'  data m(t)=11')
        elseif m_t(i)==12 
            text(x(i),y(i)-0.05,'  data m(t)=12')
        elseif m_t(i)==13 
            text(x(i),y(i)-0.05,'  data m(t)=13')
        elseif m_t(i)==14 
            text(x(i),y(i)-0.05,'  data m(t)=14')
        elseif m_t(i)==15 
            text(x(i),y(i)-0.05,'  data m(t)=15')
        end; 
    end;
end;

xlabel('x soufázová modulacní slozka (in-phase axis)')
ylabel('y kvadraturní modulacni slozka (quadrature axis)') 
grid on
title('Konstelacní diagram (IQ diagram) - QAM modulovaná komplexní obálka g(t)')