clc; clear all; format long;


%% ZADANI [°] [m] =========================================================
fA = dms2rad([  49  9 57.177 ]); % zadani 47
lA = dms2rad([  12 27 43.765 ]);
aA = dms2rad([ 272  4  1.093 ]);
s = 694141.67;

% vysledky z programu Elipsoid (vypocty na elipsoidu) ---------------------
fB = dms2rad([  48 59 56.06  ]);
lB = dms2rad([   2 57 21.61  ]);
% prusecik s 0 polednikem
fP = dms2rad([  48 47 08.9994 ]);
lP = 0;
% prusecik s rovnikem
fR = 0;
lR = dms2rad([ -80 04 19.1396 ]);

% NASTAVENI ZOBRAZENI =====================================================
stred = [24 -40 0]; % stred zobrazovaneho uzemi [fi lambda orientace]
% barva, styl
krivka = 'g';
bodA = 'kx';
bodB = 'rx';
bodP = 'rx';
bodR = 'rx';
polednik = 'k--'; % a rovnik
% rozmer
krivka_w = 2;
polednik_w = 2; % a rovnik
sit_w = 1;
bod = [10 2]; % velikost, sire cary ('MarkerSize', 'LineWidth')

%% RESENI =================================================================
% vsechny hodnoty behem reseni jsou v radianech a metrech
[a b e2] = bessel; % parametry Besselova elipsoidu:

%% PPØÍKLAD 2.1 ===========================================================
R = R(a, e2, fA); % stredni polomer v bode A
[fB, lB] = reckon(fA, lA, s/R, aA, 'rad'); % souradnice bodu B
aBA = azimuth(fB, lB, fA, lA, 'rad'); % azimut z B na A
    fprintf('PR 1. reseni na kouli:\n  R = % 10.4f m\n', R);
    fprintf(' fB ='); rad2dms(fB);
    fprintf(' lB ='); rad2dms(lB);
    fprintf('aBA ='); rad2dms(aBA);

%% PPØÍKLAD 2.2 ===========================================================
% bod B
[fB, lB] = reckon(fA, lA, s, aA, [a sqrt(e2)], 'rad'); % souradnice bodu B
aBA = azimuth(fB, lB, fA, lA, [a sqrt(e2)], 'rad'); % azimut z B na A
    fprintf('\nPR 2. reseni na elipsoidu:\n');
    fprintf(' fB ='); rad2dms(fB);
    fprintf(' lB ='); rad2dms(lB);
    fprintf('aBA ='); rad2dms(aBA);

%% VIZUALIZACE %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
cla
axesm ('gnomonic', 'Grid','on', 'Origin',stred, 'GLineWidth', sit_w );
axis off; hold on % vypnuti zobrazeni os; zakaz prepisovani novymi objekty

% prevod na stupne
ro = 180/pi;

% rovnik
[lat,long] = track1('gc',0,0,90);
plotm(lat,long,polednik,'LineWidth',polednik_w)
% nulty polednik
[lat,long] = track1('gc',0,0,0,[-90 90]);
plotm(lat,long,polednik,'LineWidth',polednik_w)

% ORTODROMA -------------------------------------------------------------
[lat,long] = track1('gc', fA*ro, lA*ro, aA*ro);
plotm(lat,long,krivka,'LineWidth',krivka_w)
% body
plotm(fA*ro, lA*ro, bodA, 'MarkerSize',bod(1),'LineWidth',bod(2))
plotm(fB*ro, lB*ro, bodB, 'MarkerSize',bod(1),'LineWidth',bod(2))
plotm(fP*ro, lP*ro, bodP, 'MarkerSize',bod(1),'LineWidth',bod(2))
plotm(fR*ro, lR*ro, bodR, 'MarkerSize',bod(1),'LineWidth',bod(2))

