close all;
clear all;

h1 = 0.8;
P = 100;

% Merici vzdalenost
r = 3;

% Izotropni zaric (Smerovky)
D = ones(1,360);



Xplot =  30:10:1000;



%% Rozmitani vysky [h2] s skokem 0.05 m [reference]

i = 1;
for h2 = 1:0.05:4
    
    r1 = sqrt((h2-h1)^2 + r^2);
    r2 = sqrt((h1+h2)^2 + r^2);
    del_r = r2 - r1;

    q = 1;
    for freq = 30e6:10e6:1000e6
    
        lam = 3e8/freq;
        k = 2*pi/lam;

        alfaD = (180/pi)*asin((h2-h1)/r1);
        betaD = (180/pi)*asin((h2+h1)/r2);

        Up = round(alfaD);
        Uo = round(360-betaD);

        E1(i,q) = (sqrt(30*P*D(Up))/r1);
        E2(i,q) = (sqrt(30*P*D(Uo))/r2);

        Ev_005(i,q) = abs( E1(i,q) + E2(i,q)*exp(j*k*del_r));
        Eh_005(i,q) = abs( E1(i,q) - E2(i,q)*exp(j*k*del_r));
  
 

    q = q+1;
    end

i = i+1;
end

 % Nalezení maxima ze vsech vysek pro kazdou frekvenci.
 for i = 1:1:(q-1)
     
    Evmax_005(i) = max(Ev_005(:,i));   
    Ehmax_005(i) = max(Eh_005(:,i));

 end
 
Evmax_005_dB = 20*LOG10(Evmax_005*1e6);
Ehmax_005_dB = 20*LOG10(Ehmax_005*1e6);


%% Rozmitani vysky [h2] s skokem 0.25 m

i = 1;
for h2 = 1:0.25:4
    
    r1 = sqrt((h2-h1)^2 + r^2);
    r2 = sqrt((h1+h2)^2 + r^2);
    del_r = r2 - r1;

    q = 1;
    for freq = 30e6:10e6:1000e6
    
        lam = 3e8/freq;
        k = 2*pi/lam;

        alfaD = (180/pi)*asin((h2-h1)/r1);
        betaD = (180/pi)*asin((h2+h1)/r2);

        Up = round(alfaD);
        Uo = round(360-betaD);

        E1(i,q) = (sqrt(30*P*D(Up))/r1);
        E2(i,q) = (sqrt(30*P*D(Uo))/r2);

        Ev_025(i,q) = abs( E1(i,q) + E2(i,q)*exp(j*k*del_r));
        Eh_025(i,q) = abs( E1(i,q) - E2(i,q)*exp(j*k*del_r));
  
 

    q = q+1;
    end

i = i+1;
end

 % Nalezení maxima ze vsech vysek pro kazdou frekvenci.
 for i = 1:1:(q-1)
     
    Evmax_025(i) = max(Ev_025(:,i));   
    Ehmax_025(i) = max(Eh_025(:,i));

 end
 
Evmax_025_dB = 20*LOG10(Evmax_025*1e6);
Ehmax_025_dB = 20*LOG10(Ehmax_025*1e6);


%% Rozmitani vysky [h2] s skokem 0.5 m

i = 1;
for h2 = 1:0.50:4
    
    r1 = sqrt((h2-h1)^2 + r^2);
    r2 = sqrt((h1+h2)^2 + r^2);
    del_r = r2 - r1;

    q = 1;
    for freq = 30e6:10e6:1000e6
    
        lam = 3e8/freq;
        k = 2*pi/lam;

        alfaD = (180/pi)*asin((h2-h1)/r1);
        betaD = (180/pi)*asin((h2+h1)/r2);

        Up = round(alfaD);
        Uo = round(360-betaD);

        E1(i,q) = (sqrt(30*P*D(Up))/r1);
        E2(i,q) = (sqrt(30*P*D(Uo))/r2);

        Ev_05(i,q) = abs( E1(i,q) + E2(i,q)*exp(j*k*del_r));
        Eh_05(i,q) = abs( E1(i,q) - E2(i,q)*exp(j*k*del_r));
  
 

    q = q+1;
    end

i = i+1;
end

 % Nalezení maxima ze vsech vysek pro kazdou frekvenci.
 for i = 1:1:(q-1)
     
    Evmax_05(i) = max(Ev_05(:,i));   
    Ehmax_05(i) = max(Eh_05(:,i));

 end
 
Evmax_05_dB = 20*LOG10(Evmax_05*1e6);
Ehmax_05_dB = 20*LOG10(Ehmax_05*1e6); 
      

%% Rozmitani vysky [h2] s skokem 1 m

i = 1;
for h2 = 1:1:4
    
    r1 = sqrt((h2-h1)^2 + r^2);
    r2 = sqrt((h1+h2)^2 + r^2);
    del_r = r2 - r1;

    q = 1;
    for freq = 30e6:10e6:1000e6
    
        lam = 3e8/freq;
        k = 2*pi/lam;

        alfaD = (180/pi)*asin((h2-h1)/r1);
        betaD = (180/pi)*asin((h2+h1)/r2);

        Up = round(alfaD);
        Uo = round(360-betaD);

        E1(i,q) = (sqrt(30*P*D(Up))/r1);
        E2(i,q) = (sqrt(30*P*D(Uo))/r2);

        Ev_10(i,q) = abs( E1(i,q) + E2(i,q)*exp(j*k*del_r));
        Eh_10(i,q) = abs( E1(i,q) - E2(i,q)*exp(j*k*del_r));
  
 

    q = q+1;
    end

i = i+1;
end

 % Nalezení maxima ze vsech vysek pro kazdou frekvenci.
 for i = 1:1:(q-1)
     
    Evmax_10(i) = max(Ev_10(:,i));   
    Ehmax_10(i) = max(Eh_10(:,i));

 end
 
Evmax_10_dB = 20*LOG10(Evmax_10*1e6);
Ehmax_10_dB = 20*LOG10(Ehmax_10*1e6); 
      
%% Vykresleni prubehu


figure(1)     
subplot(2,1,1)
semilogx(Xplot, Evmax_10, '--k')
hold on;
semilogx(Xplot, Evmax_005, '--g')
hold on;
semilogx(Xplot, Evmax_025, '--r')
hold on;
semilogx(Xplot, Evmax_05, '--b')
hold on;
title('Vertikální polarizace');
legend('PRE-SCAN', 'krok 5 cm', 'krok 25 cm', 'krok 50 cm')
xlabel('Frekvence [MHz]');
ylabel('Intenzita [V/m]');

subplot(2,1,2)
semilogx(Xplot, Ehmax_10, 'k')
hold on;
semilogx(Xplot, Ehmax_005, '--g')
hold on;
semilogx(Xplot, Ehmax_025, '--r')
hold on;
semilogx(Xplot, Ehmax_05, '--b')
hold on;
title('Horizontální polarizace');
legend('PRE-SCAN', 'krok 5 cm', 'krok 25 cm', 'krok 50 cm')
xlabel('Frekvence [MHz]');
ylabel('Intenzita [V/m]');


figure(2)     
subplot(2,1,1)
semilogx(Xplot, Evmax_005_dB, '--g')
hold on;
semilogx(Xplot, Evmax_025_dB, '--r')
hold on;
semilogx(Xplot, Evmax_05_dB, '--b')
hold on;
title('Vertikální polarizace');
legend('PRE-SCAN', 'krok 5 cm', 'krok 25 cm', 'krok 50 cm')
xlabel('Frekvence [MHz]');
ylabel('Intenzita [dBuV/m]');

subplot(2,1,2)
semilogx(Xplot, Ehmax_005_dB, '--g')
hold on;
semilogx(Xplot, Ehmax_025_dB, '--r')
hold on;
semilogx(Xplot, Ehmax_05_dB, '--b')
hold on;
title('Horizontální polarizace');
legend('PRE-SCAN', 'krok 5 cm', 'krok 25 cm', 'krok 50 cm')
xlabel('Frekvence [MHz]');
ylabel('Intenzita [dBuV/m]');

%% Srovnani

figure(3)
set(gcf, 'color', [1 1 1])

subplot(2,1,1)
semilogx(Xplot, Evmax_005_dB - Evmax_10_dB, 'b','LineWidth',2)
hold on;
semilogx(Xplot, Evmax_005_dB - Evmax_05_dB, 'r','LineWidth',2)
hold on;
semilogx(Xplot, Evmax_005_dB - Evmax_025_dB, 'g','LineWidth',2)
hold on;
set(gca,'FontSize',25)
title('Vertikální polarizace -S');
legend( 'krok 100 cm', 'krok 50 cm', 'krok 25 cm')
xlabel('Frekvence [MHz]');
ylabel('Odchylka [dB]');
set(gca,'FontSize',20)

subplot(2,1,2)
semilogx(Xplot, Ehmax_005_dB - Ehmax_10_dB, 'b','LineWidth',2)
hold on;
semilogx(Xplot, Ehmax_005_dB - Ehmax_05_dB, 'r','LineWidth',2)
hold on;
semilogx(Xplot, Ehmax_005_dB - Ehmax_025_dB, 'g','LineWidth',2)
hold on;
set(gca,'FontSize',25)
title('Horizontální polarizace -S');
legend( 'krok 100 cm', 'krok 50 cm', 'krok 25 cm')
xlabel('Frekvence [MHz]');
ylabel('Odchylka [dB]');
set(gca,'FontSize',20)




