Vyrovnání oboustraně připojeného a orientovaného polygonového pořadu pomocí MNČ






zpracoval: Zdeněk Nejedlý
dne: 27.12.2002


Ovládání programu    |    Příklad    |    Zdrojový kód

Soubor "plg.m" je třeba umístit do aktuálního pracovního adresáře (Current Directory) nebo do některého z adresářů uvedených v [File / Set Path...] jinak ho Matlab nenajde a nebude s ním počítat. Do stejného adresáře je vhodné umístit i seznam souřadnic a zápisník s měřenými daty (tedy pokud nedáš přednost ručnímu zadání hodnot)

Ovládání programu
Program se spustí příkazem plg('seznam souřadnic','zápisník měřených dat'), konkrétně třeba plg('ss.txt','zap.txt'). Jméno soubodu je třeba zadat i s apostrofy ('). Pokud jsou vstupní soubory v jiném adresáři než program, je třeba zadat jejich úplnou cestu, například plg('g:/ss.txt','c:\data/zap.txt')

Pokud se program spustí bez parametrů, objeví se nabídka s volbou vstupu

Nacist data ze souboru - 1
      Zadat data rucne - 2

              Typ vstupu :

Po stisknutí 1 můžeš zadat vstupní soubory, ...

Nacist data ze souboru - 1
      Zadat data rucne - 2

              Typ vstupu : 1
 
 
Nacteni dat ze souboru
(jmeno souboru zadavej ve tvaru jmeno.pripona (bez'))
    seznam souradnic = ss.txt
      merene hodnoty = zap.txt

... nebo po stisknutí 2 můžeš zadat vstupní hodnoty (čísla bodů, souřadnice y a x, úhly a délky) ručně

Nacist data ze souboru - 1
      Zadat data rucne - 2

              Typ vstupu : 2
 
 
Pocatecni pripojovaci bod
   cislo = 101
   y [m] = 726982.256
   x [m] = 1026352.842
 
Orientace na bod
   cislo =

Po zadání vstupních dat (ať již ze souboru, nebo ručně) je ještě třeba zvolit, zda se budou počítat apriorní nebo aposteriorní střední chyby a parametry středních elips chyb

    Vypocet APRIORNICH strednich chyb a parametru elips - 1
Vypocet APOSTERIORNICH strednich chyb a parametru elips - 2

              Typ vypoctu :

Poté je proveden vlastní výpočet (objeví se při něm výzva k zadání apriorní střední chyby jednotkové m0, střední chyby délek md a úhlů mw) jehož výsledkem je níže uvedený protokol a obrázek polygonu se středními elipsami chyb


Příklad
Vstupní soubory:
seznam souřadnic (číslo bodu, souřadnice y, souřadnice x)

Čísla bodů mohou obsahovat čísla i písmena, nesmí tam ale být mezery
Desetinným oddělovačem je tečka

Q	733762.80	1030209.84
1	733590.48	1029780.16
8	733105.50	1028711.45
L	733313.52	1028743.33

zápisník měřených hodnot (číslo bodu, úhel na bodě, délka vpřed)

1	249.4886	397.209
2	244.4452	283.242
3	089.1833	379.817
4	240.3127	267.978
5	150.8930	300.268
6	140.1824	254.409
7	149.7205	164.782
8	201.8150

Uvedené hodnoty jsou z příkladu ve skriptech Geodézie 50 - Příklady a návody na cvičení na straně 86 až 91 (výsledné střední chyby a délky poloos elips chyb se liší, protože v příkladu počítají aposteriorní hodnoty)

Výsledkem výpočtu je takovýto protokol ...

 ========================================================================
   VYROVNANI OBOUSTRANE PRIPOJENEHO A ORIENTOVANEHO POLYGONOVEHO PORADU    
 ========================================================================

 
 apriorni stredni chyba jednotkova m0 [mm/cc] = 15
                  stredni chyba delek md [mm] = 15
                   stredni chyba uhlu mw [cc] = 20
 
  
 Dane a merene hodnoty
 ---------------------
 bod            y[m]           x[m]            uhel[g]    delka[m]
 -----------------------------------------------------------------
 Q          733762.8000    1030209.8400
 1          733590.4800    1029780.1600        249.4886   397.209
 2                                             244.4452   283.242
 3                                              89.1833   379.817
 4                                             240.3127   267.978
 5                                             150.8930   300.268
 6                                             140.1824   254.409
 7                                             149.7205   164.782
 8          733105.5000    1028711.4500        201.8150
 L          733313.5200    1028743.3300
 
 
 
 Uzavery (matice podminek):
 --------------------------
  Ux = -48.59 mm
  Uy = +16.03 mm
  Uw = +28.48 cc
 
 

 Kontroly vypoctu
 ----------------
  vTPv = 481.159304238
  -UTK = 481.159304238
 
 Aposteriorni stredni chyba jednotkova m0
 ----------------------------------------
    m0 = 12.66 [mm/cc]
 
 Dosazeni do pretvorenych podminkovych rovnic:
 ---------------------------------------------
    Ux = +0.000000000000007
    Uy = -0.000000000000011
    Uw = +0.000000000000021
 
 
Mezni uzavery na hladine vyznamnosti 95%:
-----------------------------------------
      uzaver            mezni uzaver          vysledek
------------------------------------------------------
   Ux = -48.59 mm     Mx = ±90.52 mm          splneno
   Uy = +16.03 mm     My = ±164.67 mm 	      splneno
   Uw = +28.48 cc     Mw = ±141.42 cc 	      splneno
 
 Vyrovnane merene hodnoty a jejich apriorni stredni chyby
 --------------------------------------------------------
 bod    uhel[g]       m[cc]     delka[m]       m[mm]
 ---------------------------------------------------
 1      249.48979     14.94     397.2085       13.89
 2      244.44554     17.79     283.2465       14.28
 3      89.18315      17.39     379.8101       13.45
 4      240.31219     18.51     267.9745       13.56
 5      150.89198     17.96     300.2607       13.52
 6      140.18118     17.52     254.4028       14.22
 7      149.71959     17.31     164.7808       14.08
 8      201.81444     17.18
 
 
 Vyrovnane souradnice urcovanych bodu, jejich stredni chyby a parametry elips chyb
 ---------------------------------------------------------------------------------
 bod       y[m]           x[m]       my[mm]    mx[mm]    mp[mm]    mxy[mm]   omega[g]       a[mm]     b[mm]
 ----------------------------------------------------------------------------------------------------------
 2    733226.5105    1029621.0979    13.26     10.20     16.73     11.83     73.6064        13.89     9.32
 3    732954.7807    1029701.0449    18.22     14.66     23.39     16.54     77.4553        18.72     14.02
 4    732910.7346    1029323.7974    16.66     16.26     23.28     16.46     64.9828        16.90     16.01
 5    732728.1801    1029127.6236    16.64     17.75     24.33     17.20     194.4874       17.76     16.63
 6    732734.7656    1028827.4351    16.81     13.39     21.49     15.20     128.8071       17.74     12.13
 7    732943.3984    1028681.8566    14.00     4.69      14.77     10.44     92.0260        14.10     4.39


... a obrázek s náčrtem polygonového pořadu se středními elipsami chyb


Polygonový pořad s elipsami chyb


Zdrojový kód programu

function plg(ss,zap)
%
%+------------------------------------------------------------------------+
%|  VYROVNANI OBOUSTRANE PRIPOJENEHO A ORIENTOVANEHO POLYGONOVEHO PORADU  |
%+------------------------------------------------------------------------+
%
%funkci spust ve tvaru PLG('SEZNAM SOURADNIC' , 'ZAPISNIK MERENYCH HODNOT')
%
%jmena souboru zadat vcetne apostrofu '
%
%Pokud funkci spustis bez parametru, zobrazi se nabidka s moznostmi
%zadat vstupni soubory (seznam souradnic a zapisnik merenych hodnot)
%nebo zadat vstupni hodnoty (cisla bodu, souradnice y a x a merene
%uhly a delky) rucne
%
%Vstupni soubory musi byt ve stejnem adresari jako tato funkce (soubor plg.m),
%jinak je treba zadat uplnou cestu, napriklad PLG('G:/SS.TXT','ZAP.TXT')
%
%Pri prvnim pouziti je asi nejlepsi zadat hodnoty rucne a pote se podivat
%na struktutu souboru 'ss.txt' a 'zap.txt'
%
%
%Napsano a odladeno v Matlabu 6 (R12), nezarucuji tudiz, ze to ve starsich
%verzich bude chodit
%
%
%rsc@email.cz



if (nargin == 2)
%nacteni vstupnich hodnot ze zadanych souboru (pokud existuji)
    
%------------------------------------------

elseif (nargin ~= 2)

clc;
disp(
'Nacist data ze souboru - 1');
disp(
'      Zadat data rucne - 2');
c = input(
'\n              Typ vstupu : ');

if c==1
   
disp(
' ');
disp(
' ');
disp(
'Nacteni dat ze souboru');
disp(
'(jmeno souboru zadavej ve tvaru jmeno.pripona (bez''))');
ss  = input(
'    seznam souradnic = ','s');
zap = input(
'      merene hodnoty = ','s');

%------------------------------------------

%rucni zadani vstupnich hodnot
elseif c==2
disp(
' ');
disp(
' ');
%---------------------------------
%pripojovaci body
disp('Pocatecni pripojovaci bod');

ss_c(2,:) = input(
'   cislo = ','s');
d = length(ss_c(2,:));

ss(2,1) = input(
'   y [m] = ');
ss(2,2) = input(
'   x [m] = ');
%---------------------------------
disp(' ');
disp(
'Orientace na bod');

cislo = input(
'   cislo = ','s');

if(length(cislo) > d)
    ss_c(1,length(cislo)) = 
' ';
    d = length(cislo);
end

if
(length(cislo) < d)
    cislo(1,d) = 
' ';
end

ss_c(1,:) = cislo;
ss(1,1) = input(
'   y [m] = ');
ss(1,2) = input(
'   x [m] = ');
%---------------------------------
disp(' ');
disp(
' ');
disp(
'Koncovy pripojovaci bod');

cislo = input(
'   cislo = ','s');

if(length(cislo) > d)
    ss_c(3,length(cislo)) = 
' ';
    d = length(cislo);
end

if
(length(cislo) < d)
    cislo(1,d) = 
' ';
end

ss_c(3,:) = cislo;
ss(3,1) = input(
'   y [m] = ');
ss(3,2) = input(
'   x [m] = ');
%---------------------------------
disp(' ');
disp(
'Orientace na bod');

cislo = input(
'   cislo = ','s');

if(length(cislo) > d)
    ss_c(4,length(cislo)) = 
' ';
    d = length(cislo);
end

if
(length(cislo) < d)
    cislo(1,d) = 
' ';
end

ss_c(4,:) = cislo;
ss(4,1) = input(
'   y [m] = ');
ss(4,2) = input(
'   x [m] = ');


%zapis do souboru "ss.txt"
fid = fopen('ss.txt','w');

for(i=1:4)
    fprintf(fid,
'%s \t %f \t %f\n',ss_c(i,:),ss(i,1),ss(i,2));
end

fclose(fid);
disp(
' ');
disp(
'Zadane hodnoty byly zapsany do souboru "ss.txt"');



%merene uhly a delky
disp(' ');
disp(
' ');
p = input(
'Pocet vrcholu polygonoveho poradu = ');
%---------------------------------
disp(' ');
disp(
'Vrchol 1 (pocatecni pripojovaci bod)');

zap_c(1,:) = ss_c(2,:); 
%cislo prvniho bodu
d = length(zap_c(1,:));

zap(1,1) = input(
'          uhel [g] = ');
zap(1,2) = input(
'   delka vpred [m] = ');
%---------------------------------
for i=2:p-1
    fprintf(1,
'\nVrchol %i\n',i);

    cislo = input(
'             cislo = ','s');

    
if(length(cislo) > d)
        zap_c(1,length(cislo)) = 
' ';
        d = length(cislo);
    
end

    if
(length(cislo) < d)
        cislo(1,d) = 
' ';
    
end
    
    
zap_c(i,:) = cislo;

    zap(i,1) = input(
'          uhel [g] = ');
    zap(i,2) = input(
'   delka vpred [m] = ');
end
%---------------------------------
fprintf(1,'\nVrchol %i (koncovy pripojovaci bod)\n',p);

cislo = ss_c(3,:); 
%cislo posledniho bodu

    
if(length(cislo) > d)
        ss_c(1,length(cislo)) = 
' ';
        d = length(cislo);
    
end

    if
(length(cislo) < d)
        cislo(1,d) = 
' ';
    
end

zap_c(p,:) = cislo;

zap(p,1) = input(
'          uhel [g] = ');


%zapis do souboru "zap.txt"
fid = fopen('zap.txt','w');

for(i=1:p-1)
    fprintf(fid,
'%s \t %f \t %f\n',zap_c(i,:),zap(i,1),zap(i,2));
end

    
fprintf(fid,'%s \t %f\n',zap_c(p,:),zap(p,1));


fclose(fid);
disp(
' ');
disp(
'Zadane hodnoty byly zapsany do souboru "zap.txt"');


clear ss;
clear zap;
clear ss_c;
clear zap_c;

ss = [];
zap = [];
ss_c = [];
zap_c = [];


end    
end

%------------------------------------------
%vlastni nacteni vstupnich dat

if(isempty(ss) == 1) ss = 'ss.txt';
end

if
(isempty(zap) == 1) zap = 'zap.txt';
end


%body


[fid,mes]=fopen(ss,'r');

if (fid < 0)
   fprintf(2,
'PRI NACTENI SEZNAMU SOURADNIC DOSLO K CHYBE\n',ss,mes);
   error(
'')
end

[p1,p2]=fscanf(fid,'%s');
fseek(fid,0,-1);
poc_b = p2/3;

if (rem(p2,3)~=0)
  error(
'Chyba v seznamu souradnic');
end

for 
i=1:poc_b
   pom = fscanf(fid,
'%s',1);
   [m,n] = size(pom); 
   cis(i,1:n) = pom;
   y_dane(i,1) = fscanf(fid,
'%f',1);
   x_dane(i,1) = fscanf(fid,
'%f',1);
end

%merene hodnoty
[f,mes]=fopen(zap,'r');

if f<0
   fprintf(2,
'Soubor se nepodarilo otevrit\n',soubor,mes);
   error(
'')
end

[p1,p2]=fscanf(f,'%s');
fseek(f,0,-1);
poc_mereni = ((p2+1)/3);

if (rem((p2+1),3) ~= 0)
  error(
'Chyba v zapisniku merenych dat');
end

for 
i=1:poc_mereni
   pom = fscanf(f,
'%s',1);
   [m,n] = size(pom); 
   cis_mer(i,1:n) = pom;
   mer_uhly(i,1) = fscanf(f,
'%f',1);
   
if(i ~= poc_mereni)
       mer_delky(i,1) = fscanf(f,
'%f',1);
   
end
end


p = poc_mereni; %p = pocet vrcholu poradu


%vektor cisel bodu

for i=1:p
    cisla_bodu(i+1,:) = cis_mer(i,:);
end

[dd d] = size(cisla_bodu);

    cislo = cis(1,:);

    
if(length(cislo) > d)
        cisla_bodu(1,length(cislo)) = 
' ';
    
end

    if
(length(cislo) < d)
        cislo(1,d) = 
' ';
    
end
    
    
cisla_bodu(1,:) = cislo;


    cislo = cis(4,:);

    
if(length(cislo) > d)
        cisla_bodu(1,length(cislo)) = 
' ';
    
end

    if
(length(cislo) < d)
        cislo(1,d) = 
' ';
    
end
    
    
cisla_bodu(p+2,:) = cislo;


%******************************************************
typ_m0 = 1; %implicitni hodnota
disp(
' ');
disp(
'    Vypocet APRIORNICH strednich chyb a parametru elips - 1');
disp(
'Vypocet APOSTERIORNICH strednich chyb a parametru elips - 2');
typ_m0 = input(
'\n              Typ vypoctu : ');

%******************************************************


clc;
disp(
'========================================================================');
disp(
'  VYROVNANI OBOUSTRANE PRIPOJENEHO A ORIENTOVANEHO POLYGONOVEHO PORADU  ');
disp(
'========================================================================');
disp(
' ')
disp(
' ')

ro = 200/pi;
roq = 2000/pi;

%apriorni stredni chyby
  
m0 = input('apriorni stredni chyba jednotkova m0 [mm/cc] = ');
 m_d = input(
'                 stredni chyba delek md [mm] = ');
 m_w = input(
'                  stredni chyba uhlu mw [cc] = ');

%m0 = 15;
%m_d = 15;
%m_w = 20;


disp(' ')
disp(
' ')


disp(
'Dane a merene hodnoty');
disp(
'---------------------');
disp(
'bod       y[m]           x[m]            uhel[g]         delka[m]');
disp(
'-----------------------------------------------------------------');
fprintf(1,
'%s \t %.4f \t %.4f\n',cisla_bodu(1,:),y_dane(1,1),x_dane(1,1));
    
if(mer_uhly(1,1) <= 100)
        fprintf(1,
'%s \t %.4f \t %.4f \t \t  %.4f \t %.3f \n',cisla_bodu(2,:),y_dane(2,1),x_dane(2,1),mer_uhly(1,1),mer_delky(1,1));
    
end
    if
(mer_uhly(1,1) >= 100)
        fprintf(1,
'%s \t %.4f \t %.4f \t \t %.4f \t %.3f \n',cisla_bodu(2,:),y_dane(2,1),x_dane(2,1),mer_uhly(1,1),mer_delky(1,1));
    
end
for 
i=2:p-1
    
if(mer_uhly(i,1) <= 100)
        fprintf(1,
'%s \t \t \t \t \t \t \t \t  %.4f \t %.3f\n',cisla_bodu(i+1,:),mer_uhly(i,1),mer_delky(i,1));
    
end
    if
(mer_uhly(i,1) >= 100)
        fprintf(1,
'%s \t \t \t \t \t \t \t \t %.4f \t %.3f\n',cisla_bodu(i+1,:),mer_uhly(i,1),mer_delky(i,1));
    
end
end
    if
(mer_uhly(p,1) <= 100)
        fprintf(1,
'%s \t %.4f \t %.4f \t \t  %.4f\n',cisla_bodu(p+1,:),y_dane(3,1),x_dane(3,1),mer_uhly(p,1));
    
end
    if
(mer_uhly(p,1) >= 100)
        fprintf(1,
'%s \t %.4f \t %.4f \t \t %.4f\n',cisla_bodu(p+1,:),y_dane(3,1),x_dane(3,1),mer_uhly(p,1));
    
end


fprintf(1,'%s \t %.4f \t %.4f\n',cisla_bodu(p+2,:),y_dane(4,1),x_dane(4,1));

disp(
' ');
disp(
' ');
disp(
' ');

disp(
'Uzavery (matice podminek):');
disp(
'--------------------------');

%smernik PQ
sigma_pq = ro*atan2(y_dane(1,1) - y_dane(2,1),x_dane(1,1) - x_dane(2,1));
if(sigma_pq < 0)
    sigma_pq = sigma_pq + 400;
end

%smernik KL    
sigma_kl = ro*atan2(y_dane(4,1) - y_dane(3,1),x_dane(4,1) - x_dane(3,1));
if(sigma_kl < 0)
    sigma_kl = sigma_kl + 400;
end

%priblizne smerniky a souradnicove rozdily
pribl_smer(1,1) = sigma_pq + mer_uhly(1,1);

for i=2:p
    pribl_y(i-1,1) = mer_delky(i-1,1) * sin(pribl_smer(i-1,1)/ro);
    pribl_x(i-1,1) = mer_delky(i-1,1) * cos(pribl_smer(i-1,1)/ro);    

    pribl_smer(i,1) = pribl_smer(i-1,1) + mer_uhly(i,1) - 200;
end

%uzaver sour x
U_x = 1000*(x_dane(2,1) + sum(pribl_x) - x_dane(3,1));
U(1,1) = U_x;
fprintf(1,
'\tUx = %+.2f mm\n', U_x);

%uzaver sour y
U_y = 1000*(y_dane(2,1) + sum(pribl_y) - y_dane(3,1));
U(2,1) = U_y;
fprintf(1,
'\tUy = %+.2f mm\n', U_y);

%uhlovy uzaver
U_w = 10000*(sigma_pq + sum(mer_uhly) - (200*round((sigma_pq + sum(mer_uhly) - sigma_kl)/200)) - sigma_kl);
U(3,1) = U_w;

fprintf(1,
'\tUw = %+.2f cc\n\n\n', U_w);


%matice D
%delky
for i=1:p-1
    D(1,i) = cos(pribl_smer(i,1)/ro);
    D(2,i) = sin(pribl_smer(i,1)/ro);    
end

%uhly
dy = sum(pribl_y);
dx = sum(pribl_x);
for i=p:2*p-2
    
for j=1:i-p+1
        D(1,i) = -dy/roq;
        D(2,i) =  dx/roq;
        D(3,i) = 1;
    
end
dy = dy - pribl_y(j,1);
dx = dx - pribl_x(j,1);
end

D(3,2*p-1) = 1;



%matice vah P
%delky
for i=1:p-1
    P(i,i) = (m0^2)/(m_d^2);
end

%uhly
for i=p:2*p-1
    P(i,i) = (m0^2)/(m_w^2);
end

%matice normalnich rovnic
N = D*inv(P)*D';

%korelaty
K = -inv(N)*U;

%vektor oprav merenzch hodnot
v = inv(P)*D'*K;

%kontroly
disp(' ');
disp(
'Kontroly vypoctu');
disp(
'----------------');
%vTPv
vTPv = v'*P*v;
fprintf(1,
'\tvTPv = %.9f \n', vTPv);

%-UTK
UTK = -U'*K;
fprintf(1,
'\t-UTK = %.9f \n\n', UTK);

%aposteriorni m0
disp('Aposteriorni stredni chyba jednotkova m0');
disp(
'----------------------------------------');
apost_m0 = (vTPv/3)^0.5;
fprintf(1,
'\t  m0 = %.2f [mm/cc]\n\n', apost_m0);

if(typ_m0 == 2) m0 = apost_m0;
end




%zpetne dosazeni do PPR (=0)
disp('Dosazeni do pretvorenych podminkovych rovnic:');
disp(
'---------------------------------------------');
PPR = D*v+U;
fprintf(1,
'\t  Ux = %+.15f \n', PPR(1,1));
fprintf(1,
'\t  Uy = %+.15f \n', PPR(2,1));
fprintf(1,
'\t  Uw = %+.15f \n\n\n', PPR(3,1));

%vyrovnane delky
for i=1:p-1
    vyr_delky(i,1) = mer_delky(i,1) + v(i,1)/1000;
end

%vyrovnane uhly
for i=1:p
    vyr_uhly(i,1) = mer_uhly(i,1) + v(i+p-1,1)/10000;
end

%mezni uzavery
disp('Mezni uzavery na hladine vyznamnosti 95%:');
disp(
'-----------------------------------------');
Mu = 2.5*sqrt(m0^2*N);
disp(
'           uzaver          mezni uzaver       vysledek');
disp(
'------------------------------------------------------');
if(U_x <= Mu(1,1))
    fprintf(1,
'\t Ux = %+.2f mm \t Mx = ±%.2f mm \t splneno\n',U_x,Mu(1,1));
end
if
(U_x > Mu(1,1))
    fprintf(1,
'\t Ux = %+.2f mm \t Mx = ±%.2f mm \t NESPLNENO\n',U_x,Mu(1,1));
end

if
(U_y <= Mu(2,2))
    fprintf(1,
'\t Uy = %+.2f mm \t My = ±%.2f mm \t splneno\n',U_y,Mu(2,2));
end
if
(U_y > Mu(2,2))
    fprintf(1,
'\t Uy = %+.2f mm \t My = ±%.2f mm \t NESPLNENO\n',U_y,Mu(2,2));
end

if
(U_w <= Mu(3,3))
    fprintf(1,
'\t Uw = %+.2f cc \t Mw = ±%.2f cc \t splneno\n\n',U_w,Mu(3,3));
end
if
(U_w > Mu(3,3))
    fprintf(1,
'\t Uw = %+.2f cc \t Mw = ±%.2f cc \t NESPLNENO\n\n',U_w,Mu(3,3));
end


%kovariancni matice vyrovnanych mereni
Qt = inv(P)-inv(P)*D'*inv(N)*D*inv(P);
Mt = (m0^2)*Qt;

%vypis

if(typ_m0 == 2)
    disp(
'Vyrovnane merene hodnoty a jejich aposteriorni stredni chyby');
    disp(
'------------------------------------------------------------');
end    

if 
(typ_m0 ~= 2)
    disp(
'Vyrovnane merene hodnoty a jejich apriorni stredni chyby');
    disp(
'--------------------------------------------------------');
end

disp('bod    uhel[g]       m[cc]     delka[m]       m[mm]');
disp(
'---------------------------------------------------');

for i=1:p-1
    fprintf(1,
'%s \t %.5f \t %.2f \t %.4f \t %.2f\n',cisla_bodu(i+1,:),vyr_uhly(i,1),(Mt(i+p-1,i+p-1)^0.5),vyr_delky(i,1),(Mt(i,i))^0.5);
end

fprintf(1,'%s \t %.5f \t %.2f\n\n\n',cisla_bodu(p+1,:),vyr_uhly(p,1),(Mt(2*p-1,2*p-1)^0.5));


%vyrovnane smerniky a souradnicove rozdily
vyr_smer(1,1) = sigma_pq + vyr_uhly(1,1);

for i=2:p
    vyr_smer(i,1) = vyr_smer(i-1,1) + vyr_uhly(i,1) - 200;
    
    vyr_y(i-1,1) = vyr_delky(i-1,1) * sin(vyr_smer(i-1,1)/ro);
    vyr_x(i-1,1) = vyr_delky(i-1,1) * cos(vyr_smer(i-1,1)/ro);    
end

%vyrovnane souradnice urcovanych bodu
y_vyr(1,1) = y_dane(2,1) + vyr_y(1,1);
x_vyr(1,1) = x_dane(2,1) + vyr_x(1,1);
for i=2:p-2
    y_vyr(i,1) = y_vyr(i-1,1) + vyr_y(i,1);
    x_vyr(i,1) = x_vyr(i-1,1) + vyr_x(i,1);    
end

%matice F
%delky
for i=1:p-2
    
for j=i:p-2
        F(2*j-1,i) = cos(vyr_smer(i,1)/ro);
        F(2*j,i) = sin(vyr_smer(i,1)/ro);        
    
end
end

%uhly
for i=1:p-2
    dy = 0;
    dx = 0;
    
for j=i:p-2
        dy = dy + vyr_y(j,1);
        F(2*j-1,i+p-1) = -dy/roq;
        
        dx = dx + vyr_x(j,1);
        F(2*j,i+p-1) = dx/roq;        
    
end
end

F(3,2*p-1) = 0;

%kovariancni matice souradnic urcovanych bodu
Mx = F*Mt*F';

%elipsy chyb
d=diag(Mx);
mx=d(1:2:
end);
my=d(2:2:
end);
d2=diag(Mx,1);
mxy=d2(1:2:
end);

w=0.5*atan2((2*(mxy)),(mx-my));

a=sqrt(mx.*cos(w).^2+mxy.*sin(2.*w)+my.*sin(w).^2);

b=sqrt(mx.*sin(w).^2-mxy.*sin(2.*w)+my.*cos(w).^2);

for i=1:p-2
    w_g(i,1) = ro*w(i,1);
    
if(w(i,1) < 0)
        w_g(i,1) = 200+ro*w(i,1);
    
end
end











%vypis
if(typ_m0 == 2)
    disp(
'Vyrovnane souradnice urcovanych bodu, jejich aposteriorni stredni chyby a parametry elips chyb');
    disp(
'----------------------------------------------------------------------------------------------');
end    

if 
(typ_m0 ~= 2)
    disp(
'Vyrovnane souradnice urcovanych bodu, jejich apriorni stredni chyby a parametry elips chyb');
    disp(
'------------------------------------------------------------------------------------------');
end

disp('bod       y[m]           x[m]       my[mm]    mx[mm]    mp[mm]    mxy[mm]   omega[g]       a[mm]     b[mm]');
disp(
'----------------------------------------------------------------------------------------------------------');

for i=1:p-2
    fprintf(1,
'%s \t %.4f \t %.4f \t %.2f \t %.2f \t %.2f  \t %.2f  \t %.4f  \t %.2f  \t %.2f \n',cisla_bodu(i+2,:),y_vyr(i,1),x_vyr(i,1),Mx(2*i,2*i)^0.5,Mx(2*i-1,2*i-1)^0.5,(Mx(2*i,2*i) + Mx(2*i-1,2*i-1))^0.5,((Mx(2*i,2*i) + Mx(2*i-1,2*i-1))^0.5)/(2^0.5),w_g(i,1),a(i,1),b(i,1));
end

disp(' ');
disp(
' ');
disp(
' ');
disp(
' ');
disp(
' ');



%graf
%body
x_P(1,1) = x_dane(1,1);
x_P(2,1) = x_dane(2,1);
y_P(1,1) = y_dane(1,1);
y_P(2,1) = y_dane(2,1);

y_T(1,1) = y_dane(2,1);
x_T(1,1) = x_dane(2,1);
for i=1:p-2
    x_T(i+1,1) = x_vyr(i);
    y_T(i+1,1) = y_vyr(i);    
end
y_T(p,1) = y_dane(3,1);
x_T(p,1) = x_dane(3,1);

x_K(1,1) = x_dane(3,1);
x_K(2,1) = x_dane(4,1);
y_K(1,1) = y_dane(3,1);
y_K(2,1) = y_dane(4,1);


hold on;
axis equal;

%polygon
plot(-y_P,-x_P,'b:');
plot(-y_T,-x_T,
'k-.');
plot(-y_K,-x_K,
'b:');

%popis bodu
x = [];
y = [];

x(1,1) = x_dane(1,1);
x(2,1) = x_dane(2,1);
y(1,1) = y_dane(1,1);
y(2,1) = y_dane(2,1);
for i=1:p-2
    x(i+2,1) = x_vyr(i);
    y(i+2,1) = y_vyr(i);    
end
x(p+1,1) = x_dane(3,1);
x(p+2,1) = x_dane(4,1);
y(p+1,1) = y_dane(3,1);
y(p+2,1) = y_dane(4,1);

for i=1:p+2
    text(-y(i,1),-x(i,1),cisla_bodu(i,:),
'color','b','spanWeight','bold','spanAngle','italic')
end

x = [];
y = [];

%elipsy

%zvetseni eli
M = ((sum(mer_delky)/(p-1))/2)/((max(a)+max(b))/2);

%poloosy v meritku
a_M = M*a;
b_M = M*b;
w_M = pi-w; 
%osy GEO -> MAT

%vykresleni
for i=1:p-2
    t=0:0.1:2*pi+0.1;
    x1 = a_M(i,1) * cos(t);
    y1 = b_M(i,1) * sin(t);

    x2 =  x1 * cos(w_M(i,1)) + y1 * sin(w_M(i,1)) + x_vyr(i,1);
    y2 = -x1 * sin(w_M(i,1)) + y1 * cos(w_M(i,1)) + y_vyr(i,1);

    
%vykresleni eli
    
plot(-y2,-x2,'r');
    
    
%osy eli -> a
    
xA = x_vyr(i,1) + a_M(i,1) * cos(w(i,1));
    xB = x_vyr(i,1) - a_M(i,1) * cos(w(i,1));

    yA = y_vyr(i,1) + a_M(i,1) * sin(w(i,1));
    yB = y_vyr(i,1) - a_M(i,1) * sin(w(i,1));
    
    y=[yA,yB];
    x=[xA,xB];
    plot(-y,-x,
'r:');

    
%osy eli -> b
    
xA = x_vyr(i,1) + b_M(i,1) * cos(w(i,1) + pi/2);
    xB = x_vyr(i,1) - b_M(i,1) * cos(w(i,1) + pi/2);

    yA = y_vyr(i,1) + b_M(i,1) * sin(w(i,1) + pi/2);
    yB = y_vyr(i,1) - b_M(i,1) * sin(w(i,1) + pi/2);
    
    y=[yA,yB];
    x=[xA,xB];
    plot(-y,-x,
'r:');
end


title('Polygonový porad a elipsy chyb','spanWeight','bold','spanSize',12,'spanAngle','italic','VerticalAlignment','bottom');

axis 
off;

hold off;




Exapolis, 26.12.2002
Zdeněk Nejedlý
rsc@email.cz