Wednesday, April 30, 2014

Initial Region Masking map

There are many ways to break up the global into regions.
Already spent hours trying to find ways of classification of create regions
Already spent hours trying to locate shape files, geoTiffs, or any GIS compatible sources

Waste of time, the data exists, others have it. Need to work out who to ask and gain the confident to ask.

Need to focus on the process

Josh O'Brien
So from http://stackoverflow.com/questions/20146809/how-can-i-plot-a-continents-map-with-r
SRES are IPCC Emission Scenario Regions


Same code tweaked.

library(rworldmap)
library(rgeos)

sPDF <- getMap()
sres <-
    sapply(levels(sPDF$SRES),
           FUN = function(i) {
               ## Merge polygons within a continent
               poly <- gUnionCascaded(subset(sPDF, SRES==i))
               ## Give each polygon a unique ID
               poly <- spChFIDs(poly, i)
               ## Make SPDF from SpatialPolygons object
               SpatialPolygonsDataFrame(poly,
                                        data.frame(SRES=i, row.names=i))
           },
           USE.NAMES=TRUE)

## Bind the 11 SRES-level SPDFs into a single SPDF
sres <- Reduce(spRbind, sres)

## Plot to check that it worked
plot(sres, col=heat.colors(nrow(sres)))

## Check that it worked by looking at the SPDF's data.frame
## (to which you can add attributes you really want to plot on)
data.frame(sres)

Save for future use. Use  .RData for faster loading in the future.
> class(sres)
[1] "SpatialPolygonsDataFrame"
attr(,"package")
[1] "sp"
> proj4string(sres)
[1] "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"


> setwd("~/Projects/R/test")
> getwd()
[1] "/Volumes/Data/Users/fmacgill/Projects/R/test"

> sres <- data_name
> plot(sres, col=heat.colors(nrow(sres)))
> data_name <- sres
> save(data_name, file="SRESGlobal.RData" )

Clear Variables and check retrieval

> sres <- 1
> data_name <- 1
> load("SRESGlobal.RData")
> sres <- data_name
> plot(sres, col=heat.colors(nrow(sres)))

No comments:

Post a Comment