clc; clear all; format long g;

% Rozcestník pro výběr způsobu čtení dat při spuštění programu
option = questdlg('Jak chcete načíst data?', 'Volba způsobu čtení dat', 'Port', 'Soubor', 'Soubor');
switch option
    case 'Port'
        data = seriak(); % Volání funkce pro čtení dat ze sériového portu
    case 'Soubor'
        data = soubor(); % Volání funkce pro načtení dat ze souboru
    otherwise
        disp('Neplatná volba. Zvolte "port" nebo "soubor".');
        data = []; % Prázdná data v případě neplatné volby
end

%% výběr dat
% Inicializace dvou buňkových polí pro skupiny dat
gpgsvData = {};
glgsvData = {};
gnggaData = {};
gngsaData = {};

% Smyčka pro rozdělení dat do skupin podle počátečního řetězce
for i = 1:length(data)

    % Aktuální řádek dat
    line = data{i};
    
    % Rozdělení na skupiny podle počátečního řetězce
    if startsWith(line, '$GPGSV')
        gpgsvData = [gpgsvData; line];
    elseif startsWith(line, '$GLGSV')
        glgsvData = [glgsvData; line];
    elseif startsWith(line, '$GNGGA')
        gnggaData = [gnggaData; line];
    elseif startsWith(line, '$GNGSA')
        gngsaData = [gngsaData; line];
    end
end

%% konverze
% Inicializace pole pro ukládání extrahovaných čísel
GPGSV_vektor = cell(size(gpgsvData));
GLGSV_vektor = cell(size(glgsvData));

% % Iterace přes každý textový řetězec
 for i = 1:length(gpgsvData)

   gpgsvData{i} = strrep(gpgsvData{i}, ',,', ',0,');
   gpgsvData{i} = strrep(gpgsvData{i}, ',*', ',0*');

    % Nalezení číselných hodnot pomocí regulárního výrazu
    cisla = regexp(gpgsvData, '-?\d+(\.\d+)?', 'match');
    
    % Převod nalezených číselných řetězců na čísla
    GPGSV_vektor{i} = str2double(cisla{i});
    GPGSV_vektor{i} = GPGSV_vektor{i}(1:end-1);
 end

for i=1:length(glgsvData)

   glgsvData{i} = strrep(glgsvData{i}, ',,', ',0,');
   glgsvData{i} = strrep(glgsvData{i}, ',*', ',0*');

    % Nalezení číselných hodnot pomocí regulárního výrazu
    cisla = regexp(glgsvData, '-?\d+(\.\d+)?', 'match');
    
    % Převod nalezených číselných řetězců na čísla
    GLGSV_vektor{i} = str2double(cisla{i});
    GLGSV_vektor{i} = GLGSV_vektor{i}(1:end-1);

end

% Požadované pozice
pozice = [4, 8, 12, 16, 5, 9, 13, 17, 6, 10, 14, 18, 7, 11, 15, 19];

% Příprava výsledného pole
vysledek_GP = NaN(length(GPGSV_vektor), length(pozice));
vysledek_GL = NaN(length(GLGSV_vektor), length(pozice));

% Kopírování hodnot z polí do výsledného pole na základě pozic
for j = 1:length(GPGSV_vektor)
    for i = 1:length(pozice)
        if pozice(i) <= length(GPGSV_vektor{j})
            vysledek_GP(j, i) = GPGSV_vektor{j}(pozice(i));
        end
    end
end

GP_PRN_pre=vysledek_GP(:, 1:4);
GP_ELE_pre=vysledek_GP(:, 5:8);
GP_AZI_pre=vysledek_GP(:, 9:12);
GP_SNR_pre=vysledek_GP(:, 13:16);

GP_PRN_vektor=GP_PRN_pre(:);
GP_ELE_vektor=GP_ELE_pre(:);
GP_AZI_vektor=GP_AZI_pre(:);
GP_SNR_vektor=GP_SNR_pre(:);

M=[GP_PRN_vektor, GP_ELE_vektor, GP_AZI_vektor, GP_SNR_vektor];
M=sortrows(M, 1);


A=M(:,1); B=M(:,2); C=M(:,3); D=M(:,4);
% Vyřazení hodnot NaN z vektorů
valid_indices = ~isnan(A);
A_clean = A(valid_indices);
B_clean = B(valid_indices);
C_clean = C(valid_indices);
D_clean = D(valid_indices);

% Najděte unikátní platné hodnoty v A
unique_values = unique(A_clean);

% Určete maximální počet prvků podle nejčetnější hodnoty v A
max_count = max(histcounts(A_clean, unique_values));

% Vytvořte matice pro výsledky
GP_PRN = NaN(max_count, length(unique_values));
GP_ELE = NaN(max_count, length(unique_values));
GP_AZI = NaN(max_count, length(unique_values));
GP_SNR = NaN(max_count, length(unique_values));

% Vyplňte matice podle unikátních hodnot
for i = 1:length(unique_values)
    value = unique_values(i);
    indices = find(A_clean == value);
    
    % Vyplňte celý sloupec hodnotou 'value' v matici result_A
    GP_PRN(:, i) = value;
    
    % Vyplňte odpovídající hodnoty z B_clean do matice result_B
    for j = 1:length(indices)
        GP_ELE(j, i) = B_clean(indices(j));
        GP_AZI(j, i) = C_clean(indices(j));
        GP_SNR(j, i) = D_clean(indices(j));
    end
end

for j = 1:length(GLGSV_vektor)
    for i = 1:length(pozice)
        if pozice(i) <= length(GLGSV_vektor{j})
            vysledek_GL(j, i) = GLGSV_vektor{j}(pozice(i));
        end
    end
end

GL_PRN_pre=vysledek_GL(:, 1:4);
GL_ELE_pre=vysledek_GL(:, 5:8);
GL_AZI_pre=vysledek_GL(:, 9:12);
GL_SNR_pre=vysledek_GL(:, 13:16);

GL_PRN_vektor=GL_PRN_pre(:);
GL_ELE_vektor=GL_ELE_pre(:);
GL_AZI_vektor=GL_AZI_pre(:);
GL_SNR_vektor=GL_SNR_pre(:);

M=[GL_PRN_vektor, GL_ELE_vektor, GL_AZI_vektor, GL_SNR_vektor];
M=sortrows(M, 1);


A=M(:,1); B=M(:,2); C=M(:,3); D=M(:,4);
% Vyřazení hodnot NaN z vektorů
valid_indices = ~isnan(A);
A_clean = A(valid_indices);
B_clean = B(valid_indices);
C_clean = C(valid_indices);
D_clean = D(valid_indices);

% Najděte unikátní platné hodnoty v A
unique_values = unique(A_clean);

% Určete maximální počet prvků podle nejčetnější hodnoty v A
max_count = max(histcounts(A_clean, unique_values));

% Vytvořte matice pro výsledky
GL_PRN = NaN(max_count, length(unique_values));
GL_ELE = NaN(max_count, length(unique_values));
GL_AZI = NaN(max_count, length(unique_values));
GL_SNR = NaN(max_count, length(unique_values));

% Vyplňte matice podle unikátních hodnot
for i = 1:length(unique_values)
    value = unique_values(i);
    indices = find(A_clean == value);
    
    % Vyplňte celý sloupec hodnotou 'value' v matici result_A
    GL_PRN(:, i) = value;
    
    % Vyplňte odpovídající hodnoty z B_clean do matice result_B
    for j = 1:length(indices)
        GL_ELE(j, i) = B_clean(indices(j));
        GL_AZI(j, i) = C_clean(indices(j));
        GL_SNR(j, i) = D_clean(indices(j));
    end
end
%% Skyplot 
m=size(GP_PRN);
figure ('Name','Skyplot GPS satelitů') 
pax = polaraxes;
hold(pax, 'on');
for i = (1:m(2))
    polarplot(pax, deg2rad(GP_AZI(:,i)),GP_ELE(:,i), 'o');  
end

% Přizpůsobení vizuálního stylu grafu
%pax = gca; % Získání aktuální osy
pax.RAxis.Direction = 'reverse'; % Obrácení směru osy r
pax.ThetaZeroLocation = 'top'; % Nastavení nulového azimutu nahoru
pax.ThetaDir = 'clockwise'; % Nastavení směru azimutu
% ve směru hodinových ručiček
rlim([0,90]) % Nastavení limitu na radiální ose
pax.RTick = 0:10:90; % Nastavení kroku na radiální ose
pan on;
zoom on;
hold(pax, 'off');

% Inicializace prázdného pole pro legendu
    legenda = cell(1, m(2));

% Procházení každého prvku vektoru
for i = 1:m(2)
    % Vytvoření textu pro jednotlivé prvky legendy
    legenda{i} = ['GPS ', num2str(GP_PRN(1,i))];
end

legend(legenda, 'Location', 'eastoutside', 'Orientation', 'vertical'); 
% Přidání legendy mimo graf

m=size(GL_PRN);
figure ('Name','Skyplot satelitů GLONASS')
pax = polaraxes;
hold(pax, 'on');
for i = (1:m(2))
    polarplot(pax, deg2rad(GL_AZI(:,i)),GL_ELE(:,i), 'o');  
end

pax.RAxis.Direction = 'reverse'; % Obrácení směru osy r
pax.ThetaZeroLocation = 'top'; % Nastavení nulového azimutu nahoru
pax.ThetaDir = 'clockwise'; % Nastavení směru azimutu ve směru hodinových ručiček
rlim([0,90]) % Nastavení limitu na radiální ose
pax.RTick = 0:10:90; % Nastavení kroku na radiální ose
pan on;
zoom on;
hold(pax, 'off');

% Inicializace prázdného pole pro legendu
    legenda = cell(1, m(2));

% Procházení každého prvku vektoru
for i = 1:m(2)
    % Vytvoření textu pro jednotlivé prvky legendy
    legenda{i} = ['GLO ', num2str(GL_PRN(1,i))];
end
legend(legenda, 'Location', 'eastoutside', 'Orientation', 'vertical'); % Přidání legendy mimo graf

%% SNR grafy
% GPS
m=size(GP_PRN);
figure ('Name','SNR satelitů GPS')
    plot (1:(length(GP_SNR)), GP_SNR)
    grid on
    % Inicializace prázdného pole pro legendu
    legenda = cell(1, m(2));

% Procházení každého prvku vektoru
for i = 1:m(2)
    % Vytvoření textu pro jednotlivé prvky legendy
    legenda{i} = ['GPS ', num2str(GP_PRN(1,i))];
end

legend(legenda, 'Location', 'eastoutside', 'Orientation', 'vertical'); 
% Přidání legendy mimo graf

% GLONASS
m=size(GL_PRN);
figure ('Name','SNR satelitů GLONASS')
    plot (1:length(GL_SNR), GL_SNR)
    grid on
    % Inicializace prázdného pole pro legendu
    legenda = cell(1, m(2));

% Procházení každého prvku vektoru
for i = 1:m(2)
    % Vytvoření textu pro jednotlivé prvky legendy
    legenda{i} = ['GLO ', num2str(GL_PRN(1,i))];
end

legend(legenda, 'Location', 'eastoutside', 'Orientation', 'vertical'); % Přidání legendy mimo graf

%% konverze GGA
% Inicializace pole pro ukládání extrahovaných čísel
GNGGA_vektor = cell(size(gnggaData));

% % Iterace přes každý textový řetězec
 for i = 1:length(gnggaData)

   gnggaData{i} = strrep(gnggaData{i}, ',,', ',0,');
   gnggaData{i} = strrep(gnggaData{i}, ',,,', ',0,0,');
   gnggaData{i} = strrep(gnggaData{i}, ',*', ',0*');
 end
 
 for i = 1:length(gnggaData)

   gnggaData{i} = strrep(gnggaData{i}, ',,', ',');
    
    % Nalezení číselných hodnot pomocí regulárního výrazu
    cisla = regexp(gnggaData, '-?\d+(\.\d+)?', 'match');

    % Převod nalezených číselných řetězců na čísla
    GNGGA_vektor{i} = str2double(cisla{i});
    GNGGA_vektor{i} = GNGGA_vektor{i}(1:end-1);
 end

 for i=1:length(GNGGA_vektor)
     CAS(i)=GNGGA_vektor{i}(1);
     N(i)=GNGGA_vektor{i}(2)/100;
     E(i)=GNGGA_vektor{i}(3)/100;
     H(i)=GNGGA_vektor{i}(7);
 end

% Výpočet kvartilů a interkvartilového rozpětí (IQR)
Q1 = prctile(N, 25);
Q3 = prctile(N, 75);
IQR = Q3 - Q1;

% Definice dolního a horního limitu
lowerLimit = Q1 - 1.5 * IQR;
upperLimit = Q3 + 1.5 * IQR;

% Identifikace a filtrování odlehlých hodnot
N = N(N>= lowerLimit & N <= upperLimit);
N=fix(N)+((N-fix(N))*(10/6));

% Výpočet kvartilů a interkvartilového rozpětí (IQR)
Q1 = prctile(E, 25);
Q3 = prctile(E, 75);
IQR = Q3 - Q1;

% Definice dolního a horního limitu
lowerLimit = Q1 - 1.5 * IQR;
upperLimit = Q3 + 1.5 * IQR;

% Identifikace a filtrování odlehlých hodnot
E = E(E>= lowerLimit & E <= upperLimit);
E=fix(E)+((E-fix(E))*(10/6));

% Výpočet kvartilů a interkvartilového rozpětí (IQR)
Q1 = prctile(H, 25);
Q3 = prctile(H, 75);
IQR = Q3 - Q1;

% Definice dolního a horního limitu
lowerLimit = Q1 - 1.5 * IQR;
upperLimit = Q3 + 1.5 * IQR;

% Identifikace a filtrování odlehlých hodnot
H = H(H>= lowerLimit & H <= upperLimit);

NP=mean(N); stredni_chybaN = std(N);
EP=mean(E); stredni_chybaE = std(E);
VYP=mean(H); stredni_chybaH = std(H);

NP=sprintf('%.8f',NP);
EP=sprintf('%.8f',EP);
VYP=sprintf('%.4f',VYP);
disp(['N: ', NP, ' +- ', num2str(stredni_chybaN)]);
disp(['E: ', EP, ' +- ', num2str(stredni_chybaE)]);
disp(['H: ', VYP, ' +- ', num2str(stredni_chybaH)]);

% Výpočet četností pro jednotlivé intervaly s automatickým určením intervalů
[counts, edges] = histcounts(E, 'BinMethod', 'auto');

% Vytvoření sloupcového grafu
figure ('Name','četnost E')
bar(edges(1:end-1) + diff(edges)/2, counts,  'histc'); % Středové hodnoty intervalů
grid on


% Výpočet četností pro jednotlivé intervaly s automatickým určením intervalů
[counts, edges] = histcounts(N, 'BinMethod', 'auto');

% Vytvoření sloupcového grafu
figure ('Name','četnost N')
bar(edges(1:end-1) + diff(edges)/2, counts,  'histc'); % Středové hodnoty intervalů
grid on


% Výpočet četností pro jednotlivé intervaly s automatickým určením intervalů
[counts, edges] = histcounts(H, 'BinMethod', 'auto');

% Vytvoření sloupcového grafu
figure ('Name','četnost H')
bar(edges(1:end-1) + diff(edges)/2, counts, 'histc'); % Středové hodnoty intervalů
grid on


%% konverze GSA
% Inicializace pole pro ukládání extrahovaných čísel
GNGSA_vektor = cell(size(gngsaData));

for i = 1:length(gngsaData)
    
    % Nalezení číselných hodnot pomocí regulárního výrazu
    cisla = regexp(gngsaData, '-?\d+(\.\d+)?', 'match');

    % Převod nalezených číselných řetězců na čísla
    GNGSA_vektor{i} = str2double(cisla{i});
    GNGSA_vektor{i} = GNGSA_vektor{i}(1:end-1);
 end

 for i=1:length(GNGSA_vektor)
    PDOP(i)=GNGSA_vektor{i}(:, end-2);
    HDOP(i)=GNGSA_vektor{i}(:, end-1);
    VDOP(i)=GNGSA_vektor{i}(:, end);
 end

% Výpočet kvartilů a interkvartilového rozpětí (IQR)
Q1 = prctile(PDOP, 25);
Q3 = prctile(PDOP, 75);
IQR = Q3 - Q1;

% Definice dolního a horního limitu
lowerLimit = Q1 - 1.5 * IQR;
upperLimit = Q3 + 1.5 * IQR;

% Identifikace a filtrování odlehlých hodnot
PDOP = PDOP(PDOP>= lowerLimit & PDOP <= upperLimit);

% Výpočet kvartilů a interkvartilového rozpětí (IQR)
Q1 = prctile(HDOP, 25);
Q3 = prctile(HDOP, 75);
IQR = Q3 - Q1;

% Definice dolního a horního limitu
lowerLimit = Q1 - 1.5 * IQR;
upperLimit = Q3 + 1.5 * IQR;

% Identifikace a filtrování odlehlých hodnot
HDOP = HDOP(HDOP>= lowerLimit & HDOP <= upperLimit);

% Výpočet kvartilů a interkvartilového rozpětí (IQR)
Q1 = prctile(VDOP, 25);
Q3 = prctile(VDOP, 75);
IQR = Q3 - Q1;

% Definice dolního a horního limitu
lowerLimit = Q1 - 1.5 * IQR;
upperLimit = Q3 + 1.5 * IQR;

% Identifikace a filtrování odlehlých hodnot
VDOP = VDOP(VDOP>= lowerLimit & VDOP <= upperLimit);
 
PP=mean(PDOP);
HP=mean(HDOP);
VP=mean(VDOP);

disp(['PDOP: ', num2str(PP)]);
disp(['HDOP: ', num2str(HP)]);
disp(['VDOP: ', num2str(VP)]);


% graf PDOP
figure('Name','PDOP')
plot(PDOP, 'b')
grid on
hold on
yline(PP,'r--')
hold off

% graf HDOP
figure('Name','HDOP')
plot(HDOP, 'b')
grid on
hold on
yline(HP,'r--')
hold off

% graf VDOP
figure('Name','VDOP')
plot(VDOP, 'b')
grid on
hold on
yline(VP,'r--')
hold off

%% funkce soubor
function data = soubor()
    [filename, filepath] = uigetfile('*.*', 'Vyberte soubor'); % Dialogové okno pro výběr souboru
    if isequal(filename, 0) % Pokud uživatel nevybral soubor
        disp('Soubor nebyl vybrán!');
        data = []; % Prázdná data
    else
        disp(['Vybraný soubor: ', fullfile(filepath, filename)]);
        data = importdata(fullfile(filepath, filename)); % Načtení dat ze souboru
    end
end
%% funkce port
function data = seriak()

%% zadání vstupních hodnot
% Dialog pro výběr sériového portu
portPrompt = {'Zadejte název sériového portu:'};
portTitle = 'Výběr sériového portu';
portDims = [1 50];
portDefinput = {'COM3'}; % Výchozí hodnota pro sériový port
portAnswer = inputdlg(portPrompt, portTitle, portDims, portDefinput);

% Dialog pro výběr baudového přenosového pásma
baudPrompt = {'Zadejte baudrate:'};
baudTitle = 'Výběr baudrate';
baudDims = [1 50];
baudDefinput = {'9600'}; % Výchozí hodnota pro baudrate
baudAnswer = inputdlg(baudPrompt, baudTitle, baudDims, baudDefinput);

% Dialog pro zadání časového limitu měření
timePrompt = {'Zadejte časový limit měření (v sekundách):'};
timeTitle = 'Časový limit měření';
timeDims = [1 50];
timeDefinput = {'90'}; % Výchozí hodnota pro časový limit (90 sekund)
timeAnswer = inputdlg(timePrompt, timeTitle, timeDims, timeDefinput);

% Získání hodnot z dialogových oken
selectedPort = portAnswer{1};
selectedBaudrate = str2double(baudAnswer{1});
timeLimitInSeconds = str2double(timeAnswer{1});

% Pokračování s vybranými hodnotami
disp("Vybraný port: " + selectedPort);
disp("Vybrané baudrate: " + selectedBaudrate);
disp("Časový limit měření: " + timeLimitInSeconds + " sekund");

%% čtení dat z portu
try        
% Vytvoření objektu sériového portu s vybraným portem a baudrate
device = serialport(selectedPort, selectedBaudrate);

    startTime = tic; % Nastavení časového počátku
    startDateTime = datetime('now','Format','yyyy-MM-dd_HH-mm-ss'); 
    % Získání data a času začátku měření
    
    % Vytvoření názvu souboru s datem a časem začátku měření
    fileName = sprintf('data_%s.txt', char(startDateTime));
    
    % Otevření souboru pro zápis
    fileID = fopen(fileName, 'w');

    % Inicializace buňkového pole pro ukládání dat
    data = {};
    
    while true
        % Kontrola, zda uplynula určitá doba od začátku
        if toc(startTime) >= timeLimitInSeconds
            disp("Časový limit byl dosažen.");
            break; % Ukončení smyčky, pokud uplynul časový limit
        end
        
        % Přečtení řádku dat ze sériového portu
        line = readline(device);

        % Zápis řádku do souboru
        fprintf(fileID, '%s\n', line);

        % Přidání načteného řádku do buňkového pole
        data = [data; {line}];

        % Zpracování přečteného řádku dat
        disp(line); % Zobrazení přečteného řádku
    end

catch ME
    disp("Chyba při čtení ze sériového portu:");
    disp(ME.message);
end

clear device

end