S-JTSK
Článek popisuje základní práci s daty v souřadnicovém systému S-JTSK Systém Jednotné trigonometrické sítě katastrální ve FOSS GIS. Na příkladech ukazuje transformaci běžných geodat mezi vybranými (WGS84, S-JTSK, UTM 33N) souřadnicovými referenčními systémy angl. coordinate reference system (CRS) obvykle používané na území České, případně Slovenské republiky s využitím nástroje z balíků Proj.4 (cs2cs) a GDAL (ogr2ogr, gdalwarp). Dále je součástí článku ukázka nastavení souřadnicových referenčních systémů v GIS aplikacích GRASS a QGIS.
Druhy transformací podle přesnosti (WGS84/ERTS89->S-JTSK) ve 2D:
- 100-110 metrů na severozápad - základní transformace bez transformačního klíče; nepoužívá se
- < 1 metr - sedmiprvková transformace s klíčem "+towgs"; viz níže globální transformační klíč pro ČR/SR
- < 0,1 metr - pomocí gridu (ČR); zvláštní článek
Obsah |
Teorie Křovákova zobrazení
Souřadnicový referenční systém S-JTSK je mimojiné definován vlastním kartografickým zobrazením – tzv. Křovákovým zobrazením.
Křovákovo zobrazení je konformní kuželové zobrazení v obecné poloze, které v roce 1922 navrhl Ing. Josef Křovák. Transformace souřadnic <math>\varphi, \lambda</math> na pravoúhlé <math>X, Y\,</math> se provádí v několika krocích. Nejprve je provedeno Gaussovo konformní zobrazení Besselova elipsoidu na kouli a poté konformní zobrazení na kuželovou plochu obecně položenou.
Schéma zobrazení:
<math>(\varphi, \lambda) \rightarrow (U, V) \rightarrow (\check{S}, D) \rightarrow (\rho, \epsilon) \rightarrow (X, Y)</math>
V první fázi je zobrazen Besselův elipsoid na kouli. Vzhledem k tomu, že tento krok předchází zobrazení na kuželovou plochu, nazývá se Křovákovo zobrazení dvojité. Podmínkou zobrazení je minimální délkové zkreslení kolem základní rovnoběžky, která byla zvolena <math>\varphi_{o}=49^{o}30'</math> (Pro <math>\varphi_{o}</math> je hodnota délkového zkreslení <math>m_{o}=1\,</math>).
Křovák zvolil za zobrazovací plochu kužel v obecné poloze. Empiricky zjistil, že ČSR lze ohraničit horizontálními kružnicemi do pásu o šířce <math>2^{o}31'\,</math>, přičemž zkreslení na okrajích dosahovalo hodnot +24cm/km. V případě normální polohy kužele by byl pás široký <math>3^{o}20'\,</math> a zkreslení v maximální vzdálenosti od střední rovnoběžky až +42 cm/km. Při definitivní úpravě byla zvolena jako základní rovnoběžka kartografická rovnoběžka <math> \check{S}_{o}=78^{o}30'</math> a okrajové <math>\check{S}_{1}=77^{o}13'</math> a <math>\check{S}_{2}=79^{o}44'</math>. Základní rovnoběžka je kolmá na zeměpisný poledník <math>\lambda=42^{o}30'\,</math> východně od Ferra a jejich průsečík <math>A</math> má šířku <math>\varphi=48^{o}15'</math>.
Referenční koule je dále konformně zobrazena na kužel v obecné poloze.
Křovák umístil osu X do obrazu základního poledníku <math>\lambda=42^{o}30'\,</math> v.F. a počátek souřadnic do vrcholu kužele Q. Tím byla celá ČSR umístěna do jediného kvadrantu (viz obrázek níže).
Reference: [1]
Definice CRS v Proj.4/GDAL
Pro práci s CRS používat:
- EPSG a ESRI kódy odkazující do vnitřní knihovny programu
- kombinace kódů a doplňujících parametrů (jen Proj.4)
- definici CRS ze souboru
| název CRS | oficiální EPSG/ESRI kódy | odkazy v Proj.4 | odkazy v GDAL | správná definice dle EPSG/ESRI (formát Proj) |
| S-JTSK East North (Greenwitch) | ESRI:102067 EPSG:5514 |
ESRI:102067 | EPSG:102067 | viz níže |
| S-JTSK South West (Ferro) | EPSG:2065 | chybně EPSG:2065 (bug report) | chybně EPSG:2065 (report) | viz níže |
| WGS84 | EPSG:4326 | EPSG:4326 | EPSG:4326 | +proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs |
| WGS84 / UTM zone 33N | EPSG:32633 | EPSG:32633 | EPSG:32633 | +proj=utm +zone=30 +ellps=WGS84 +datum=WGS84 +units=m +no_defs |
Pozn: Pro méně přesné submetrové transformace je možno WGS84 považovat za identický s ETRS 89 z důvodu jejich velmi malých odlišností.
EPSG:5514/ESRI:102067
Definice S-JTSK v typické verzi pro GIS se zápornými souřadnicemi (osy: East-North, poledník: Greenwich) vč. transformačního klíče pro ČR ve formátu OGC WKT:
PROJCS["Czech GIS S-JTSK (Greenwich) / Krovak",
GEOGCS["Czech S-JTSK (Greenwich)",
DATUM["Czech S-JTSK",
SPHEROID["Bessel 1841",6377397.155,299.1528128,
AUTHORITY["EPSG","7004"]],
TOWGS84[570.8,85.7,462.8,4.998,1.587,5.261,3.56]],
PRIMEM["Greenwich", 0.0, AUTHORITY["EPSG","8901"]],
UNIT["degree",0.0174532925199433]],
PROJECTION["Krovak",
AUTHORITY["EPSG","9819"]],
PARAMETER["latitude_of_center",49.5],
PARAMETER["longitude_of_center",24.83333333333333],
PARAMETER["azimuth",0],
PARAMETER["pseudo_standard_parallel_1",0],
PARAMETER["scale_factor",0.9999],
PARAMETER["false_easting",0],
PARAMETER["false_northing",0],
UNIT["Meter",1]]
Formát Proj (další formáty viz Spatial Reference/102067):
+proj=krovak +lat_0=49.5 +lon_0=24.83333333333333 +alpha=30.28813972222222 +k=0.9999 +x_0=0 +y_0=0 +ellps=bessel +pm=greenwich +units=m +no_defs +towgs84=570.8,85.7,462.8,4.998,1.587,5.261,3.56
EPSG:2065
Definice méně užívané verze S-JTSK s kladnými souřadnicemi (osy: South-West, poledník: Ferro) vč. transformačního klíče pro ČR ve formátu OGC WKT:
PROJCS["Czech S-JTSK (Ferro) / Krovak",
GEOGCS["Czech S-JTSK (Ferro)",
DATUM["Czech S-JTSK",
SPHEROID["Bessel 1841",6377397.155,299.1528128,
AUTHORITY["EPSG","7004"]],
TOWGS84[570.8, 85.7, 462.8, 4.998, 1.587, 5.261, 3.56]],
PRIMEM["Ferro",-17.666666666666668,
AUTHORITY["EPSG","8909"]],
UNIT["degree",0.017453292519943295],
AXIS["Geodetic longitude",EAST],
AXIS["Geodetic latitude",NORTH]],
PROJECTION["Krovak",
AUTHORITY["EPSG","9819"]],
PARAMETER["latitude_of_center",49.5],
PARAMETER["longitude_of_center",42.5],
PARAMETER["azimuth",30.288139722222223],
PARAMETER["pseudo_standard_parallel_1",78.5],
PARAMETER["scale_factor",0.9999],
PARAMETER["false_easting",0.0],
PARAMETER["false_northing",0.0],
UNIT["Meter",-1]]
Formát Proj (další formáty viz Spatial Reference/2065):
+proj=krovak +lat_0=49.5 +lon_0=42.5 +alpha=30.28813972222222 +k=0.9999 +x_0=-0 +y_0=-0 +ellps=bessel +pm=ferro +to_meter=-1 +no_defs +towgs84=570.8,85.7,462.8,4.998,1.587,5.261,3.56
Transformace mezi CRS
V případě potřeby transformace geodat do jiných SRS (např. WGS84 → S-JTSK) je nutné doplnit dvojici CRS o zpřesňující transformační klíč [1]. Tyto parametry byly určeny pro celé území ČR, resp. SR. V případě menšího zájmového území je možné vypočítat lokální transformační klíč, nebo rovnou použít S-JTSK-Grid případná chyba se tak až řádově zmenší. Nepoužití klíče vede ke zhoršení přesnosti o dva řády, přibližně < 100 metrů do S-JTSK, 3 vteřiny do WGS84 v horizontálních souřadnicích. Ve výšce Z versus HAE je chyba právě separace geoidu. Použití zpřesňujícího transformačního klíče nemá vliv na rychlost transformace.
- ČR [2]:
+towgs84=570.8,85.7,462.8,4.998,1.587,5.261,3.56
- SR [3]:
+towgs84=485.0,169.5,483.8,7.786,4.398,4.103,0
Jednotlivé prvky globálního klíče lze definovat přesněji, přínos ve vyšší přesnosti je však minimální: ČR [4], SR [5] Od roku 2011 je v databázi EPSG:2065 díky ČÚZK modifikovaný globální klíč, ani jeho přínos není významný[6].
Testovací bod č. 109 sítě DOPNUL [7]
| WGS84 | S-JTSK | ||||
| zem. délka (lat) | zem. šířka (lon) | výška nad elipsoidem (HAE) | Y | X | Z |
| 12d48'25.15992 12,8069888666667 |
49d27'8.14571 49,4522626972222 |
559.417 | -868208.52 | -1095793.96 | 512.30 |
Testované verze GDAL, Proj.4
Níže uvedené příklady byly testovány:
- Ubuntu 9.04, Proj.4 distribuční verze 4.6.0 (2007/12/21) a GDAL distribuční verze 1.5.2 (2008/05/29) dne 8.-10. května 2009
- Ubuntu 9.04, Proj.4 Les-ejk [8] verze 4.6.1 (2008/09/21) dne 11. května 2009
- Ubuntu 8.04, Proj.4 distribuční verze 4.6.0 (2007/12/21) a GDAL distribuční verze 1.4.4 (2007/11/23) dne 11. května 2009
Transformace textových geodat
Příklady transformace textových geodat za pomoci sady Proj.4 nástroje cs2cs. Výchozí souřadnici Z nebo HAE lze pro výpočty v rovině vynechat.
WGS84 → S-JTSK
echo "12d48'25.15992 49d27'8.14571 559.417" | cs2cs +init=epsg:4326 +to +init=esri:102067 \ +towgs84=570.8,85.7,462.8,4.998,1.587,5.261,3.56 -868208.54 -1095793.58 512.46
Variantní zápis:
echo "12d48'25.15992 49d27'8.14571 559.417" | cs2cs +proj=latlong +datum=WGS84 \ +to +proj=krovak +ellps=bessel +towgs84=570.8,85.7,462.8,4.998,1.587,5.261,3.56 -868208.54 -1095793.58 512.46
Diference: Δ Y = 1cm; Δ X = 38cm; Δ Z = 16cm
S-JTSK → WGS84
echo "-868208.53 -1095793.57 512.30" | cs2cs +init=esri:102067 \ +towgs84=570.8,85.7,462.8,4.998,1.587,5.261,3.56 +to +init=epsg:4326 12d48'25.16"E 49d27'8.146"N 559.261
Variantní zápis:
echo "-868208.53 -1095793.57 512.30" | cs2cs +proj=krovak +ellps=bessel \ +towgs84=570.8,85.7,462.8,4.998,1.587,5.261,3.56 +to +proj=latlong +datum=WGS84 12d48'25.16"E 49d27'8.146"N 559.261
Diference: Δ lat = -1.0E-007 stupně; Δ long = -0.7E-007 stupně; Δ Z = 21cm
S-JTSK → UTM 33N
echo "-868208.53 -1095793.57 512.30" | cs2cs +init=esri:102067 \ +towgs84=570.8,85.7,462.8,4.998,1.587,5.261,3.56 +to +init=epsg:32633 341060.93 5480045.49 559.26
Variantní zápis:
echo "-868208.53 -1095793.57 512.30" | cs2cs +proj=krovak +ellps=bessel \ +towgs84=570.8,85.7,462.8,4.998,1.587,5.261,3.56 \ +to +proj=utm +zone=33N +datum=WGS84 341060.93 5480045.49 559.26
UTM 33N → S-JTSK
echo "341060.93 5480045.49 559.26" | cs2cs +init=epsg:32633 +to +init=esri:102067 \ +towgs84=570.8,85.7,462.8,4.998,1.587,5.261,3.56 -868208.53 -1095793.57 512.30
Variantní zápis je obdobný jako inverzní transformace S-JTSK → UTM 33N.
WGS84 → UTM 33N
echo "12d48'25.15992 49d27'8.14571 559.417" | cs2cs +init=epsg:4326 +to +init=esri:32633 341060.93 5480045.47 559.42
Variantní zápis:
echo "12d48'25.15992 49d27'8.14571 559.417" | cs2cs +proj=latlong +datum=WGS84 \ +to +proj=utm +zone=33N +datum=WGS84 341060.93 5480045.47 559.42
UTM 33N → WGS84
echo "341060.93 5480045.47 559.42" | cs2cs +init=esri:32633 +to +init=epsg:4326 12d48'25.16"E 49d27'8.146"N 559.420
Variantní zápis je obdobný jako inverzní transformace WGS84 → UTM 33N.
Transformace vektorových geodat
Příklady převodu vektorových geodat za pomoci sady GDAL nástroje ogr2ogr. Příklad prezentován na vrstvě typu ESRI Shapefile z datové sady FreeGeodataCZ archiv [2] nebo [3].
S-JTSK → WGS84
ogr2ogr -s_srs "+proj=krovak +lat_0=49.5 +lon_0=24.83333333333333 +alpha=30.28813972222222 \ +k=0.9999 +x_0=0 +y_0=0 +ellps=bessel +pm=greenwich +units=m +no_defs \ +towgs84=570.8,85.7,462.8,4.998,1.587,5.261,3.56" -t_srs "epsg:4326" -f "ESRI Shapefile" cfm_wgs84.shp cfm.shp
Variantní zápis využívající uživatelem vytvořený soubor krovak.prf[9]:
ogr2ogr -s_srs krovak.prf -t_srs "epsg:4326" -f "ESRI Shapefile" cfm_wgs84.shp cfm.shp
WGS84 → S-JTSK
ogr2ogr -s_srs "epsg:4326" -t_srs "+proj=krovak +lat_0=49.5 +lon_0=24.83333333333333 \ +alpha=30.28813972222222 +k=0.9999 +x_0=0 +y_0=0 +ellps=bessel +pm=greenwich +units=m +no_defs \ +towgs84=570.8,85.7,462.8,4.998,1.587,5.261,3.56" -f "ESRI Shapefile" cfm_jtsk.shp cfm.shp
Variantní zápis využívající uživatelem vytvořený soubor krovak.prf[9]:
ogr2ogr -s_srs "epsg:4326" -t_srs krovak.prf -f "ESRI Shapefile" cfm_jtsk.shp cfm.shp
Transformace rastrových geodat
Příklady transformace rastrových geodat za pomoci sady GDAL nástroje gdalwarp. Příklad prezentován na rastru typu tiff + world file z datové sady FreeGeodataCZ archiv [4] nebo [5].
S-JTSK → WGS84
gdalwarp -s_srs "+proj=krovak +lat_0=49.5 +lon_0=24.83333333333333 +alpha=30.28813972222222 \ +k=0.9999 +x_0=0 +y_0=0 +ellps=bessel +pm=greenwich +units=m +no_defs \ +towgs84=570.8,85.7,462.8,4.998,1.587,5.261,3.56" -t_srs "epsg:4326" dem_gtopo30.tif dem_gtopo30_wgs84.tif
Variantní zápis využívající uživatelem vytvořený soubor krovak.prf[9]:
gdalwarp -s_srs krovak.prf -t_srs "epsg:4326" dem_gtopo30.tif dem_gtopo30_wgs84.tif
WGS84 → S-JTSK
gdalwarp -s_srs "epsg:4326" -t_srs "+proj=krovak +lat_0=49.5 +lon_0=24.83333333333333 \ +alpha=30.28813972222222 +k=0.9999 +x_0=0 +y_0=0 +ellps=bessel +pm=greenwich +units=m +no_defs \ +towgs84=570.8,85.7,462.8,4.998,1.587,5.261,3.56" dem_gtopo30.tif dem_gtopo30_jtsk.tif
Variantní zápis využívající uživatelem vytvořený soubor krovak.prf[9]:
gdalwarp -s_srs "epsg:4326" -t_srs krovak.prf dem_gtopo30.tif dem_gtopo30_jtsk.tif
Lokace (location) s CRS S-JTSK v GRASS
Platné pro verze Grass 6.4.0 a vyšší, případně 6.2.3.
Definice projekce GRASS lokace v S-JTSK by měla vypadat přibližně takto:
GRASS> g.proj -p -PROJ_INFO------------------------------------------------- name : Krovak proj : krovak ellps : bessel -PROJ_UNITS------------------------------------------------ unit : meter units : meters meters : 1.0
nebo (EPSG 2065)
-PROJ_INFO------------------------------------------------- name : Krovak proj : krovak datum : hermannskogel ellps : bessel lat_0 : 49.5 lon_0 : 42.5 alpha : 30.28813972222222 k : 0.9999 x_0 : 0 y_0 : 0 pm : ferro no_defs : defined towgs84 : 570.8,85.7,462.8,4.998,1.587,5.261,3.56 -PROJ_UNITS------------------------------------------------ unit : metre units : metres meters : 1
V případě potřeby transformace dat do jiných souřadnicových systémů (např. WGS-84, UTM, S-42) je nutné doplnit transformační klíč (towgs84).
Založení nové GRASS lokace
V GRASS 6.4.0 wxGUI pomocí Location wizard
Pomocí textové rozhraní
Would you like to create location <t1> ? (y/n) [y] <ENTER>
To create a new LOCATION, you will need the following information:
1. The coordinate system for the database
x,y (for imagery and other unreferenced data)
Latitude-Longitude
UTM
Other Projection
2. The zone for the UTM database
and all the necessary parameters for projections other than
Latitude-Longitude, x,y, and UTM
3. The coordinates of the area to become the default region
and the grid resolution of this region
4. A short, one-line description or title for the location
Do you have all this information? (y/n) [y] <ENTER>
Please specify the coordinate system for location <t1>
A x,y
B Latitude-Longitude
C UTM
D Other Projection
RETURN to cancel
D <ENTER>
Other Projection coordinate system? (y/n) [y] <ENTER>
Please enter a one line description for location <t1>
> popis grass lokace
=====================================================
popis grass lokace
=====================================================
ok? (y/n) [y] <ENTER>
Please specify projection name
Enter 'list' for the list of available projections
Hit RETURN to cancel request
> krovak <ENTER>
Do you wish to specify a geodetic datum for this location?(y/n) [y] <ENTER>
Please specify datum name
Enter 'list' for the list of available datums
or 'custom' if you wish to enter custom parameters
Hit RETURN to cancel request
>hermannskogel <ENTER>
Now select Datum Transformation Parameters
Please think carefully about the area covered by your data
and the accuracy you require before making your selection.
Enter 'list' to see the list of available Parameter sets
Enter the corresponding number, or <RETURN> to cancel request
>2 (3 pro SR)
Plural form of units: meters <ENTER>
<ESC><ENTER>
Do you accept this region? (y/n) [y] > <ENTER>
Hit RETURN <ENTER>
<ESC><ENTER>
GRASS> g.proj -p
-PROJ_INFO-------------------------------------------------
name : Krovak
datum : hermannskogel
towgs84 : 570.8,85.7,462.8,4.998,1.587,5.261,3.56
proj : krovak
ellps : bessel
a : 6377397.1550000003
es : 0.0066743722
f : 299.1528128000
-PROJ_UNITS------------------------------------------------
unit : meter
units : meters
meters : 1.0
Volba projekce S-JTSK v QGIS
QGIS od verze 1.0.1 do 1.7.3 obsahuje korektní CRS pro S-JTSK S-JTSK (Greenwich) / Krovak EPSG:102067, ale ve verzi 1.7.4 a 1.8.0 byl dočasně vypuštěn
+proj=krovak +lat_0=49.5 +lon_0=24.83333333333333 +alpha=30.28813972222222 +k=0.9999 +x_0=0 +y_0=0 +ellps=bessel +pm=greenwich +units=m +no_defs
Pro potřeby použití více druhů CRS v jednom mapovém projektu je třeba nadefinovat S-JTSK s konkrétním globálním transformačním klíčem S-JTSK +towgs84:
+proj=krovak +lat_0=49.5 +lon_0=24.83333333333333 +alpha=30.28813972222222 +k=0.9999 +x_0=0 +y_0=0 +ellps=bessel +pm=greenwich +units=m +no_defs +towgs84=570.8,85.7,462.8,4.998,1.587,5.261,3.56
Mapový projekt S-JTSK s vrstvou WGS84
Pokud chcete použít transformaci on-the-fly např. pro projekt v S-JTSK s vloženou vrstvou ve WGS84, musíte definovat mapovému projektu nový CRS identický k ESRI:102067 s transformačním klíčem +towgs84.
Mapový projekt WGS84 s vrstvou S-JTSK
Pokud chcete použít transformaci on-the-fly pro projekt v WGS84 s vloženou vrstvou ve S-JTSK, musíte této vrstvě definovat CRS identický k ESRI:102067 s transformačním klíčem +towgs84.
Odkazy
Související články
Reference
- ↑ Jan Ježek: Vývoj programového modulu pro převod souřadnic mezi kartografickými zobrazeními, diplomová práce ČVUT, Fakulta stavební, Praha 2003
- ↑ Bundesamt für Kartographie und Geodäsie: Description of Transformation - CZ_S-JTSK to ETRS89, 6. prosinec 2004
- ↑ Bundesamt für Kartographie und Geodäsie: Description of Transformation - SK_S-JTSK to ETRS89, 6. prosinec 2004
- ↑ Pavel Mátl: Využití systému WGS84 pro katastrální mapování, str. 54, diplomová práce ZČÚ, Plzeň 2006
- ↑ Geodetický a kartografický ústav Bratislava: Transformácia medzi S-JTSK a ETRS89, SKPOS - Slovenská Priestorová Observacná Služba yužitia signálov GNSS, 7. dubna 2009
- ↑ Coordinate Transformation: S-JTSK to WGS 84 (5) in GeoRepository
- ↑ Zdeněk Hrdina: Transformace souřadnic ze systému WGS-84 do systému S-JTSK, str. 18-21, ČVUT, Praha 1997
- ↑ Jáchym Čepický: Les-ejk UbuntuGIS repository. 2009
- ↑ 9,0 9,1 9,2 9,3 ESRI:102607 in Spatial Reference
Externí odkazy
- Romana Kubátová: Databáze EPSG, semestrální práce z předmětu KMA/PDB, Plzeň 2007
- Romana Kubátová: Systém JTSK a WGS-84, jejich charakteristika a vzájemná transformace, bakalářská práce, Plzeň 2007
- Jan Ježek: Coordinate Reference System Transformations.
- Petr Přidal: Maptilet: CRS Czech Republic
- EPSG registry
