Technická příručka k OSM seriálu

Seriál 32. ročníku KSP byl věnován práci s velkými daty a některé jeho díly pracovaly s daty z projektu OpenStreetMap (neboli OSM). Tato stránka poskytuje doplňuje některé technické detaily, na které v textu seriálu nebylo místo, a nabízí odkazy ke stažení dat.

Stažení dat

OSM mapy

Data z OpenStreetMap si lze opatřit dvěma způsoby:
  • Dotazem na API, kterému se zadá obdélník a OSM server pošle jako odpověď výřez. Toto se hodí, když chceme malý kousek, větší server odmítne vydat. Používá se běžně na editaci kousků map.
  • Z projektu Planet OSM. Tam se dají stáhnout mapová data pro celý svět nebo výřezy pro celé státy. Pokud by tě zajímala jedna ulice, tak budeš zbytečně stahovat gigabajty. Ale my si teď chceme hrát s mapou celé Evropy, takže zvolíme tuto možnost.

API i Planet OSM vám dají data ve stejném formátu, takže vás o moc neochudíme, když budeme řešit jen Planet OSM. Na Planet OSM si můžete stáhnout mapu celého světa, zajímá vás Latest Weekly Planet XML File. Pravděpodobně ale nechcete celý svět, v takovém případě se podívejte na Planet.osm na OSM Wiki. Tam najdete seznam serverů, na kterých je Planet OSM zrcadlený a kde mají i výřezy. Mezi zajímavé patří:

  • Stabilní německý Geofabrik poskytuje výřez asi pro každou zemi s docela pěkným rozhraním.
  • Nově existující francouzský OSM download nemá tolik formátů, ale zase má pro některé státy i podregiony.
  • Pro Českou republiku může být zajímavý zdroj OSM archiv VUT v Brně včetně historie.

SRTM výšková data

Výšková data z projektu SRTM (více jsme ho popsali v textu druhé série) se dají opatřit z více zdrojů, tím hlavním jsou asi stránky americké USGC, ale není vůbec jednoduché tam ta data najít. Proto může být dobrým prvním krokem podívat se na aktuální seznam odkazů na stránce SRTM na OSM Wiki. Data jsou typicky distribuovaná v souborech po "čtvercích" o straně jednoho úhlového stupně pojmenovaných podle svého jihovýchodního rohu (například N50E014 obsahuje data od bodu 50°N 14°E do 51°N 15°E).

Můžete se setkat s několika verzemi a několika přesnostmi SRTM dat. Verze 1 jsou prakticky surová data z původního snímání, ve kterých jsou občas scházející údaje a různé chyby měření (když se například nepovedlo změřit výšku konkrétního bodu v horách kvůli divným odrazům nebo oblačnosti, případně se radar nechal zmást rovnými plochami s nízkou odrazivostí v pouštích nebo naopak mořské vlny započítal jako pevninu). Verze 2 jsou vyčištěná data z verze 1 oříznutá podle map pobřeží, ale stále v něm na pevnině schází některé oblasti, které se nepovedlo nasnímat. Nakonec nejpoužívanější verze 3 obsahuje data z verze 2 doplněná nějakou extrapolací o scházející data, takže by všechny body na pevnině měly mít udanou výšku. Dnes se prakticky pro mapové účely používají jen data verze 3 a my nebudeme výjimkou.

Data byla původně dostupná jen v rozlišení 3 úhlových vteřin (což je asi 90 metrů na rovníku), až později byla uvolněna data v rozlišení 1 úhlové vteřiny (asi 30 metrů na rovníku). My budeme používat ta přesnější.

Nakonec zmíníme ještě několik externích agregátorů, které buď odkazují na správné soubory na USGC, nebo na nějaké jejich kopie v jiném formátu:

  • Na stránkách http://srtm.csi.cgiar.org/srtmdata/ lze najít dobrou klikací mapu, výběr z několika formátů (GeoTIFF a Esri ASCII) a různě velké oblasti ke stáhnutí. Nově mají data s rozlišením jedné úhlové vteřiny.
  • Na stránkách https://dwtkns.com/srtm30m/ je dobře udělané stahovátko, které vám nabízí odkazy na správné soubory ze stránek earthdata.nasa.gov. Ke stažení si jenom musíte na earthdata.nasa.gov vytvořit účet.

Data použitá v seriálu jsme postahovali přes druhý zmíněný agregátor. Jsou to data v binárním HGT formátu, který jsme popsali na záložce Parsování SRTM.

Testovací data k úlohám

Pro účely řešení našich úloh jsme pro v době 32. ročníku připravili konkrétní výřezy, abychom mohli verifikovat výsledky proti stejným datům. Ty nyní již zastaraly a tak vám doporučíme spíše stažení dat pro Evropu a Českou republiku z výše zmíněných serverů (například z Geofabrik.de). Výřezy pro Hrochův Týnec nebo Brno si pak musíte z dat pro Českou republiku pořídit sami, například pomocí nástroje osmium-extract. Všechna data, která jsme používali, byla stažena k datu 1. 9. 2019.

SRTM binární data:

Stáhněte si níže zmíněné čtverce z jednoho ze zdrojů SRTM dat s přesností jedné úhlové vteřiny, HGT soubory jsou pojmenované podle svého jihozápadního rohu. Více detailů k formátu souborů je na záložce Parsování SRTM.

  • Česká republika – čtverce N48E012 až N51E018
  • Praha – čtverce N49E014 a N50E014
  • Brno – čtverec N49E016
  • Hrochův Týnec – čtverec N49E015

OSM XML

Z Planet OSM si můžete stáhnout data ve dvou formátech - zakódovaná v XML a nebo v PBF. XML je poměrně často používaný textový formát, zkratka znamená "eXtensible Markup Language", česky něco jako "rozšiřitelný značkovací jazyk". PBF je binární formát a jeho název je zkratkou za "protocol buffer". Oba jsou to obecné jazyky pro zápis dat a potkáme se s nimi nejen při práci s OSM ale i na mnoha dalších místech. I když je PBF formát úspornější na velikost dat, tak se v seriálu budeme držet formátu XML, protože si myslíme, že s ním bude méně starostí a ukážeme si na něm více zajímavých věcí.

XML dokument je poskládaný z elementů a atributů – pokud jste někdy viděli HTML, tak je to velmi podobné. Elementy a jejich atributy mají vždy jméno, atributy mají libovolnou textovou hodnotu, elementy mají atributy a obsah. Jednoduchý XML dokument může vypadat nějak takto:


	<?xml version='1.0' encoding='UTF-8'?>+
	<osm version="0.6" generator="osmconvert 0.8.10" timestamp="2019-08-17T20:15:02Z">
		<node id="344942189" lat="49.9743305" lon="15.8992473" version="2" timestamp="2009-10-01T22:39:31Z" changeset="0"/>
		<node id="344942429" lat="49.9676795" lon="15.9035589" version="4" timestamp="2015-08-08T14:59:30Z" changeset="0">
			<tag k="name" v="Stíčany"/>
			<tag k="place" v="village"/>
			<tag k="source" v="csu:uir-zsj"/>
			<tag k="name:cs" v="Stíčany"/>
			<tag k="ref:zsj" v="048313"/>
			<tag k="ref:cobe" v="048313"/>
			<tag k="population" v="300"/>
		</node>
	</osm>
	

Soubor začíná hlavičkou, která říká akorát verzi XML a kódování. Nic moc zajímavého. Pak je v XML souboru právě jeden element, v našem případě se jmenuje osm (špičaté závory jej ohraničují, ty nejsou součástí jména). Za jménem má element atributy version, generator a timestamp s hodnotami za znakem = v uvozovkách. Uvnitř elementu osm (mezi značkami <osm ...> a </osm> vidíte další tagy, tentokrát node. Ty tvoří obsah elementu osm, elementy takto můžou další elementy a celý XML dokument je tak vlastně jenom strom.

Asi už snadno vidíte, že každý node má atributy id, lat, lon, version, timestamp, changeset a může mít v sobě libovolný počet elementů tag.

Za zmínku ještě stojí syntaxe <node ... /> – značka, která obsahuje na konci lomítko, definuje tag, který nemá obsah. Je to tedy vlastně jenom zkratka za <node ...></node>.

XML nám ale vůbec neříká jaké elementy smíme používat, jaké mají mít atributy a co vlastně znamenají. To je užitečné, protože se tak dá jeden formát použít skoro na cokoliv, ale je potřeba jeho použití ještě doplnit další dokumentací. V případě OSM jí najdete na OSM Wiki.

Gzip a bzip2

Protože XML je velmi "ukecaný" formát (obsahuje pořád dokola názvy elementů a atributů), tak se mapová data alespoň komprimují algoritmem gzip nebo bzip2. To bude drobná komplikace při zpracování, ale pro většinu jazyků najdete snadno knihovnu, která bude umět dekomprimovat gzip streamově (důvody a princip streamového zpracování jsme popsali v textu prvního dílu seriálu).

Volba jazyka a ukázky parsování

K řešení můžete zvolit libovolný svůj oblíbený programovací jazyk. Jen si ověřte, že má alespoň nějaké rozumné knihovny pro práci s XML soubory a pro jejich streamové parsování – bez toho se vám úlohy na velkých datech asi vyřešit nepovede.

My jsme úlohy ozkoušeli na třech jazycích a to na Pythonu 3, Go a C#. V Pythonu využíváme knihovnu lxml, což je vlastně jen Pythoní napojení na C knihovnu libxml2, a díky tomu je i Python rozumně rychlý na zpracování velkých XML. I tak ale Python dosahuje (podle našich hrubých měření) sotva třetinového výkonu třeba oproti řešení v Go. Pokud se tedy rozhodnete řešit úlohy v něm, tak vám doporučujeme si vyhradit ještě více času na to, aby vám doběhl výpočet i na velkém vstupu.

Níže můžete najít ukázky streamového parsování ve všech třech jmenovaných jazycích. Pro ukázku zde zkoušíme najít všechny vrcholy i cesty, které označují nějakou knihovnu (mají tag amenity=library):


		

		

		

Parsování SRTM

SRTM data se dají získat v mnoha formátech, ale ten nejvíce klasický (a ten, který vám dáváme ke stažení na záložce Stažení OSM/SRTM dat) je v podobě binárních HGT souborů (typicky ještě zabalených v zipu).

Každý HGT soubor udává výšky pro "čtverec" o straně jednoho úhlového stupně a první část formátu je již název samotného souboru – soubory jsou pojmenované podle svého jihozápadního (levého dolního) rohu, například N50E014 obsahuje data od bodu 50°N 14°E do 51°N 15°E.

Když využíváme data s přesností jedné úhlové vteřiny, tak takový soubor obsahuje 3601x3601 výškových bodů. Jak vám právě asi došlo, tak data v jednotlivých HGT souborech jsou včetně okrajů – což také znamená, že dva sousedící soubory se budou překrývat právě v jednom řádku/sloupci, myslete na to při parsování!

Samotný soubor obsahuje výšky zapsané jako 16-bitová čísla, hezky jedno číslo za druhým (postupně od západu k východu a ze severu na jih, takže prvních 3601 hodnot souboru N50E014 bude udávat výšky bodů od 51°N 14°E až do 51°N 15°E). Mimo výšek v souboru nic jiného není (a velikost souboru je tak přesně 3601·3601·2 neboli 25934402 bytů).

Pro ukázku načítání jsme si pro vás připravili malé prográmky, které dostanou na vstupu dvě celá čísla, načtou čtverec odpovídající těmto číslům a vypíší z něj souřadnice nejvyššího bodu. Zkuste si programy spustit například na čtverci N50E015.