Prevenzione della mappatura parziale delle linee costiere

Question

Come evitare che i poligoni dei paesi vengano tagliati con diverse proiezioni?

Nell'esempio seguente, mi piacerebbe realizzare una mappa di proiezione stereografica dell'Antartide comprese le latitudini < -45 ° S. Impostando i miei limiti y su questo intervallo, l'area del tracciato è corretta, ma i poligoni del paese vengono anche ritagliati a questi limiti. C'è un modo per fare in modo che le coste siano tracciate ai margini dell'area del tracciato?

Grazie per qualsiasi consiglio.

library(maps) 
library(mapproj) 

ylim=c(-90,-45) 
orientation=c(-90, 0, 0) 

x11() 
par(mar=c(1,1,1,1)) 
m <- map("world", plot=FALSE) 
map("world",project="stereographic", orientation=orientation, ylim=ylim) 
map.grid(m, nx=18,ny=18, col=8) 
box() 

Answer 1

Il problema è dovuto al fatto che si sta specificando un massimo di -45 ylim, ei dati sono tagliati con precisione a quella latitudine. Chiaramente i limiti effettivi della trama sono calcolati in un modo separato, che probabilmente ha a che fare con qualche tipo di assunto prudente basato sui limiti proiettati risultanti (ma di cui non so nulla senza esplorare il codice sorgente).

È possibile evitare il problema attraverso la creazione di una trama fittizia che non attira la costa, quindi aggiungere i dati sulla parte superiore con una maggiore ylim:

library(maps) 
library(mapproj) 

ylim = c(-90,-45) 
orientation=c(-90, 0, 0) 

x11() 
par(mar=c(1,1,1,1)) 
m <- map("world", plot=FALSE) 
map("world",project="stereographic", orientation=orientation, ylim=ylim, col = "transparent") 
map("world",project="stereographic", orientation=orientation, ylim=ylim + c(-5, 5), add = TRUE) 
map.grid(m, nx=18,ny=18, col=8) 
box() 

BTW, questa strategia non funziona per tutte le proiezioni dal momento che alcune implicheranno una sovrapposizione per alcuni limiti di dati, ma la stereografica polare si estende appena fuori dal centro quindi non c'è alcun problema di questo tipo.

Inoltre, ci sarebbe un modo più aggiornato per farlo con maptools, rgdal e rgeos che ti permetterebbe di usare i dati condividendo una PROJ.4 Proiezione stereografica appropriata (ma non ho ancora capito il ritaglio passo). Le proiezioni in mapproj sono "di sola forma" e sarebbe un po 'di lavoro ottenere altri dati in questi, afaik.

Answer 2

Mi sono reso conto che un'altra opzione consisteva nel definire i limiti dell'area di stampa e non nella mappa stessa. Ciò potrebbe dare una certa flessibilità nella definizione dell'area di trama esatta (ad esempio xaxs = "i" e yaxs = "i"). Si garantisce anche che tutti i poligoni vengano tracciati - nell'esempio di @mdsumner, l'Australia è mancante e i limiti y dovrebbero essere estesi per tracciare correttamente il poligono.

orientation=c(-90, 0, 0) 

ylim <- c(mapproject(x=-180,y=-45, project="stereographic", orientation=orientation)$y, mapproject(x=0,y=-45, project="stereographic", orientation=orientation)$y) 
xlim <- c(mapproject(x=-90,y=-45, project="stereographic", orientation=orientation)$x, mapproject(x=90,y=-45, project="stereographic", orientation=orientation)$x) 

x11(width=6,height=6) 
par(mar=c(1,1,1,1)) 
plot(0,0, t="n", ylim=ylim, xlim=xlim, xaxs="i", yaxs="i", xlab="", ylab="", xaxt="n", yaxt="n") 
map("world",project="stereographic", orientation=orientation, add=TRUE) 
map.grid(nx=18,ny=18, col=8) 
box() 

Answer 3

Un altro modo è quello di utilizzare una proiezione PROJ.4 reale e utilizzare i dati impostati nel pacchetto maptools. Ecco il codice, c'è ancora un problema nel continente Antartico dove la linea costiera oceanica incontra il polo a -90 (che è un po 'fastidioso, e dovremmo trovare quella coordinata e rimuoverla dal risultato trasformato, o correre attraverso uno strumento di topologia per trovare lo scheggia).

library(maptools) 
data(wrld_simpl) 
xrange <- c(-180, 180, 0, 0) 
yrange <- c(-90, -45, -90, -45) 

stere <- "+proj=stere +lat_0=-90 +lat_ts=-71 +lon_0=0 +k=1 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs" 

library(rgdal) 

w <- spTransform(wrld_simpl, CRS(stere)) 

xy <- project(cbind(xrange, yrange), stere) 
plot(w, xlim = range(xy[,1]), ylim = range(xy[,2])) 
box()