18. Geografi

Det är mycket vanligt med data där koordinaterna är ”geographics” eller ”latitude/longitude”.

Till skillnad från koordinater i Mercator, UTM eller Stateplane är geografiska koordinater inte kartesiska koordinater. Geografiska koordinater representerar inte ett linjärt avstånd från ett ursprung som plottas på ett plan. Snarare beskriver dessa sfäriska koordinater vinkelkoordinater på en jordglob. I sfäriska koordinater specificeras en punkt av rotationsvinkeln från en referensmeridian (longitud) och vinkeln från ekvatorn (latitud).

_images/cartesian_spherical.jpg

Du kan behandla geografiska koordinater som ungefärliga kartesiska koordinater och fortsätta att göra spatiala beräkningar. Mätningar av avstånd, längd och area kommer dock att vara nonsensiska. Eftersom sfäriska koordinater mäter vinkel avstånd, är enheterna i ”grader” Dessutom kan de ungefärliga resultaten från index och sant/falskt-test som skär och innehåller bli väldigt fel. Avståndet mellan punkterna blir större när man närmar sig problemområden som polerna eller den internationella datalinjen.

Här är till exempel koordinaterna för Los Angeles och Paris.

  • Los Angeles: POINT(-118.4079 33.9434)

  • Paris: POINT(2.3490 48.8533)

Följande beräknar avståndet mellan Los Angeles och Paris med hjälp av PostGIS kartesiska standardkommando: ST_Distance(geometry, geometry). Observera att SRID:en 4326 deklarerar ett geografiskt spatialt referenssystem.

SELECT ST_Distance(
  'SRID=4326;POINT(-118.4079 33.9434)'::geometry, -- Los Angeles (LAX)
  'SRID=4326;POINT(2.5559 49.0083)'::geometry     -- Paris (CDG)
  );
121.898285970107

Aha! 122! Men, vad betyder det?

Enheterna för spatial referens 4326 är grader. Så vårt svar är 122 grader. Men (återigen), vad betyder det?

På en sfär är storleken på en ”gradkvadrat” ganska varierande och blir mindre ju längre bort från ekvatorn man kommer. Tänk på att meridianerna (de vertikala linjerna) på jordklotet kommer närmare varandra när du rör dig mot polerna. Så ett avstånd på 122 grader betyder ingenting. Det är ett nonsensnummer.

För att kunna beräkna ett meningsfullt avstånd måste vi behandla geografiska koordinater inte som ungefärliga kartesiska koordinater utan snarare som verkliga sfäriska koordinater. Vi måste mäta avstånden mellan punkterna som verkliga banor över en sfär - en del av en storcirkel.

PostGIS tillhandahåller denna funktionalitet genom typen geography.

Observera

Olika spatiala databaser har olika metoder för att ”hantera geografi”

  • Oracle försöker dölja skillnaderna genom att på ett transparent sätt göra geografiska beräkningar när SRID:en är geografisk.

  • SQL Server använder två spatiala typer, ”STGeometry” för kartesiska data och ”STGeography” för geografiska data.

  • Informix Spatial är ett rent kartesiskt tillägg till Informix, medan Informix Geodetic är ett rent geografiskt tillägg.

  • I likhet med SQL Server använder PostGIS två typer, ”geometry” och ”geography”.

Med hjälp av typen geography istället för geometry, låt oss försöka igen att mäta avståndet mellan Los Angeles och Paris.

SELECT ST_Distance(
  'SRID=4326;POINT(-118.4079 33.9434)'::geography, -- Los Angeles (LAX)
  'SRID=4326;POINT(2.5559 49.0083)'::geography     -- Paris (CDG)
  );
9124665.27317673

Ett stort tal! Alla returvärden från geography-beräkningar är i meter, så vårt svar är 9125km.

Äldre versioner av PostGIS stödde mycket grundläggande beräkningar över sfären med hjälp av funktionen ST_Distance_Spheroid(point, point, measurement). Men ST_Distance_Spheroid är väsentligt begränsad. Funktionen fungerar bara på punkter och ger inget stöd för indexering över polerna eller den internationella datalinjen.

Behovet av att stödja icke-punktgeometrier blir mycket tydligt när man ställer en fråga som ”Hur nära kommer ett flyg från Los Angeles till Paris att komma Island?”

_images/lax_cdg.jpg

Att arbeta med geografiska koordinater på ett kartesiskt plan (den lila linjen) ger ett mycket felaktigt svar! Att använda storcirkelrutter (de röda linjerna) ger rätt svar. Om vi konverterar vår flygning LAX-CDG till en linjestring och beräknar avståndet till en punkt på Island med hjälp av geography får vi rätt svar (recall) i meter.

SELECT ST_Distance(
  ST_GeographyFromText('LINESTRING(-118.4079 33.9434, 2.5559 49.0083)'), -- LAX-CDG
  ST_GeographyFromText('POINT(-22.6056 63.9850)')                        -- Iceland (KEF)
);
502454.906643729

Så den närmaste inflygningen till Island (mätt från dess internationella flygplats) på LAX-CDG-rutten är relativt små 502 km.

Den kartesiska metoden för att hantera geografiska koordinater bryts ned helt för objekt som korsar den internationella datalinjen. Den kortaste storcirkelrutten från Los Angeles till Tokyo korsar Stilla havet. Den kortaste kartesiska rutten korsar Atlanten och Indiska oceanen.

_images/lax_nrt.png
SELECT ST_Distance(
  ST_GeometryFromText('Point(-118.4079 33.9434)'),  -- LAX
  ST_GeometryFromText('Point(139.733 35.567)'))     -- NRT (Tokyo/Narita)
    AS geometry_distance,
ST_Distance(
  ST_GeographyFromText('Point(-118.4079 33.9434)'), -- LAX
  ST_GeographyFromText('Point(139.733 35.567)'))    -- NRT (Tokyo/Narita)
    AS geography_distance;
 geometry_distance | geography_distance
-------------------+--------------------
  258.146005837336 |   8833954.76996256

18.1. Använda geografi

För att läsa in geometridata i en geograftabell måste geometrin först projiceras till EPSG:4326 (longitud/latitud), sedan måste den ändras till geografi. Funktionen ST_Transform(geometry,srid) konverterar koordinater till geografi och funktionen Geography(geometry) eller suffixet ::geography ”casts” till geografi.

CREATE TABLE nyc_subway_stations_geog AS
SELECT
  ST_Transform(geom,4326)::geography AS geog,
  name,
  routes
FROM nyc_subway_stations;

Att bygga ett spatialt index på en geograftabell är exakt samma sak som för geometri:

CREATE INDEX nyc_subway_stations_geog_gix
ON nyc_subway_stations_geog USING GIST (geog);

Skillnaden finns under täckmanteln: geografiindexet kommer att hantera frågor som täcker polerna eller den internationella datumlinjen korrekt, medan geometriindexet inte gör det.

Här är en fråga för att hitta alla tunnelbanestationer inom 500 meter från Empire State Building.

WITH empire_state_building AS (
  SELECT 'POINT(-73.98501 40.74812)'::geography AS geog
)
SELECT name,
  ST_Distance(esb.geog, ss.geog) AS distance,
  degrees(ST_Azimuth(esb.geog, ss.geog)) AS direction
FROM nyc_subway_stations_geog ss,
     empire_state_building esb
WHERE ST_DWithin(ss.geog, esb.geog, 500);

Det finns bara ett litet antal inbyggda funktioner för geografin:

  • ST_AsText(geografi) returnerar text

  • ST_GeographyFromText(text) returnerar geography

  • ST_AsBinary(geography) returnerar bytea

  • ST_GeogFromWKB(bytea) returnerar geography

  • ST_AsSVG(geografi) returnerar text

  • ST_AsGML(geography) returnerar text

  • ST_AsKML(geography) returnerar text

  • ST_AsGeoJson(geography) returnerar text

  • ST_Distance(geography, geography) returnerar double

  • ST_DWithin(geography, geography, float8) returnerar boolean

  • ST_Area(geography) returnerar double

  • ST_Length(geography) returnerar double

  • ST_Covers(geography, geography) returnerar boolean

  • ST_CoveredBy(geography, geography) returnerar boolean

  • ST_Intersects(geography, geography) returnerar boolean

  • ST_Buffer(geography, float8) returnerar geography [1]

  • ST_Intersection(geography, geography) returnerar geography [1]

18.2. Skapa en geografitabell

SQL för att skapa en ny tabell med en geografikolumn är ungefär som för att skapa en geometritabell. Geografi inkluderar dock möjligheten att ange objekttypen direkt vid skapandet av tabellen. Till exempel

CREATE TABLE airports (
    code VARCHAR(3),
    geog GEOGRAPHY(Point)
  );

INSERT INTO airports
  VALUES ('LAX', 'POINT(-118.4079 33.9434)');
INSERT INTO airports
  VALUES ('CDG', 'POINT(2.5559 49.0083)');
INSERT INTO airports
  VALUES ('KEF', 'POINT(-22.6056 63.9850)');

I tabelldefinitionen specificerar GEOGRAPHY(Point) vår flygplatsdatatyp som punkter. De nya geografiska fälten registreras inte i vyn geometry_columns. Istället registreras de i en vy som heter geography_columns.

SELECT * FROM geography_columns;
          f_table_name    | f_geography_column | srid |   type
--------------------------+--------------------+------+----------
 nyc_subway_stations_geog | geog               |    0 | Geometry
 airports                 | geog               | 4326 | Point

Observera

Vissa kolumner har utelämnats från ovanstående utdata.

18.3. Casting till geometri

Även om de grundläggande funktionerna för geograftyper kan hantera många användningsfall, finns det tillfällen då du kan behöva tillgång till andra funktioner som endast stöds av geometritypen. Lyckligtvis kan du konvertera objekt fram och tillbaka från geografi till geometri.

PostgreSQL-syntaxkonventionen för casting är att lägga till ::typename i slutet av det värde du vill casta. Så, 2::text kommer att konvertera en numerisk två till en textsträng ’2’. Och 'POINT(0 0)'::geometry kommer att konvertera textrepresentationen av point till en geometripunkt.

Funktionen ST_X(point) stöder endast geometritypen. Hur kan vi läsa X-koordinaten från våra geografier?

SELECT code, ST_X(geog::geometry) AS longitude FROM airports;
 code | longitude
------+-----------
 LAX  | -118.4079
 CDG  |    2.5559
 KEF  |  -21.8628

By appending ::geometry to our geography value, we convert the object to a geometry with an SRID of 4326. From there we can use as many geometry functions as strike our fancy. But, remember – now that our object is a geometry, the coordinates will be interpreted as Cartesian coordinates, not spherical ones.

18.4. Varför (inte) använda geografi

Geografiska koordinater är allmänt accepterade koordinater - alla förstår vad latitud/longitud betyder, men väldigt få förstår vad UTM-koordinater betyder. Varför inte använda geografi hela tiden?

  • För det första finns det, som tidigare nämnts, mycket färre funktioner tillgängliga (just nu) som direkt stöder geograftypen. Du kan få ägna mycket tid åt att arbeta runt begränsningar för geograftyper.

  • För det andra är beräkningarna på en sfär beräkningsmässigt mycket dyrare än kartesiska beräkningar. Till exempel innebär den kartesiska formeln för avstånd (Pythagoras) ett anrop till sqrt(). Den sfäriska formeln för avstånd (Haversine) innebär två sqrt()-anrop, ett arctan()-anrop, fyra sin()-anrop och två cos()-anrop. Trigonometriska funktioner är mycket kostsamma och sfäriska beräkningar innehåller många sådana.

Slutsatsen?

Om dina data är geografiskt kompakta (inom en stat, ett län eller en stad), använd geometritypen med en kartesisk projektion som är logisk för dina data. Se webbplatsen http://epsg.io och skriv in namnet på din region för att få ett urval av möjliga referenssystem.

Om du behöver mäta avstånd med en dataset som är geografiskt spridd (täcker stora delar av världen), använd typen geography. Den applikationskomplexitet du sparar genom att arbeta i geography kommer att kompensera för eventuella prestandaproblem. Och casting till geometry kan kompensera för de flesta funktionsbegränsningar.

18.5. Funktionslista

ST_Distance(geometri, geometri): För geometrityp Returnerar det 2-dimensionella kartesiska minimiavståndet (baserat på spatial ref) mellan två geometrier i projicerade enheter. För geograftyp returneras som standard det sfäroida minsta avståndet mellan två geografier i meter.

ST_GeographyFromText(text): Returnerar ett specificerat geografivärde från Well-Known Text representation eller utökad (WKT).

ST_Transform(geometry, srid): Returnerar en ny geometri med dess koordinater transformerade till den SRID som refereras till av heltalsparametern.

ST_X(point): Returnerar punktens X-koordinat, eller NULL om den inte är tillgänglig. Indata måste vara en punkt.

ST_Azimuth(geography_A, geography_B): Returnerar riktningen från A till B i radianer.

ST_DWithin(geography_A, geography_B, R): Returnerar sant om A ligger inom R meter från B.

Fotnoter