Chcę wykreślić granice kraju Ameryki Północnej na obrazie rastrowym przedstawiającym pewną zmienną, a następnie nałożyć wykres na wierzch wykresu za pomocą R. Udało mi się to zrobić przy użyciu podstawowej grafiki i sieci, ale wydaje się, że proces kreślenia jest o wiele za wolno! Nie zrobiłem tego jeszcze w ggplot2, ale wątpię, czy będzie lepiej, jeśli chodzi o szybkość.
Mam dane w pliku netcdf utworzonym z pliku grib. Na razie pobrałem granice kraju dla Kanady, USA i Meksyku, które były dostępne w plikach RData z GADM, które wczytują R jako obiekty SpatialPolygonsDataFrame.
Oto kod:
# Load packages
library(raster)
#library(ncdf) # If you cannot install ncdf4
library(ncdf4)
# Read in the file, get the 13th layer
# fn <- 'path_to_file'
r <- raster(fn, band=13)
# Set the projection and extent
p4 <- "+proj=lcc +lat_1=50.0 +lat_2=50.0 +units=km +x_0=32.46341 +y_0=32.46341 +lon_0=-107 +lat_0=1.0"
projection(r) <- CRS(p4)
extent(r) <- c(-5648.71, 5680.72, 1481.40, 10430.62)
# Get the country borders
# This will download the RData files to your working directory
can<-getData('GADM', country="CAN", level=1)
usa<-getData('GADM', country="USA", level=1)
mex<-getData('GADM', country="MEX", level=1)
# Project to model grid
can_p <- spTransform(can, CRS(p4))
usa_p <- spTransform(usa, CRS(p4))
mex_p <- spTransform(mex, CRS(p4))
### USING BASE GRAPHICS
par(mar=c(0,0,0,0))
# Plot the raster
bins <- 100
plot(r, axes=FALSE, box=FALSE, legend=FALSE,
col=rev( rainbow(bins,start=0,end=1) ),
breaks=seq(4500,6000,length.out=bins))
plot(r, legend.only=TRUE, col=rev( rainbow(bins,start=0,end=1)),
legend.width=0.5, legend.shrink=0.75,
breaks=seq(4500,6000,length.out=bins),
axis.args=list(at=seq(4500,6000,length.out=11),
labels=seq(4500,6000,length.out=11),
cex.axis=0.5),
legend.args=list(text='Height (m)', side=4, font=2,
line=2, cex=0.8))
# Plot the borders
# These are so slow!!
plot(can_p, add=TRUE, border='white', lwd=2)
plot(usa_p, add=TRUE, border='white', lwd=2)
plot(mex_p, add=TRUE, border='white', lwd=2)
# Add the contours
contour(r, add=TRUE, nlevel=5)
### USING LATTICE
library(rasterVis)
# Some settings for our themes
myTheme <- RdBuTheme()
myTheme$axis.line$col<-"transparent"
myTheme$add.line$alpha <- 1
myTheme2 <- myTheme
myTheme2$regions$col <- 'transparent'
myTheme2$add.text$cex <- 0.7
myTheme2$add.line$lwd <- 1
myTheme2$add.line$alpha <- 0.8
# Get JUST the contour lines
contours <- contourplot(r, margin=FALSE, scales=list(draw=FALSE),
par.settings=myTheme2, pretty=TRUE, key=NULL, cuts=5,
labels=TRUE)
# Plot the colour
levels <- levelplot(r, contour=FALSE, margin=FALSE, scales=list(draw=FALSE),
par.settings = myTheme, cuts=100)
# Plot!
levels +
layer(sp.polygons(can_p, col='green', lwd=2)) +
layer(sp.polygons(usa_p, col='green', lwd=2)) +
layer(sp.polygons(mex_p, col='green', lwd=2)) +
contours
Czy istnieje sposób na przyspieszenie kreślenia wielokątów? W systemie, nad którym pracuję, kreślenie zajmuje kilka minut. Chcę w końcu stworzyć funkcję, która z łatwością wygeneruje pewną liczbę tych wykresów do kontroli, i sądzę, że wykreślę wiele z tych map, więc chcę zwiększyć prędkość wykresów!
Dzięki!