REF_LKPR_WGS84=[50.10083333 14.26000000 376];
gs_ref_airport_coord = REF_LKPR_WGS84;  % Referencni bod letiste

%% Definice TWYs
%            A-Zem.Sirka A-Zem.delka     A-Vyska         B-Zem.Sirka B-Zem.delka      B-Vyska           TWY  
gt_twy_wgs84 = [50.10177484 14.23359286 REF_LKPR_WGS84(1,3) 50.10046022 14.22921528 REF_LKPR_WGS84(1,3) 1   % TWY F1
               ,50.10306288 14.24290753 REF_LKPR_WGS84(1,3) 50.10221234 14.23610717 REF_LKPR_WGS84(1,3) 2   % TWY F2
               ,50.10408979 14.25112175 REF_LKPR_WGS84(1,3) 50.10377631 14.24861437 REF_LKPR_WGS84(1,3) 3   % TWY F3
               ,50.10468765 14.25590817 REF_LKPR_WGS84(1,3) 50.10453175 14.25466096 REF_LKPR_WGS84(1,3) 4   % TWY F4
               ,50.10601915 14.24359393 REF_LKPR_WGS84(1,3) 50.10516649 14.24266569 REF_LKPR_WGS84(1,3) 5   % TWY D1
               ,50.10496968 14.24964785 REF_LKPR_WGS84(1,3) 50.10492978 14.24951481 REF_LKPR_WGS84(1,3) 6   % TWY D2
               ,50.10597414 14.25459067 REF_LKPR_WGS84(1,3) 50.10570644 14.25244840 REF_LKPR_WGS84(1,3) 7   % TWY D3
               ,50.11259334 14.26852742 REF_LKPR_WGS84(1,3) 50.11217087 14.26711747 REF_LKPR_WGS84(1,3) 8   % TWY A1
               ,50.11389302 14.27286546 REF_LKPR_WGS84(1,3) 50.11296055 14.26975303 REF_LKPR_WGS84(1,3) 9]; % TWY A2
             
%% Sumarizace vsech dat do jedne tabulky gt_ouput
% Struktura tabulky gt_output vypada nasledovne
%ICAO Adress, Latitude, Longitude, Status HDG, HDG, Position Error,
%HDG Error, TWY, Type of A/C, Age of A/C,

%gt_output     ... vystupni tabulka
%lt_stat       ... tabulka pro vypocet obecnych statistik
%gv_mean_pe    ... prumer Position Error
%gv_median_pe  ... median Position Erorr
%gv_mean_hdg   ... prumer HDG Error
%gv_median_hdg ... median HDG Erorr

[m,n] = size(lt_surface_data);
[m_twy,n_twy] = size(gt_twy_wgs84);
gt_output = zeros(m,9);     % Definice velikosti vystupni tabulky
%lt_stat = zeros(m,1);       % Definice velikosti tabulky pro vypocet obecnych statistik

i_stat_pos = 1;                 % Pro pocitani hdonocenych zaznamu
i_stat_hdg = 1;
i_twy  = 1;                 % Pro pocitani Position Error a HDG Error
for i=1:m       
    gt_output(i,1) = lt_surface_data(i,1);      % ICAO ADress
    gt_output(i,2) = lt_surface_data(i,2);      % Latitude    
    gt_output(i,3) = lt_surface_data(i,3);      % Longitude
    gt_output(i,4) = lt_surface_data(i,4);      % Status HDG
    gt_output(i,5) = lt_surface_data(i,5);      % HDG
    gt_output(i,6) = 99999;
    gt_output(i,9) = lt_surface_data(i,7);      % Prijimac
    gt_output(i,10) = lt_surface_data(i,6);     % UNIX Time
    
    [ls_output_ecef(1,1),ls_output_ecef(1,2),ls_output_ecef(1,3)] = WGS2ECEF (gt_output(i,2),gt_output(i,3),gs_ref_airport_coord(1,3));       
    [ls_output_enu(1,1) ,ls_output_enu(1,2) ,ls_output_enu(1,3)]  = ECEF2ENU (gs_ref_airport_coord(1,1),gs_ref_airport_coord(1,2),gs_ref_airport_coord(1,3),ls_output_ecef(1,1),ls_output_ecef(1,2),ls_output_ecef(1,3));
    
    for i_twy=1:m_twy;        
        % Position & HDG Error
        if gt_output(i,6) == 99999;
            [lv_twy_x1,lv_twy_y1,lv_twy_z1] = WGS2ECEF (gt_twy_wgs84(i_twy,1),gt_twy_wgs84(i_twy,2),gt_twy_wgs84(i_twy,3));
            [lv_twy_x2,lv_twy_y2,lv_twy_z2] = WGS2ECEF (gt_twy_wgs84(i_twy,4),gt_twy_wgs84(i_twy,5),gt_twy_wgs84(i_twy,6));
        
            [lv_twy_e1,lv_twy_n1,lv_twy_u1] = ECEF2ENU (REF_LKPR_WGS84(1,1),REF_LKPR_WGS84(1,2),REF_LKPR_WGS84(1,3),lv_twy_x1,lv_twy_y1,lv_twy_z1);
            [lv_twy_e2,lv_twy_n2,lv_twy_u2] = ECEF2ENU (REF_LKPR_WGS84(1,1),REF_LKPR_WGS84(1,2),REF_LKPR_WGS84(1,3),lv_twy_x2,lv_twy_y2,lv_twy_z2);

            [gt_output(i,6),gt_output(i,7)] = Position_HDG_Error (ls_output_enu(1,1),ls_output_enu(1,2),gt_output(i,4),gt_output(i,5),lv_twy_e1,lv_twy_n1,lv_twy_e2,lv_twy_n2);
            if gt_output(i,6) == 99999;
                gt_output(i,8) = 0;
            else;
                gt_output(i,8) = gt_twy_wgs84(i_twy,7);
            end;
        end;
    end;  
    
    % Obecne statisticke udaje (Prumer, Median ...)
    if gt_output(i,6) ~= 99999;
        lt_stat_pos(i_stat_pos,1) = gt_output(i,6); % Position Error
        if gt_output(i,7) ~= 999;
            lt_stat_hdg(i_stat_hdg,1) = gt_output(i,7); % HDG Error
            i_stat_hdg = i_stat_hdg + 1;
        end;
        i_stat_pos = i_stat_pos + 1;
    end;   
end;

% Vypocet obecnych statistik
gv_mean_pe = mean(abs(lt_stat_pos(:,1)));
gv_median_pe = median(abs(lt_stat_pos(:,1)));
gv_mean_hdg = mean(abs(lt_stat_hdg(:,1)));
gv_median_hdg = median(abs(lt_stat_hdg(:,1)));
