Készítette: Laky Piroska, BME Általános és Felsőgeodézia Tanszék Frissítve: 2014.04.14.
MÉRŐÁLLOMÁSSAL MÉRT POLÁRIS MÉRÉSEK FELDOLGOZÁSA POLÁRIS PONT SZÁMÍTÁSÁHOZ SZÜKSÉGES FÜGGVÉNYEK Mérőállomással egyszerre tudjuk mérni a ferde távolságot, a vízszintes és a magassági vagy zenitszöget. Ezekből a mérési eredményekből meghatározhatóak a részletpontok koordinátái, feltéve, hogy tudjuk az álláspont koordinátáit és a tájékozási szöget. Nézzük meg, hogy milyen függvények szükségesek ezek kiszámításához. FOK-PERC-MÁSODPERC FÜGGVÉNY BŐVÍTÉSE (FFF.PPMM ALAKKAL) Az egyik fontos tudnivaló, ha Matlab/Octave segítségével szeretnénk szögfüggvényeket használni, hogy vagy radiánban használjuk a szögeket (pl. sin, cos, tan, atan, atan2 függvények) vagy tizedfokban (ekkor sind, cosd, tand függvények). Fok-perc-másodpercet alapból nem kezel a se a Matlab, se az Octave (igaz a Matlab Mapping toolboxában vannak ilyen függvények is). Az első kiegyenlítő számítások gyakorlaton már előkerült ez a probléma, ahol írtunk is egy függvényt, ami fok-perc-másodpercből átszámol tizedfokba, akár vesszővel vannak elválasztva az értékek, akár vektorban(mátrixban) vannak tárolva. Van azonban egy másik gyakori megadási mód, amikor fff.ppmm formában adják meg az értékeket, áldecimális értékként, ahol a tizedespont előtt van a fok, az első két tizedes a perc, a második két tizedes értéke pedig a másodperc értéke. Pl: 312.0456 jelentése 312 fok 04 perc 56 másodperc. Egészítsük ki a korábbi függvényünket, hogy ilyen bemenetet is értelmezni tudjon. Ehhez egy plusz feltételt kell beiktatni, ha egy bemenet van, annak egy vagy 3 oszlopa van-e (size(f,2))? Ha egy, akkor áldecimálisként adott a fok-perc-másodperc, ha 3, akkor vektorként. function [fok rad]= fpm(f,p,m) % A fuggveny fok-perc-masodperc ertekekbol tizedfokot és radiant szamol. % A bemenet lehet vesszovel elvalasztva fok,perc, masodperc vagy % lehet vektorban megadva [fok perc másodperc] alakban is, vagy % aldecimalis fok ertekben fff.ppmm alakban. % Matrixban/vektorban adott bemenetet is elfogad a fuggveny. % 1. kimenet: tizedfok (skalar vagy vektor bemenettol fuggoen) % 2. kimenet: a szog erteke radianban switch nargin case 3 fok = f + p/60 + m/3600; rad = fok * pi / 180; case 1 if size(f,2)==1 fok = fix(f); f = (f-fok)*100; perc = fix(f); f = (f-perc)*100; mp = fix(f); fok = fok + perc/60 + mp/3600; rad = fok * pi / 180; end if size(f,2)==3 fok = f(:,1) + f(:,2)/60 + f(:,3)/3600; rad = fok *pi / 180;
end end end
MEGJELENÍTÉS FFF-PP-MM FORMÁTUMBAN Másik hasznos program lehet, ha az eredményt meg tudjuk jeleníteni a geodéziában szokásos fff-pp-mm formában. Ehhez alakítsuk át a szintén korábban megírt függvényünket, ami tizedfokból fok-perc-másodpercbe alakított át. Most a kimenet egy formázott string lesz. function str = fok2fpmstr(x); % A fuggveny tizedfokbol fok-perc-masodperc ertekekbe szamol at. % A kimenet egy formázott string: fff-pp-mm alakban f = fix(x); p = fix((x-f) * 60); m = fix(((x-f)*60-p)*60); str = sprintf('%03d-%02d-%02d', f, p, m); end
A poláris pontok számításához szükség lesz a tájékozási szög ismeretére, ehhez írjunk egy külön függvényt. A tájékozáshoz szükséges irányszöget számítani a második geodéziai főfeladattal, ezt célszerű külön függvényben megírni. Két kimenete legyen mindkét függvénynek, az irányszög és a távolság. Középtájékozási szögeknél a távolság ismerete is szükséges. IRÁNYSZÖG, TÁVOLSÁG SZÁMÍTÁS Az irányszög számításához használhatjuk a hagyományos atan, vagyis arcus tangens függvényt, ekkor figyelni kell a különböző szögnegyedeket, vagy használhatjuk az atan2 függvényt, ami már figyelembe veszi a szögnegyedeket, és az eredményt –π és +π között adja vissza radiánban. Geodéziában a 0-360 fok (vagyis 0-2π ) között tartományban dolgozunk, ezért ha negatív az eredmény, adjunk hozzá 2π-t, majd váltsuk át tizedfokba. A függvénynek két bemenete lesz, P1(y1,x1) és P2(y2,x2) pont és két kimenete d - irányszög, t - távolság. function [d t] = iranyszog(P1,P2) % 2. geodéziai főfeladat, az irányszöget tizedfokban kapjuk y1 = P1(1); x1 = P1(2); y2 = P2(1); x2 = P2(2); dy = y2 - y1; dx = x2 - x1; d = atan2(dy,dx); if sign(d)==-1 d = d + 2*pi; end d = d * 180 / pi; t = sqrt(dy^2+dx^2); end
2
TÁJÉKOZÁS A tájékozáshoz szükségesek az álláspont koordinátái (A), a tájékozó irány koordinátái (T) és a leolvasott irányérték a tájékozó irányra (l). Kimenet a tájékozási szög (z) és a távolság(t). function [z t] = tajekozas(A, T, l) % tajekozas szamitasa % bemenet: % A - allaspont koordinatai % T - tajekozoirany koordinatai % l - iranyertek (tizedfokban) % kimenet: % z - tajekozasi szog % t - tavolsag sulyozashoz [d t] = iranyszog(A,T); z = d - l; if sign(z)==-1 z = z + 360; end end
POLÁRIS PONT SZÁMÍTÁSA A poláris pont számításához szükséges az álláspont koordinátáinak ismerete, vízszintes távolság és irányszög ismerete. Ezekből megkaphatjuk a P pont koordinátáit (yP,xP) function yA = dy = dx = P(1) P(2) end
P = polaris(A,t,d) A(1); xA = A(2); t * sind(d); t * cosd(d); = yA + dy; = xA + dx;
TRIMBLE M3 MÉRŐÁLLOMÁS ADATAINAK FELDOLGOZÁSA A fenti függvények általános függvények, bármilyen formátumú mérések feldolgozásához használhatóak. Nézzük meg, hogy például a tanszéken általánosan használt Trimble M3 műszerből hogyan tudjuk kinyerni a mérési adatokat és feldolgozni! Először a műszer fromátumát kell megismerni. Ehhez használhatjuk az interneten fellelhető leírásokat. A műszer egyik kimeneti formátuma az M5 formátum. most ezzel fogunk kicsit ismerkedni, beolvasni, feldolgozni.
3
TRIMBLE M5 FORMÁTUM A Trimble M5 formátum leírása benne van a Trimble M3 műszer angol nyelvű használati utasításban. Az idevágó részt feltettem a tárgy honlapjára, onnan le lehet tölteni. Íme egy példa egy részletre egy M5 formátumú fájlból: For For For For For For For For For For For For For For For For
M5|Adr M5|Adr M5|Adr M5|Adr M5|Adr M5|Adr M5|Adr M5|Adr M5|Adr M5|Adr M5|Adr M5|Adr M5|Adr M5|Adr M5|Adr M5|Adr
00080|TI 00081|PI1 00082|PI1 00083|TI 00084|PI1 00085|TI 00086|TI 00087|TI 00088|PI1 00089|PI1 00090|PI1 00091|PI1 00092|TI 00093|TI 00094|PI1 00095|PI1
KN STAT A A S INPUT INPUT ST-OUT !
INPUT POLAR
típus info 1. blokk típusa 1. blokk értéke 2. blokk típusa 2. blokk értéke 3. blokk típusa 3. blokk értéke
| 107|Y 107| | 3|Y |th |th | 550|Y 550|SD 550|dl 550|dy |th | 951|SD 951|Y
576676.580 m
577210.871 m 1.500 m 2.000 m 577188.620 50.218 0.003 -0.011 1.600
m m m m m
26.349 m 577192.302 m
| |X |Hz |Om |X |ih |ih | |X |Hz |dc |dx | | |Hz |X
| |Z |V1 | |Z | | | 188926.115 m | 333.4248 DMS |V1 -0.011 m |dr -0.002 m | | | 229.3221 DMS |V1 188865.255 m |Z 190128.530 217.0159 119.4651 188881.092 1.485 1.485
m DMS DMS m m m
karakter példa szöveges típusra példa pont típusra 18-20 ’TI ’ ’PI1’ 22-48 ’KN STAT’, ’POLAR’, ’ST OUT’, ’INPUT’ pontkód, pontszám 50-51 ’th’ ’Y ’, ’SD’ 53-66 2.000 577060.939, 24.659 73-74 ’ih’ ’X ’, ’Hz’ 76-89 1.485 188805.047, 349.3453 96-97 ’Z ’, ’V1’ 99-112 143.709, 94.0623
4
| 187.285 m | 88.5541 DMS | | 179.471 m | | | | | 90.0951 DMS | 0.012 m | | | | 112.0856 DMS | 169.422 m |
Látszik, hogy oszlopok szerint vannak az adatok (egy sor 121 karakter), a leírásban pontosan meg van határozva, hogy egy sorban hányadik karaktertől hányadik karakterig mi található. A fontosabb részek | vonallal vannak elválasztva, így 6 oszlopot különböztethetünk meg. A fájl beolvasásánál használhatjuk ezeket az | elválasztókat pl. textscan paranccsal: fid = fopen('reszletm.m5'); s=textscan(fid,'%s%s%s%s%s%s','Delimiter','|'); fclose(fid);
Vagy használhatjuk a megadott karaktereket a sorok beolvasásánál. Nézzük meg részletesebben ez utóbbit! MÉRÉSEK BEOLVASÁSA Nézzünk egy olyan sort, ahol mérések vannak tárolva (lásd például előző oldalon az utolsó előtti, For M5|Adr 00094 -es sort! Itt a típusnál (18-21. karakter) PI1 szerepel (point identification), utána az infóban pontszám (951), majd az első blokkban ferde távolság (SD – slope distance), a következő blokkban vízszintes szög (áldecimális fff.ppmm alakban), majd az utolsó blokkban zenitszög (szintén fff.ppmm alakban). Írjunk egy függvényt, ami ebből a sorból kiveszi a szükséges mérési eredményeket, számmá alakítja őket és a szögeket pedig fok-pec-másodperc értékből tizedfokká! Hogy tesztelni tudjuk töltsük le a tárgy oldaláról az M5 mérési állományt. Kezdjünk el egy új m fájlt és egy s1 stringbe másoljunk be egy olyan sort, amiben csak mérések vannak! Egy s2 stringbe pedig egy méréseket tartalmazó sort egy ismert álláspont megadásából ’KN STAT’, ahol a pontszámon kívül pontkód is van! clc; clear all; page_screen_output(0); s1 = 'For M5|Adr 00205|PI1 |Hz 204.1534 DMS |V1 s2 = 'For M5|Adr 00082|PI1 |Hz 217.0159 DMS |V1
1005|SD 103.5359 DMS | ' A 88.5541 DMS | '
20.121 m
107|
Írjuk meg, vagy töltsük le az általános dokumentumok közül a kiegészített fpm.m függvényt (fok-perc-másodpercből tizedfok és radián) és a fok2fpmstr.m fájlt, tizedfokból fok-perc-másodperc függvény megjelenítéshez (string)! Ezután írhatjuk a függvényt, figyelembe véve a dokumentációban megadott karakterszámokat! Figyeljünk, hogy nincs mindig minden blokk kitöltve! Ahol nincs érték, ott a függvény NaN-t ad vissza (Not a Number). Az ismert álláspontnál tájékozást végeztünk, itt például nincs távolságmérés, de van a pontszám előtt egy pontkód (A). A strtrim paranccsal levághatjuk egy string elejéről és végéről a felesleges nullákat, a strtok paranccsal két részre vághatjuk a szóközöknél a szöveget. A str2num szöveget számmá alakít. Az isempty eldönti, hogy egy változó üres-e vagy sem, ezt használhatjuk arra, hogy megnézzük a pontszámon kívül pontkód is tartozik-e az adott ponthoz a strtok használatával. A strcmp összehasonlít két szöveget, hogy egyezik-e, ezzel tudjuk megnézni, hogy van-e ferde távolság, vízszintes szög, zenitszög rögzítve. 5
function [SD Hz V psz] = meresek(l) SD=NaN; Hz=NaN; V=NaN; pid=NaN; type = l(18:20); psz = strtrim(l(22:48)); [a b] = strtok(psz); if isempty(b) psz = str2num(a); else psz = str2num(b); end type1 = l(50:51); type2 = l(73:74); type3 = l(96:97); if strcmp(type1,'SD') SD=str2num(l(53:66)); end if strcmp(type2,'Hz') n = str2num(l(76:89)); Hz = fpm(n); end if strcmp(type3,'V1') n = str2num(l(99:112)); V = fpm(n); end end
Teszteljük az eredményt! [SD Hz V psz] = meresek(s1) [SD Hz V psz] = meresek(s2)
Az eredmény: SD = 20.120999999999999 Hz = 2.042594444444445e+002 V = 1.038994444444445e+002 psz = 1005 SD = NaN Hz = 2.170327777777778e+002 V = 88.928055555555559 psz = 107
KOORDINÁTÁK BEOLVASÁSA Most olvassunk be koordinátákat. A teszteléshez most is vegyünk ki két sort, amiben koordináták szerepelnek, legyen olyan, amiben van pontkód és olyan, amiben nincs! s3 = 'For M5|Adr 00088|PI1 |X 188926.115 m | s4 = 'For M5|Adr 00095|PI1 |X 188865.255 m |Z
!
550|Y
577188.620 m
951|Y
577192.302 m
| ' 169.422 m
6
| '
A függvény: function [yxz psz] = koord(l) y=NaN; x=NaN; z=NaN; psz = strtrim(l(22:48)); [a b] = strtok(psz); if isempty(b) psz = str2num(a); else psz = str2num(b); end type1 = l(50:51); type2 = l(73:74); type3 = l(96:97); if strcmp(type1,'Y ') y = str2num(l(53:66)); end if strcmp(type2,'X ') x = str2num(l(76:89)); end if strcmp(type3,'Z ') z = str2num(l(99:112)); end yxz = [y x z]; end
Teszteljük az eredményt! [yxz psz] = koord(s3) [yxz psz] = koord(s4)
Az eredmény: yxz = 1.0e+005 * 5.771886200000000 psz = 550 yxz = 1.0e+005 * 5.771923020000000 psz = 951
1.889261150000000
NaN
1.888652550000000
0.001694220000000
A MÉRÉSI ÁLLOMÁNY BEOLVASÁSA ÉS A FELDOLGOZÁS Most egy olyan mérési állományt fogunk beolvasni, ahol csak egy álláspontból történtek mérések. Egy ismert álláspont lesz (’KN STAT’) és utána a mérések (’POLAR’). Közben figyelni kell, hogy változik-e a jelmagasság (’INPUT’, ’th’)! Ha több álláspontunk lenne, akkor értelemszerűen folyamatosan figyelnünk kell, hogy mikor van új álláspont rekord, új mérések. Az álláspont rekord beolvasása után számítsunk tájékozást! A Trimble M3 műszer egyébként számolja és rögzíti is a tájékozást (’0m’), de más műszerek nem. Most legalább lehetőségünk van ellenőrizni is a számításainkat.
7
A Trimble műszer másik érdekessége, hogy álláspont rekord után nem a sima irányértékeket rögzíti, hanem a tájékozott irányértékeket, ezért nincs szükség arra, hogy hozzáadjuk a tájékozási szöget az irányértékekhez, lehet rögtön a ’Hz’-ben tárolt értékeket felhasználni a poláris pont számításához. A ferde távolságokat viszont ne felejtsük el vízszintesre redukálni! A műszer egyébként rögzíti a koordinátákat is, így elvileg nem is kellene nekünk semmit számítani, de a gyakorlatban előfordulhat, hogy utólag jövünk rá, hogy rossz prizmaállandó volt beállítva, rossz jelmagasság, és akkor vissza kell nyúlnunk a mérésekig, hogy ki tudjuk javítani a hibát. így most a rögzített koordinátákat szintén csak ellenőrzésre használjuk, hogy jó-e a számításunk. clc; clear all; close all; fid = fopen('meres.m5'); page_screen_output(0); %% Álláspont % Megkeressük az első álláspont rekordot while feof(fid)==0 l=fgetl(fid); tipus = l(18:20); info = strtrim(l(22:48)); if strcmp(tipus,'TI ') && strcmp(info,'KN STAT') break end end % beolvasunk 4 sort, amiben az álláspont rekordot tárolják l1 = fgetl(fid); % tájékozó irány koordinátái l2 = fgetl(fid); % mérések a tájékozó irányra l3 = fgetl(fid); % tájékozási szög l4 = fgetl(fid); % álláspont koordinátái [T T_psz] = koord(l1); [SD Hz V psz] = meresek(l2); [A A_psz] = koord(l4); z = tajekozas(A,T,Hz); fprintf('\nTajekozas\n'); fprintf('Allaspont: %d (Y=%.2f; X=%.2f)\n', A_psz, A(1), A(2)); fprintf('Tajekozo irany: %d (Y=%.2f; X=%.2f)\n', T_psz, T(1), T(2)); fprintf('Iranyertek: %s\n', fok2fpmstr(Hz)); fprintf('Iranyszog: %s\n', fok2fpmstr(iranyszog(A,T))); fprintf('Tajekozasi szog: %s\n', fok2fpmstr(z)); %% Mérések % Közben jelmagasság, műszermagasság változás figyelembe vétele % Megkeressük, hogy hol kezdődnek a mérések, közben ha van beállítjuk a % jelmagasságot (th), műszermagasságot(ih) th=NaN; ih=NaN; pontok = []; while feof(fid)==0 l=fgetl(fid); tipus = l(18:20); info = strtrim(l(22:48)); if strcmp(tipus,'TI ') && strcmp(info,'POLAR') break end if strcmp(tipus,'TI ') && strcmp(info,'INPUT')
8
type1 = l(50:51); type2 = l(73:74); if strcmp(type1,'th') th = str2num(l(53:66)); end if strcmp(type2,'ih') ih = str2num(l(76:89)); end end end fprintf('Muszermagassag: %.2f; jelmagassag: %.2f\n', ih, th) while feof(fid)==0 l=fgetl(fid); tipus = l(18:20); info = strtrim(l(22:48)); type1 = l(50:51); type2 = l(73:74); if strcmp(tipus,'TI ') && strcmp(info,'INPUT') if strcmp(type1,'th') th = str2num(l(53:66)); end if strcmp(type2,'ih') ih = str2num(l(76:89)); end end if strcmp(tipus,'PI1') && strcmp(type1,'SD') [SD Hz V psz] = meresek(l); % A tajekozasi szoget Trimble M3 esetében nem kell hozzaadni az % iranyertekekhez, mert mar eleve tajekozott iranyertekeket tarol % (ha volt elotte Ismert allaspont program) HD = SD * sind(V); P = polaris(A,HD,Hz); Z = A(3) + ih + SD*cosd(V) - th; pont = [psz P Z]; pontok = [pontok; pont]; end end % Eredmenyek kiiratasa fprintf('\nReszletpontok\n'); fprintf(' P Y fprintf('%8d %9.2f %9.2f
X Z\n'); %6.2f\n', pontok');
% Plottolas figure(1); hold off; plot3(pontok(:,2), pontok(:,3), pontok(:,4), 'b+'); fclose(fid); % Fajlba iras fid = fopen('koord_m5.txt','w'); fprintf(fid, ' P Y fprintf(fid, '%5d %9.2f %9.2f fclose(fid);
X Z \r\n'); %6.2f\r\n', pontok');
9
