2014-06-13 15 views
6

Ho scoperto che quando provo a sovrapporre più raster usando la trama (..., add = T) se provo a sovrapporre più di 3 raster insieme il grafico successivo non allinea correttamente i raster.Sovrapposizione trama raster utilizzando trama (..., add = T) porta a disallineamento arbitrario della trama finale

Il mio intento originale era quello di creare una mappa categoriale di landcover modellato in cui l'oscurità del colore che rappresenta una classe di copertura variava in base alla certezza nella nostra proiezione del modello. Per fare ciò, ho creato uno script semplice che passasse in rassegna ogni classe di copertina e lo tracciasse (ad esempio, foresta, colore verde sulla mappa) usando un gradiente di colore da grigio (previsione della foresta a bassa certezza) a colore di copertura completo (ad esempio, verde scuro per le aree sono fortemente previsti). Quello che ho trovato è che usando questo approccio, dopo che la terza copertina è stata aggiunta alla trama, tutti i raster successivi che sono sovrapposti sul tracciato sono arbitrariamente disallineati. Ho invertito l'ordine di pianificazione delle classi di copertura e viene mostrato lo stesso comportamento, il che significa che non è un problema con i raster di classe di copertura individuale. Ancora più sconcertante in Rstudio, quando uso il pulsante dello zoom per ispezionare da vicino la trama finale, il disallineamento peggiora.

Avete qualche idea del perché questo comportamento esiste? Soprattutto, hai qualche soluzione suggerita o soluzioni alternative?

Il codice e i dati sul link sottostante riportano tutti i comportamenti descritti. https://dl.dropboxusercontent.com/u/332961/r%20plot%20raster%20add%20issue.zip Turn plot_gradient = F per vedere come se semplicemente si sommi semplicemente uno stesso raster e si aggiungano i sottoinsiemi in modo sequenziale allo stesso grafico è possibile replicare il problema. Ho già provato a impostare l'estensione della trama del dispositivo (..., ext) ma non ha funzionato. Ho anche controllato e l'estensione di ogni copertina raster è la stessa.

Di seguito è riportata la figura delle classi di copertura disallineate. il tracciamento sul dispositivo jpeg avrà come risultato un'immagine simile (ad esempio, questo non è un problema del rendering Rstudio). enter image description here Stranamente, se zoom della immagine utilizzando Rstudio, il disallineamento è diverso enter image description here Per fare un confronto, questo è come le coperture devono essere allineati correttamente nel paesaggio enter image description here

library(raster) 
library(colorRamps) 
raster_of_classes=raster("C:/r plot raster add issue/raster_of_classes.tif") 
raster_of_certainty_of_classes=raster("C:/r plot raster add issue/raster_of_certainty_of_classes.tif") 
endCols=c("darkorchid4", "darkorange3", "red3", "green4", "dodgerblue4") #colors to be used in gradients for each class 
classes=unique(raster_of_classes) 
minVal=cellStats(raster_of_certainty_of_classes, min) 
tmp_i=1 
addPlot=F 
plot_gradient=F #this is for debug only 
#classes=rev(classes) #turn this off and on to see how last 2 classes are mis aligned, regardless of plotting order 
for (class in classes){ 
    raster_class=raster_of_classes==class #create mask for individual class 
    raster_class[raster_class==0]=NA #remove 0s from mask so they to do not get plotted 
    if (plot_gradient){ 
    raster_of_certainty_of_class=raster_of_certainty_of_classes*raster_class #apply class mask to certainty map 
    }else{ 
    raster_of_certainty_of_class=raster_class #apply class mask to certainty map 
    } 
    endCol=endCols[tmp_i] #pick color for gradient 
    col5 <- colorRampPalette(c('grey50', endCol)) 
    if (plot_gradient){ 
    plot(raster_of_certainty_of_class, 
     col=col5(n=49), breaks=seq(minVal,1,length.out=50), #as uncertainty values range from 0 to 1 plot them with fixed range 
     useRaster=T, axes=FALSE, box=FALSE, add=addPlot, legend=F)  
    }else{ 
    plot(raster_of_certainty_of_class, 
     col=endCol,  
     useRaster=T, axes=FALSE, box=FALSE, add=addPlot, legend=F)  
    } 
    tmp_i=tmp_i+1 
    addPlot=T #after plotting first class, all other classes are added 
} 
+0

Non ho l'energia in questo momento al lavoro attraverso tutta la logica che si utilizza per produrre i molteplici raster che stai tramando, ma una soluzione ragionevole al problema sarà quello di istituire un unico raster codifica le diverse classi e per tracciare * that *. Se devi produrre più raster intermedi, uno per ogni classe, puoi rimetterli insieme con una serie di chiamate a 'raster :: cover (x, y)'. –

+0

In effetti, è facile tracciare una mappa di copertura del suolo con un colore monotono specificato per ciascuna copertina. Tuttavia, poiché mi interessa che ogni copertina sia mappata con un gradiente di colore indicativo della certezza del modello, non sembra possibile un semplice approccio alla mappatura dei colori –

+0

Aha. Io vedo. In teoria dovrebbe essere ancora possibile applicare un insieme di valori raster e una scala di colori corrispondente per farlo. Posso vedere come in pratica, però, sarà molto più difficile che dipingere uno strato dopo l'altro sopra a quelli precedenti. –

risposta

2

Nel dicembre 2013 ho ha postato una domanda su esattamente questo comportamento to the R-sig-geo mailing list, e non ha avuto alcuna risposta utile (oltre alla conferma che si verifica anche con versioni R e sistemi operativi diversi dal mio).

Qui, per la cronaca, è l'esempio riproducibile che ho utilizzato per illustrare il problema. (Si veda la questione legata per qualche spiegazione in più.)

library(maptools) ## Only needs to be installed for example data 
library(raster) 
library(rgeos) 

## Create an example raster 
p <- shapefile(system.file("shapes/co37_d90.shp", package="maptools")) 
p <- p[31,] ## A tall narrow county polygon 
pr <- gDifference(gBuffer(p, width=.01), p) 
r <- rasterize(pr, raster(extent(pr), ncol=100, nrow=100)) 

## These three are properly registered on one another 
plot(r, col="yellow", legend=FALSE) 
plot(r, col="green", legend=FALSE, add=TRUE) 
plot(r, col="grey", legend=FALSE, add=TRUE) 
## All subsequent "layers" are improperly shifted/skewed to right 
plot(r, col="yellow", legend=FALSE, add=TRUE) 
plot(r, col="blue", legend=FALSE, add=TRUE) 
plot(r, col="red", legend=FALSE, add=TRUE) 
plot(r, col="grey20", legend=FALSE, add=TRUE) 
## Following the above, SpatialPolygons are also shifted/skewed 
plot(p, border="red", lwd=2, add=TRUE) 

enter image description here

+0

Mentre sono sollevato nel sentire che non sono l'unico ad avere il problema (cioè, non sono pazzo), non è incoraggiante vedere che altri hanno fallito nel trovare una soluzione. So che in SO la parola 'bug' non viene lanciata leggermente, ma questo certamente odora come uno ... –

+2

Concordato. È sicuramente un bug. Un giorno, se nessun altro lo affronta, ho intenzione di farlo. Sarà qualcosa di basso livello, e la mia ipotesi è che abbia a che fare con un'ottimizzazione volta a permettere al dispositivo di tracciamento di non dover "ricordare" tutto ciò che riguarda il layer dopo il quale viene tracciato su di esso. In qualche modo, durante quel passo, un pezzo del relativo metainfo non viene passato abbastanza a destra: le coordinate che dovrebbero essere espresse rispetto all'intero dispositivo sono espresse in relazione a 'par (" mar ")' o 'par (" omar ") ', o ... qualcosa del genere. (Solo un'intuizione, ma molte cose a riguardo indicano in questo modo.) –

3

ho avuto anche questo problema e risolto chiamando la funzione dei parametri grafica, par(), con un sottoinsieme di parametri, e la maggior parte importante, mettere il new=TRUE nella chiamata par(), non la chiamata plot(), prima di ogni chiamata plot addizionale(). Ad esempio:

png(fullname, 
     width = 3000, 
     height= 3000) 

    # original par() call 
    par(mfrow=c(1,1), cex=3, mar=c(3,3,3,7), bg=bgcol, col=txtcol) 

    # first plot 
    plot(zreate, 
     maxpixels=ncell(zreate), 
     col=qcol, 
     colNA=mapbg, 
     xaxt='n', 
     yaxt='n', 
     ext=map_extent, 
     breaks=tq, 
     bty='n', 
     legend=FALSE) 

    #second plot and par() call 
    par(mfrow=c(1,1), cex=3, mar=c(3,3,3,7), bg=bgcol, col=txtcol, new=TRUE) 
    plot(rt, 
     maxpixels=ncell(rt), 
     col=dcol, 
     legend=FALSE, 
     xaxt='n', 
     yaxt='n', 
     ext=map_extent, 
     bty='n') 

    #third plot and par() call 
    par(mfrow=c(1,1), cex=3, mar=c(3,3,3,7), bg=bgcol, col=txtcol, new=TRUE) 
    plot(r0, 
     maxpixels=ncell(r0), 
     col="#9e9ac8", 
     xaxt='n', 
     yaxt='n', 
     ext=map_extent, #PRENAFILTERING fix 
     bty='n', 
     legend=FALSE) 
+0

grazie mille per la risposta, ma sfortunatamente l'ho provato e non allinea ancora i raster tracciati –

+0

Grazie. Ho avuto questo problema prima e ho provato tutto per risolverlo. La risposta sopra fornita da "the_skua" funziona! È un po 'hacker, ma è sufficiente per iniziare. Assicurati di rimuovere l'argomento 'add = TRUE' dalle successive chiamate di funzione' trama (raster) '. Inoltre, per ognuna delle tue funzioni 'par()' assicurati di impostare 'new = TRUE', ad eccezione del primo. –

+0

Anche questo funziona per me! –

0

Interessante problema. Come probabilmente saprai, image() non sembra avere lo stesso problema ma generalmente rende più brutte le mappe, giusto?

library(raster) 
    library(rgeos) 

    ## Create an example raster 
    p <- shapefile(system.file("shapes/co37_d90.shp", package="maptools")) 
    p <- p[31,] ## A tall narrow county polygon 
    pr <- gDifference(gBuffer(p, width=.01), p) 
    r <- rasterize(pr, raster(extent(pr), ncol=100, nrow=100)) 


    ## These three are properly registered on one another 
    image(r, col="yellow") 
    image(r, col="green", add=TRUE) 
    image(r, col="grey", add=TRUE) 
    ## All subsequent "layers" are also registered 
    image(r, col="yellow", add=TRUE) 
    image(r, col="blue", add=TRUE) 
    image(r, col="red", add=TRUE) 
    image(r, col="grey20", add=TRUE) 
    ## Following the above, SpatialPolygons are no longer shifted/skewed 
    plot(p, border="red", lwd=2, add=TRUE)