Czy przejście na nową projekcję, a następnie powrót, wpływa na dokładność danych?

Mam klasę obiektów (hrabstwa Karoliny Południowej, czyli dość duży obszar geograficzny) w samolocie stanowym NAD83 SC. Musi zostać przekształcony w drugą projekcję (NAD83 UTM 17), a następnie przekształcony z powrotem w oryginał. W tym celu użyję narzędzia Projekt Esri .

Czy ta podwójna transformacja może spowodować przesunięcie położenia współrzędnych wielokątów i o ile - centymetrów, metrów, kilometrów?

Nie wiem, którego silnika projekcyjnego używa ArcGis, ale bardzo interesujące pytanie dotyczy również proj.4. Dlatego próbuję przetestować silnik projekcyjny proj.4 w środowisku GNU-R. Używam narożników NAD 83 - UTM 17 i EPSG 26917 i rekursywnie ponownie 10000 i 1000000 razy i obliczam różnicę do wartości początkowych.

Oto wyniki:

Wygląda na to, że błąd „ponownej projekcji” mieści się w zakresie centymetrów dla 10000 pętli.

"LON/LAT differences after  10000  loops"
           DLON          DLAT
1 -2.441464e-07 -1.341807e-07
2  2.441129e-07 -1.341807e-07
3  1.852679e-07 -1.691737e-08
4 -1.853157e-07 -1.691819e-08

"X/Y differences after  10000  loops"
            DX           DY
1 -0.025169783 -0.014338141
2  0.025166375 -0.014338208
3  0.002419045 -0.002016762
4 -0.002419690 -0.002016889

I dojdzie do błędu w zakresie metrów, jeśli uruchomisz pętlę 1000000 razy.

"LON/LAT differences after  1000000  loops"
           DLON          DLAT
1 -2.441464e-05 -1.341845e-05
2  2.441128e-05 -1.341846e-05
3  1.852621e-05 -1.691837e-06
4 -1.853105e-05 -1.691828e-06

"X/Y differences after  1000000  loops"
          DX         DY
1 -2.5172288 -1.4339977
2  2.5168869 -1.4340064
3  0.2419201 -0.2017070
4 -0.2419859 -0.2017094

Oto skrypt.

# load the package
require('proj4')

# the LON/LAT frame of NAD83 UTM 17 
lon = c(-84.00, -78.00, -84.00, -78.00 ) 
lat = c( 24.00,  24.00,  83.00,  83.00)

# build the projection conform object
ll0 = matrix(c(lon,lat),nrow=4,ncol=2)
xy0 = project(ll0,"+init=epsg:26917",ellps.default='GRS80')

# make a copy
ll1 = ll0
xy1 = xy0

# number of iterations
num = 1000000

# reproject the stuff num times
for(i in 1:num) {
 # project forward  
 xy1 = project(ll1,"+init=epsg:26917", ellps.default='GRS80')
 # project backward
 ll1 = project(xy1,"+init=epsg:26917", inverse=T, ellps.default='GRS80')
}

# build difference table ll
dll = as.data.frame(ll1-ll0)
names(dll) = c('DLON','DLAT')
# print results LON/LAT
print(paste("LON/LAT differences after ", num," loops"))
print(dll)

# build difference table xy
dxy = as.data.frame(xy1-xy0)
names(dxy) = c('DX','DY')
# print results X/Y
print(paste("X/Y differences after ", num," loops"))
print(dxy)

Dalsze testy w środowisku statystyk powinny być łatwe. Skrypty i objaśnienia kodu dla środowiska linux są dostępne na stronie github.com/bigopensky .

Esri ma własny silnik projekcyjny.

Większość metod prognozowania i transformacji geograficznej / danych jest dobrze zachowana, gdy jest stosowana w odpowiednim obszarze zainteresowania. Jeśli znajdziesz się zbyt daleko poza strefą UTM, poprzeczny Mercator nie zawsze dokładnie „odwraca” (przelicz na długość i szerokość geograficzną) dokładnie. W rzutach używanych dla całego świata mogą występować pewne problemy na biegunach lub wokół biegunów lub południka +/- 180 lub „anty-południka” (południk, który znajduje się naprzeciwko środka rzutowanego układu odniesienia za pomocą współrzędnych).

Przebiegłem 4 punkty, które wypadają poza Karoliną Południową za pośrednictwem silnika projekcyjnego Esri. Aby wykonać test warunków skrajnych na 1k, 10k lub 1M punktów, będę musiał coś zakodować, ponieważ mój istniejący podobny test po prostu wykonuje „podróż w obie strony” - od projekcji geograficznej do projekcji. 32133 to NAD 1983 State Plane South Carolina (metry). 26917 to NAD 1983 UTM zone 17 North.

C:\Users\melita>inverse 32133
382000 20000
      -83.40806392522212        31.98974518135408
382000 383000
      -83.50098893136905        35.26180827475587
839100 20000
      -78.57184097446545        31.98934439195045
839100 383000
      -78.47814111839074        35.26139222680582

C:\Users\melita>forward 26917
  -83.40806392522212        31.98974518135408
       272490.5730967618        3541832.738731374
  -83.50098893136905        35.26180827475587
       272485.6257057797         3904944.98998655
  -78.57184097446545        31.98934439195045
       729409.4734382738        3541830.781689366
  -78.47814111839074        35.26139222680582
       729414.4926270114        3904946.919009762

C:\Users\melita>inverse 26917
 272490.5730967618        3541832.738731374
      -83.40806392522212        31.98974518135408
  272485.6257057797         3904944.98998655
      -83.50098893136905        35.26180827475587
  729409.4734382738        3541830.781689366
      -78.57184097446545        31.98934439195045
  729414.4926270114        3904946.919009762
      -78.47814111839074        35.26139222680582
^Z

C:\Users\melita>forward 32133
  -83.40806392522212        31.98974518135408
                382000.0                  20000.0
  -83.50098893136905        35.26180827475587
                382000.0                 383000.0
  -78.57184097446545        31.98934439195045
                839100.0        19999.99999999814
  -78.47814111839074        35.26139222680582
                839100.0        382999.9999999981

Widać więc, że mieliśmy dwa punkty, które wróciły o 10e-09.

Obsługa w ArcGIS komplikuje fakt, że istnieje odniesienie przestrzenne. Odniesienie przestrzenne obejmuje układ współrzędnych oraz pewne wartości przechowywania i analizy. Domyślnie układy współrzędnych wykorzystujące mierniki są przechowywane z dokładnością do jednej dziesiątej milimetra, 0,0001.

Ujawnienie: Pracuję dla Esri.

Myślę, że jest to przypadek, w którym trzeba przetestować proponowany przepływ pracy pod kątem niektórych funkcji punktów testowych, do których można łatwo dodać pola współrzędnych XY.

Porównaj wartości XY swoich początkowych punktów z tymi, które rzutowałeś / przekształcałeś (jakkolwiek wiele razy), a oszacujesz różnicę.