close all;
clear all;

h1 = 0.8;
P = 100;

% Merici vzdalenost
r = 3;

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

% Vysky merici anteny - PRE-SCAN (10 m)
% h2_ver = [1 1.5 2 3.5];
% h2_hor = [1 1.5 2.5 4];

% Vysky merici anteny - PRE-SCAN (3 m)
h2_ver = [1 1.5 2 2.5];
h2_hor = [1 1.5 2 2.5];


Xplot =  30:10:1000;

%% Rozmitani vysky [h2] PRE-SCAN

X = max(size(h2_ver), size(h2_hor));
X = X(1,2);

for i = 1:1:X
    
r1_ver = sqrt((h2_ver(i)-h1)^2 + r^2);
r1_hor = sqrt((h2_hor(i)-h1)^2 + r^2);
r2_ver = sqrt((h1+h2_ver(i))^2 + r^2);
r2_hor = sqrt((h1+h2_hor(i))^2 + r^2);

del_r_ver = r2_ver - r1_ver;
del_r_hor = r2_hor - r1_hor;

q = 1;
    for freq = 30e6:10e6:1000e6
    
        lam = 3e8/freq;
        k = 2*pi/lam;
   
        alfaD_ver = (180/pi)*asin((h2_ver(i)-h1)/r1_ver);
        betaD_ver = (180/pi)*asin((h2_ver(i)+h1)/r2_ver);

        alfaD_hor = (180/pi)*asin((h2_hor(i)-h1)/r1_hor);   
        betaD_hor = (180/pi)*asin((h2_hor(i)+h1)/r2_hor);

        Up_ver = round(alfaD_ver);
        Uo_ver = round(360-betaD_ver);

        Up_hor = round(alfaD_hor);
        Uo_hor = round(360-betaD_hor);


        E1_ver(i,q) = (sqrt(30*P*D(Up_ver))/r1_ver);
        E1_hor(i,q) = (sqrt(30*P*D(Up_hor))/r1_hor);
        E2_ver(i,q) = (sqrt(30*P*D(Uo_ver))/r2_ver);
        E2_hor(i,q) = (sqrt(30*P*D(Uo_hor))/r2_hor);

        Ev_pre(i,q) = abs( E1_ver(i,q) + E2_ver(i,q)*exp(j*k*del_r_ver));
        Eh_pre(i,q) = abs( E1_hor(i,q) - E2_hor(i,q)*exp(j*k*del_r_hor));
 
     q = q+1;
    end
end

 % Nalezení maxima ze vsech vysek pro kazdou frekvenci.
 for i = 1:1:(q-1)
     
    Evmax_pre(i) = max(Ev_pre(:,i));   
    Ehmax_pre(i) = max(Eh_pre(:,i));

 end
 
Evmax_pre_dB = 20*LOG10(Evmax_pre*1e6);
Ehmax_pre_dB = 20*LOG10(Ehmax_pre*1e6);
      

%% 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_pre, 'k')
hold on;
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_pre, 'k')
hold on;
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_pre_dB, 'k')
hold on;
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_pre_dB, 'k')
hold on;
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, 'g','LineWidth',2)
hold on;
semilogx(Xplot, Evmax_005_dB - Evmax_025_dB, 'r','LineWidth',2)
hold on;
set(gca,'FontSize',25)
title('Vertikální polarizace');
legend( '100 cm', '50 cm', '25 cm', 'Location','northwest','Orientation','vertical')
xlabel('Frekvence [MHz]');
ylabel('Odchylka [dB]');
ylim([0 6])
xlim([30 1000]);
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, 'g','LineWidth',2)
hold on;
semilogx(Xplot, Ehmax_005_dB - Ehmax_025_dB, 'r','LineWidth',2)
hold on;
set(gca,'FontSize',25)
title('Horizontální polarizace');
legend( '100 cm', '50 cm', '25 cm', 'Location','northwest','Orientation','vertical')
xlabel('Frekvence [MHz]');
ylabel('Odchylka [dB]');
ylim([0 5])
xlim([30 1000]);
set(gca,'FontSize',20)



%% Measured


 %% Step measured

filename = 'matlab_data.xlsx';
sheet = 7;

freq = xlsread(filename, sheet, 'B5:B28');

%ver
VER_100    = xlsread(filename, sheet, 'W5:W28');
VER_50   = xlsread(filename, sheet, 'P5:P28');
VER_25  = xlsread(filename, sheet, 'I5:I28');

%hor
HOR_100    = xlsread(filename, sheet, 'W33:W56');
HOR_50   = xlsread(filename, sheet, 'P33:P56');
HOR_25  = xlsread(filename, sheet, 'I33:I56');

%% Make it fine

 O = 1;
for Q = 1:1:length(freq)
    for i = 1:1:5
        freqs(O) = freq(Q);
  
        VER_100s(O) = VER_100(Q);
        VER_50s(O) = VER_50(Q);
        VER_25s(O) = VER_25(Q);
        
        HOR_100s(O) = HOR_100(Q);
        HOR_50s(O) = HOR_50(Q);
        HOR_25s(O) = HOR_25(Q);
       
       
        
       O = O+1;
    end
end
   W = window(@hann, 10);
      
   freqn = filtfilt(W, 4.5, freqs);
      
   VER_100n = filtfilt(W, 4.5, VER_100s);
   VER_50n = filtfilt(W, 4.5, VER_50s);
   VER_25n = filtfilt(W, 4.5, VER_25s);
   
   HOR_100n = filtfilt(W, 4.5, HOR_100s);
   HOR_50n = filtfilt(W, 4.5, HOR_50s);
   HOR_25n = filtfilt(W, 4.5, HOR_25s);


     
%% Results 

% figure(1)
% set(gcf, 'color', [1 1 1])
% 
% subplot(2,1,1)
% semilogx(freq, VER_100, 'b','LineWidth',2)
% hold on;
% semilogx(freq, VER_50, 'g','LineWidth',2)
% hold on;
% semilogx(freq, VER_25, 'r','LineWidth',2)
% hold on;
% 
% set(gca,'FontSize',25)
% title('Vertikální polarizace');
% legend( '100 cm', '50 cm', '25 cm', 'Location','northwest','Orientation','vertical')
% xlabel('Frekvence [MHz]');
% ylabel('Odchylka [dB]');
% ylim([0 6])
% xlim([30 1000]);
% set(gca,'FontSize',20)
% 
% subplot(2,1,2)
% 
% semilogx(freq, HOR_100, 'b','LineWidth',2)
% hold on;
% semilogx(freq, HOR_50, 'g','LineWidth',2)
% hold on;
% semilogx(freq, HOR_25, 'r','LineWidth',2)
% hold on;
% set(gca,'FontSize',25)
% title('Horizontální polarizace');
% legend( '100 cm', '50 cm', '25 cm', 'Location','northwest','Orientation','vertical')
% xlabel('Frekvence [MHz]');
% ylabel('Odchylka [dB]');
% ylim([0 5])
% xlim([30 1000]);
% set(gca,'FontSize',20)

%%

% figure(2)
% set(gcf, 'color', [1 1 1])
% 
% subplot(2,1,1)
% plot(freqn, VER_100n, 'b','LineWidth',2)
% hold on;
% plot(freqn, VER_50n, 'g','LineWidth',2)
% hold on;
% plot(freqn, VER_25n, 'r','LineWidth',2)
% hold on;
% 
% set(gca,'FontSize',25)
% title('Vertikální polarizace');
% legend( '100 cm', '50 cm', '25 cm', 'Location','northwest','Orientation','vertical')
% xlabel('Frekvence [MHz]');
% ylabel('Odchylka [dB]');
% ylim([0 5])
% xlim([30 1000]);
% set(gca,'FontSize',20)
% 
% subplot(2,1,2)
% 
% plot(freqn, HOR_100n, 'b','LineWidth',2)
% hold on;
% plot(freqn, HOR_50n, 'g','LineWidth',2)
% hold on;
% plot(freqn, HOR_25n, 'r','LineWidth',2)
% hold on;
% set(gca,'FontSize',25)
% title('Horizontální polarizace');
% legend( '100 cm', '50 cm', '25 cm', 'Location','northwest','Orientation','vertical')
% xlabel('Frekvence [MHz]');
% ylabel('Odchylka [dB]');
% ylim([0 5])
% xlim([30 1000]);
% set(gca,'FontSize',20)

%%


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

subplot(2,1,1)
semilogx(freqn, VER_100n, 'b','LineWidth',2)
hold on;
semilogx(freqn, VER_50n, 'g','LineWidth',2)
hold on;
semilogx(freqn, VER_25n, 'r','LineWidth',2)
hold on;

set(gca,'FontSize',25)
title('Vertikální polarizace');
legend( '100 cm', '50 cm', '25 cm', 'Location','northwest','Orientation','vertical')
xlabel('Frekvence [MHz]');
ylabel('Odchylka [dB]');
ylim([0 6])
xlim([30 1000]);
set(gca,'FontSize',20)

subplot(2,1,2)

semilogx(freqn, HOR_100n, 'b','LineWidth',2)
hold on;
semilogx(freqn, HOR_50n, 'g','LineWidth',2)
hold on;
semilogx(freqn, HOR_25n, 'r','LineWidth',2)
hold on;
set(gca,'FontSize',25)
title('Horizontální polarizace');
legend( '100 cm', '50 cm', '25 cm', 'Location','northwest','Orientation','vertical')
xlabel('Frekvence [MHz]');
ylabel('Odchylka [dB]');
ylim([0 5])
xlim([30 1000]);
set(gca,'FontSize',20)

 %%
 
 figure(5)
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, '--g','LineWidth',2)
hold on;
semilogx(Xplot, Evmax_005_dB - Evmax_025_dB, '--r','LineWidth',2)
hold on;

semilogx(freqn, VER_100n, 'b','LineWidth',2)
hold on;
semilogx(freqn, VER_50n, 'g','LineWidth',2)
hold on;
semilogx(freqn, VER_25n, 'r','LineWidth',2)
hold on;

set(gca,'FontSize',25)
title('Vertikální polarizace');
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, 'g','LineWidth',2)
hold on;
semilogx(Xplot, Ehmax_005_dB - Ehmax_025_dB, 'r','LineWidth',2)
hold on;

semilogx(freqn, HOR_100n, '--b','LineWidth',2)
hold on;
semilogx(freqn, HOR_50n, '--g','LineWidth',2)
hold on;
semilogx(freqn, HOR_25n, '--r','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)
 
 
 
 
 