Jak ArcGIS oblicza odległość między dwoma punktami za pomocą projekcji nierównoodległej?

10

To jest pytanie uzupełniające do mojego poprzedniego. Czy możesz zasugerować dobrze napisane teksty wprowadzające na temat rzutów układu współrzędnych?


Załóżmy, że pracuję z projekcją mapy CH1903, która, o ile wiem, jest zgodna, ale nie w równej odległości. Znaczenie, kąty (kształt) zostały zachowane, ale nie obszary, odległości ani skala. (Przynajmniej nie zostały dokładnie zachowane ). Na razie w porządku.

Zastanawiam się, jakie obliczenia wykonuje ArcGIS, gdy chcę teraz obliczyć odległość między dwoma punktami. W ArcObjects mógłbym używać IProximityOperatorinterfejsu w następujący sposób:

IPoint a = ...,
       b = ...;

double distance = ((IProximityOperator)a).ReturnDistance(b);

Pytanie: Kiedy pracuję z układem odniesienia, który nie zachowuje dokładnie odległości, co ArcGIS zrobiłby, gdybym zapytał go o odległość między dwoma punktami (jak pokazano powyżej)?

  • Czy po prostu wykonuje niektóre matematykę pitagorejską (a 2 + b 2 = c 2 ), aby uzyskać odległość, co oznacza, że ​​zwrócona odległość będzie tak dokładna, jak pozwala na to rzutowanie?

  • A może zrobi coś bardziej skomplikowanego, na przykład jakąś formę ponownej projekcji, aby uzyskać dokładniejszy dystans?

( To samo pytanie, ale bardziej ogólnie: po projekcji geometrii, czy ArcGIS wykonuje wszystkie obliczenia po prostu w przestrzeni euklidesowej, czy też zastosowana projekcja mapy nadal wpływa na obliczenia odległości, kątów, obszarów itp.?)

stakx
źródło
2
Utwórz nowe pytanie zamiast modyfikować oryginał. W przeciwnym razie obalisz wszystkie mechanizmy na tej stronie: co oznaczają oceny, gdy dwa lub więcej pytań jest w grze w jednym wątku? Co oznaczałoby oznaczenie jednej odpowiedzi jako poprawnej? Itd.
whuber
1
@whuber: Chociaż wszystko, co zostało napisane w tym wątku, wciąż dotyczy tematu WRT, zadane pierwotne pytanie, zgadzam się, że teraz są tutaj naprawdę dwa pytania. Jest już za późno, aby to zmienić, ale pamiętaj o swoich radach następnym razem.
stakx

Odpowiedzi:

10

Jeśli chcesz mieć stabilną metodę obliczania odległości geodezyjnych, polecam owijkę Richiego Carmichaela do silnika projekcyjnego ESRI .

Aktualizacja: Właśnie wypróbowałem kod Richiego z ArcGIS 10.0 na Vista64 i otrzymałem wyjątek po wywołaniu LoadLibrary. Zajmę się tym później.

Na razie jednak oto kod w odpowiedzi na pytania w komentarzach do innej odpowiedzi.

Kod porównuje IProximityOperator dla punktów z odniesieniami przestrzennymi i bez nich. Następnie pokazuje, jak użyć azymutalnej projekcji równoodległej (gdzie pierwszy punkt jest punktem styczności), aby znaleźć wielką odległość okręgu.

private void Test()
{
    IPoint p1 = new PointClass();
    p1.PutCoords(-98.0, 28.0);

    IPoint p2 = new PointClass();
    p2.PutCoords(-78.0, 28.0);

    Debug.Print("Euclidian Distance {0}", EuclidianDistance(p1, p2));
    Debug.Print("Distance with no spatialref {0}", GetDistance(p1, p2));

    ISpatialReferenceFactory srf = new SpatialReferenceEnvironmentClass();
    IGeographicCoordinateSystem gcs =
    srf.CreateGeographicCoordinateSystem((int)esriSRGeoCSType.esriSRGeoCS_WGS1984);

    p1.SpatialReference = gcs;
    p2.SpatialReference = gcs;

    Debug.Print("Distance with spatialref {0}", GetDistance(p1, p2));
    Debug.Print("Great Circle Distance {0}", GreatCircleDist(p1, p2));

}
private double GetDistance(IPoint p1, IPoint p2)
{
    return ((IProximityOperator)p1).ReturnDistance(p2);
}

private double EuclidianDistance(IPoint p1, IPoint p2)
{
    return Math.Sqrt(Math.Pow((p2.X - p1.X),2.0) + Math.Pow((p2.Y - p1.Y), 2.0));
}

private double GreatCircleDist(IPoint p1, IPoint p2)
{
    ISpatialReferenceFactory srf = new SpatialReferenceEnvironmentClass();
    IProjectedCoordinateSystem pcs =
    srf.CreateProjectedCoordinateSystem((int)esriSRProjCSType.esriSRProjCS_WGS1984N_PoleAziEqui);
    pcs.set_CentralMeridian(true, p1.X);
    ((IProjectedCoordinateSystem2)pcs).LatitudeOfOrigin = p1.Y;
    p1.SpatialReference = pcs.GeographicCoordinateSystem;
    p1.Project(pcs);
    p2.SpatialReference = pcs.GeographicCoordinateSystem;
    p2.Project(pcs);
    return EuclidianDistance(p1, p2);
}

Oto wynik:

Euclidian Distance 20
Distance with no spatialref 20
Distance with spatialref 20
Great Circle Distance 1965015.61318737

Myślę, że byłoby interesujące przetestować to w stosunku do biblioteki dll silnika projekcyjnego (pe.dll). Opublikuje wyniki, jeśli kiedykolwiek uda mi się uruchomić kod Richiego.

Aktualizacja: po zmianie kodu Richiesa na kompilację dla x86 uruchomiłem go. Ciekawe ... odległość, jaką daje mi wielki okrąg, to 1960273.80162999 - znacząca różnica od tej zwróconej z powyższej metody azymutalnej w równej odległości.

Kirk Kuykendall
źródło
Przyczyną rozbieżności jest prawdopodobnie to, że odcinek linii (w PCS) łączący punkty nie jest rzutowanym geodezyjnym, który byłby krzywoliniowy podczas rzutowania. W związku z tym otrzymujesz mniejszą wartość niż powinieneś. Test tej teorii jest prosty do zrobienia: weź prosty geodezyjny (taki jak równik) i porównaj dwa obliczenia odległości między dwoma szeroko oddzielonymi punktami geodezyjnymi. Jednym z nich jest bezpośrednie obliczenie, jak w kodzie; drugi dzieli element geodezyjny na segmenty, bezpośrednio oblicza długości segmentów i dodaje je. To ostatnie powinno być dokładniejsze.
whuber
9

W ArcGIS 10 sprawdź IGeometryServer2, który ma teraz GetDistanceGeodesic (odległość geodezyjna między dwiema geometriami), GetLengthsGeodesic (zwróć długości geodezyjne każdej polilinii) i DensifyGeodesic (zagęszczenie polilinii poprzez wykreślenie punktów wzdłuż linii geodezyjnych łączących wierzchołki, używa IPolycurve4: : GeodesicDensify).

Jak wspomniano w innych odpowiedziach, ArcGIS nadal wykorzystuje głównie obliczenia płaskie.

Melita Kennedy


Trochę komentarza do innych odpowiedzi (jeszcze za mało powtórzeń, aby komentować bezpośrednio!).

Azymutalna projekcja azymutalna w równej odległości obsługuje elipsoidy. Kod GreatCircleDist tworzy PCS, który wykorzystuje GCS oparty na elipsoidzie / sferoidie, dlatego odległości od środka / punktu początkowego będą odległościami geodezyjnymi, a nie dużymi odległościami od okręgu. Można to również uprościć. Znamy rzutowane współrzędne pierwszego punktu, ponieważ jest to środek rzutu: 0,0. Tak więc tylko drugi punkt musi zostać wyświetlony. Można wtedy zastosować uproszczoną funkcję EuclidianDistance.

Sprawdziłem wyniki pod kątem funkcji geodezyjnych pe.dll i pasowało. Wygląda na to, że aplikacja Richiego używa kuli, więc zwraca duże odległości / współrzędne okręgu w swojej aplikacji testowej. Właśnie dlatego wyniki się nie zgadzają. Nie rozpoznałem wartości promienia; Myślę, że muszę z nim o tym porozmawiać!

Mkennedy
źródło
2
Melita - Świetnie cię tu widzieć!
Kirk Kuykendall
1
Zgadzam się, witam na pokładzie!
matt wilkie
8

Dokładność każdej odpowiedzi na temat ArcGIS może ulec zmianie w dowolnym momencie - o ile wiemy, nowe procedury zostaną wprowadzone w kolejnym pakiecie serwisowym bez ostrzeżenia i dokumentacji. To powiedziawszy, oprogramowanie ESRI od dawna używa obliczeń euklidesowych ( np . Pitagorasa dla odległości), ilekroć używane są współrzędne rzutowane. Często w obliczeniach takich jak ten, który ilustrujesz, oprogramowanie nie ma nawet dostępu do informacji o projekcji, więc co jeszcze można zrobić?

Twoje pytanie wydaje się sugerować, że obliczenia euklidesowe odległości dla projekcji w jednakowej odległości są prawidłowe. Nic nie może być dalej od prawdy. W przypadku jednopunktowego rzutu równoodległego odległość euklidesowa do punktu bazowego jest równa odległości geodezyjnej; w przypadku projekcji dwupunktowej w jednakowej odległości odległość euklidesowa do każdego punktu bazowego jest równa odległości geodezyjnej. W zamian za te gwarancje zniekształcenie metryczne między wszystkimi innymi parami punktów jest zwykle znacznie zwiększone w porównaniu do innych rzutów, które można wybrać.

Whuber
źródło
@whuber: Dzięki za odpowiedź. Odnośnie pierwszego akapitu: doszedłem do wniosku, że ArcGIS może zobaczyć, że używana jest mapa CH1903 (która wykorzystuje elipsoidę Bessela 1841), a następnie rzutuje punkty z powrotem na tę elipsoidę za pomocą układu odniesienia, a następnie wykonuje obliczenia odległości na elipsoidzie. Z twojej odpowiedzi wynika, że ​​ArcGIS nie zrobi tego wszystkiego i pozostanie w przestrzeni euklidesowej XY, aby wykonać obliczenia. (A co z innym oprogramowaniem GIS?) - akapit drugi: Oczywiście masz rację, dziękuję za wyjaśnienie tej kwestii.
stakx
Ukryty mechanizm ponownej rzutowania jest możliwy tylko wtedy, gdy obiekty punktowe zachowują odniesienia do rzutu. Nie wierzę, że tak.
whuber
@ whuber: Czy wystarczyłoby (dla dokładniejszych obliczeń) znać elipsoidę używaną do projekcji? AFAIK, ArcGIS przechowuje odniesienie do użytej projekcji dla każdej klasy obiektów (warstwa danych).
stakx
1
W rzeczywistości IPoint, który wywodzi się z IGeometry, ma SpatialReference jako właściwość. help.arcgis.com/en/sdk/10.0/arcobjects_net/componenthelp/… Jednak nie sądzę, że ReturnDistance go używa. Warto jednak sprawdzić, czy to się zmieniło.
Kirk Kuykendall
1
@stakx Zaktualizowałem swoją odpowiedź, aby uwzględnić kod, który pokazuje, że ustawienie parametru spatialref nie ma wpływu na ReturnDistance.
Kirk Kuykendall