sp :: proj4string setara dalam kerangka PROJ6?

Question

Saya bingung tentang pergeseran dari PROJ4 ke PROJ6 dalam versi (relatif) terbaru dari paket spasial R (dan gdal dalam hal ini). Saya biasa menggunakan sp::proj4string()untuk mendapatkan kode proyeksi / argumen CRS dari suatu spobjek dan menggunakan kode itu lebih jauh dalam skrip saya. Sekarang saya mendapat peringatan (yang memberi tahu saya bahwa proyeksi saya mungkin tidak akurat ). Saya jelas melakukan sesuatu yang salah. Fungsi mana yang harus saya gunakan untuk mendapatkan kode proyeksi tanpa peringatan? Saya pribadi lebih suka kode EPSG jika tingkat kesederhanaan itu mungkin ...

library(sp)

x <- SpatialPolygons(
            list(Polygons(
                 list(Polygon(cbind(c(1,2,3,4,1), c(3,8,10,12,3)))), 
                 ID = 1)), 
            proj4string = CRS(SRS_string = "EPSG:4326"))

proj4string(x)
#> Warning in proj4string(x): CRS object has comment, which is lost in output
#> [1] "+proj=longlat +datum=WGS84 +no_defs"

^{Dibuat pada 2020-09-18 oleh paket reprex (v0.3.0)}

Lihat juga referensi dan komentar dalam pertanyaan ini

Answer 1

Saya juga cukup bingung tentang menjauh dari PROJ Strings dan apa yang harus digunakan sebagai gantinya sampai saya menonton ini video ceramah yang diadakan tahun lalu di FOSS4G.

Inilah kesimpulan utama saya:

Mengapa PROJ dan GDAL menjauh dari string PROJ

Dengan PROJ4 setiap transformasi dilakukan menggunakan "pendekatan-hub". Jadi semuanya pertama-tama akan diubah menjadi WGS84 dan baru kemudian diubah menjadi proyeksi target.

Pendekatan ini dapat menyebabkan kesalahan yang lebih besar dalam proses transformasi daripada transformasi langsung antara dua proyeksi. Dalam kasus terburuk, transformasi yang dibutuhkan bahkan tidak mungkin dilakukan dengan pendekatan ini. Karena string PROJ bisa jadi ambigu, transformasi langsung seringkali tidak memungkinkan.

Contoh dari pembicaraan:

Datum Australia lama GDA94 dan GDA2020 baru keduanya didasarkan pada WGS84. Namun satu menggunakan kerangka referensi ITRF2014 di Epoch 2020.0 dan ITRF92 lainnya di Epoch 1994.0. Perbedaan tersebut tidak dapat ditangkap dalam PROJ-string (sejauh yang saya mengerti) dan dengan demikian semua proyeksi memiliki PROJ String yang hampir sama persis dengan WGS84. Konversi antara keduanya dengan pendekatan hub tidak akan menyebabkan perubahan koordinat. Sementara sebenarnya Australia telah bergerak 1,8 meter antara tahun 1994 dan 2020 dan konversi datum harus mencerminkan hal itu.

Dengan transformasi langsung, kesalahan seperti itu tidak akan terjadi.

Apa yang harus Anda gunakan sebagai gantinya

Sebaliknya Anda idealnya menggunakan definisi yang tidak ambigu untuk proyeksi seperti Kode EPSG atau definisi dalam Format WKT2.

Namun jika Anda tidak bermaksud untuk melakukan transformasi datum (atau tidak peduli tentang beberapa meter atau lebih ketidakakuratan) Anda masih dapat menggunakan PROJ-Strings, tetapi perlu diketahui bahwa praktik terbaik adalah menggunakan Kode EPSG atau WKT2.

Bagaimana ini diterapkan di sp

CRS-objects di sp sekarang memiliki field tambahan commentdi mana representasi WKT2 dari CRS disimpan.

Jadi, alur kerja praktik terbaik spsekarang terlihat seperti ini (diambil dan diperluas sedikit dari sini ):

# Define the CRS using an EPSG Code
x <- CRS(SRS_string='EPSG:4326')

# Display the stored CRS using comment()
cat(comment(x), "\n")

# Store the wkt in a variable
wkt <- comment(x)

# Use this to assign the CRS of another sp-object
y <- CRS(SRS_string = wkt)

Semua dalam semua sfmemiliki fungsi yang jauh lebih nyaman untuk mengakses dan menangani crs daripada sp. Anda dapat membaca lebih lanjut tentang cara sfmenanganinya di tautan yang sama yang diberikan di atas.

Answer 2

Proj4strings sedang dihapus. Berkat paket sf dan pengembangnya ada trik sederhana untuk mengonversi proj4string ke format WKT.

Abaikan saja peringatan datum yang datang dari rgdal dan, nyatanya, Anda bisa mematikannya karena memang menjengkelkan. Tambahkan options("rgdal_show_exportToProj4_warnings"="none")ke file Rprofile.site Anda atau terbitkan sebelum menambahkan pustaka apa pun.

Misalnya dengan proyeksi geografis WGS84, penggunaan sf::st_crs("+proj=longlat +datum=WGS84 +no_defs")akan menghasilkan

Coordinate Reference System:
  User input: +proj=longlat +datum=WGS84 +no_defs 
  wkt:
GEOGCRS["unknown",
    DATUM["World Geodetic System 1984",
        ELLIPSOID["WGS 84",6378137,298.257223563,
            LENGTHUNIT["metre",1]],
        ID["EPSG",6326]],
    PRIMEM["Greenwich",0,
        ANGLEUNIT["degree",0.0174532925199433],
        ID["EPSG",8901]],
    CS[ellipsoidal,2],
        AXIS["longitude",east,
            ORDER[1],
            ANGLEUNIT["degree",0.0174532925199433,
                ID["EPSG",9122]]],
        AXIS["latitude",north,
            ORDER[2],
            ANGLEUNIT["degree",0.0174532925199433,
                ID["EPSG",9122]]]]

Dari blog R-spasial R spasial mengikuti pengembangan GDAL dan PROJ sf negara pengembang:

Penggunaan yang disebut PROJ4-string (seperti + proj = longlat + datum = WGS84) tidak disarankan, mereka tidak lagi menawarkan penjelasan yang memadai tentang sistem referensi koordinat; penggunaan + init = epsg: XXXX mengarah ke peringatan

Answer 3

Peringatan itu tidak memengaruhi output Anda sama sekali. Seperti yang disarankan orang lain, paket sf dapat berguna untuk pengkodean EPSG yang lebih bersih dan sederhana. Ini kode Anda yang diadaptasi:

library(sf)
library(mapview)
coords <- list(rbind(c(1,3), c(2,8), c(3,10), c(4,12), c(1,3)))
x2 <- st_sf(st_geometry(st_polygon(x = coords)), crs="EPSG:4326")
x2$ID <- 1
mapview(x2)