Блог

У меня есть 2 разные маски, которые я хочу очертить на карте. В основном эти маски имеют размер [180 x 360], маски заполнены 1, а остальное в некоторых случаях равно 0 или NaN. Это файлы, которые я использую для масок: Mask1 и Mask2. Вот код того, что я пробовал:

 library(raster)
library(pracma)



R_Caribe <- read.delim("R_Caribe.txt", header=FALSE)
mask1 <- as.matrix(R_Caribe)

mask1[mask1 == 0] <- NA 


R_CUENCA <- read.delim("R_CUENCA.txt", header=FALSE)
mask2 <- as.matrix(R_CUENCA)

mask2[mask2 == 0] <- NA

###########################################################
#I need to rotate the masks so they can plot properly
############################################################

r_1 = raster(pracma::circshift(mask1, c(0,181)))
extent(r_1) <- c(-180,180,-90,90)

r_2 <- raster(pracma::circshift(mask2, c(0,181)))
extent(r_2) <- c(-180,180,-90,90)



###########################################################
#This is to get the coastlines
############################################################

mex <- getData("GADM", country = "MEX", level = 0)
us <- getData("GADM", country = "USA", level = 0)


###########################################################
# Plot
############################################################
plot(r_1, xlim = c(-130, -20), ylim = c(5 ,35))
plot(r_2,add = TRUE)
plot(mex, add = TRUE)
plot(us, add = TRUE)}
  

Используя этот код, я получаю следующую карту:

Я хочу получить карту, но каждая маска должна иметь определенный цвет, обвести их черной линией и, наконец, указать на них имя, то есть «mask1», «mask2».

Я решил свою проблему, я использовал ответ Роберта отсюда.

     library(raster)
    library(pracma)
    library(sp)
    library(maps)


    #READ THE MASKS

    R_Caribe <- read.delim("E:/PCT/Datos/Mascaras_humedad/R_Caribe.txt", header=FALSE)
    mask1 <- as.matrix(R_Caribe)

    mask1[mask1 == 0] <- NA


    R_CUENCA <- read.delim("E:/PCT/Datos/Mascaras_humedad/R_CUENCA.txt", header=FALSE)
    mask2 <- as.matrix(R_CUENCA)

    mask2[mask2 == 0] <- NA


    R_PACIFICO <- read.delim("E:/PCT/Datos/Mascaras_humedad/R_PACIFICO.txt", header=FALSE)
    mask3 <- as.matrix(R_PACIFICO)

    mask3[mask3 == 0] <- NA


    R_ROFL <- read.delim("E:/PCT/Datos/Mascaras_humedad/R_ROFL.txt", header=FALSE)
    mask4 <- as.matrix(R_ROFL)

    mask4[mask4 == 0] <- NA

#ROTATE
    r_1 = raster(pracma::circshift(mask1, c(0,181)))
    extent(r_1) <- c(-180,180,-90,90)

    r_2 <- raster(pracma::circshift(mask2, c(0,181)))
    extent(r_2) <- c(-180,180,-90,90)

    r_3 <- raster(pracma::circshift(mask3, c(0,181)))
    extent(r_3) <- c(-180,180,-90,90)

    r_4 <- raster(pracma::circshift(mask4, c(0,181)))
    extent(r_4) <- c(-180,180,-90,90)


#CREATE POLYGONS
    rr1 <- r_1 > -Inf
    caribe <- rasterToPolygons(rr1, dissolve = TRUE)


    rr2 <- r_2 > -Inf
    cuenca <- rasterToPolygons(rr2, dissolve = TRUE)


    rr3 <- r_3 > -Inf
    pacifico <- rasterToPolygons(rr3, dissolve = TRUE)


    rr4 <- r_4 > -Inf
    rofl <- rasterToPolygons(rr4, dissolve = TRUE)


    #..............................................................................
    #                           P L O T S
    #..............................................................................
    par(mar=c(4.9, 3.9, 4.1, 9.1), xpd=TRUE)

    plot(caribe, xlim = c(-130, -20), ylim = c(10 ,35), lwd=5, border='black',col = "plum", xlab="Lon", ylab="Lat")
    plot(pacifico,border='black', lwd=5,col = "coral3",add = TRUE)
    plot(rofl, border='black', lwd=5,col = "deepskyblue",add =TRUE)
    plot(cuenca, border='black', lwd=5,col = "yellow",add =TRUE)
    axis(side=1, labels=TRUE)
    axis(2, labels=TRUE)
    legend("topright", inset=c(-0.2,0),legend=c("Caribe", "Pacífico", "R.T.","Cuenca"), pt.bg=c("plum", "coral3","deepskyblue","yellow"), cex=1.5,box.lty=0, pch = 22, col=c("black","black"))
    map("world",add=TRUE,fill=FALSE,xlim = c(-130, -20), ylim = c(10 ,35))
  

И конечный результат таков

Я пытался использовать ggplot, но не мог понять, как его использовать, поэтому, если кто-нибудь знает, дайте мне знать.

   [1]: https://gis.stackexchange.com/questions/187798/create-polygons-of-the-extents-of-a-given-raster-in-r