Working with Globcolour data (Part 2)

Just a quick note to announce that the makeGlobcolourField and isin.convert functions have been added to the sinkr package. In addition, the makeGlobcolourField function now used the ncdf4 package to read the .nc files. Both functions are only set up to deal with the higher resolution 4 km data based on the ISIN grid ("L3b").
The following script is an example of extracting data for the Philippines, and produces a map of mean Chl1 values:

Example script:

# Required packages -------------------------------------------------------
 
library(sinkr)
library(maps)
library(mapproj)
library(ncdf4)
library(Matrix)
library(mapdata)
 
 
# Create field ------------------------------------------------------------
 
dat_path <- getwd() # or location of Globcolour .nc files
 
F1 <- makeGlobcolourField(
  folderPath=dat_path,
  yearRan=c(1997,2015),
  lonRan=c(116,127),
  latRan=c(4,17),
  greps=c(".nc", "GSM")
)
 
 
# Create grid and polygons ------------------------------------------------
 
grd <- isin.convert(grd=as.numeric(colnames(F1)))
polys <- isin.convert(grd=grd$grd, polygons = TRUE)
 
 
 
# Plot --------------------------------------------------------------------
 
breaks <- c(0.01, 0.05, 0.1, 0.2, 0.5, 1, 1.5, 2, 2.5, 3)
png("meanChla.png", res=400, units="in", width=6, height=5)
layout(matrix(1:2,1,2), widths=c(5,1), heights=5, respect = TRUE)
op <- par(cex=1, ps=10)
par(mar=c(4,4,1,1))
map("worldHires", "philippines", col=NA, ylim=range(grd$lat), xlim=range(grd$lon))
axis(1); axis(2)
box()
incl <- which(colSums(!is.na(F1)) / nrow(F1) > 0.5) # more that 50% of months are non-NAs
col <- val2col(colMeans(F1, na.rm = TRUE)[incl], col=jetPal(length(breaks)-1), breaks=breaks)
for(j in seq(incl)){
  polygon(polys[[incl[j]]], col=col[j], border=col[j], lwd=0.001)
}
map("worldHires", add=TRUE, col=8, resolution = 0, fill=TRUE)
par(mar=c(4,0,1,4))
imageScale(1, breaks=breaks, axis.pos=4, col=jetPal(length(breaks)-1), log="y", add.axis=FALSE)
axis(4, at=breaks, las=2)
par(op)
dev.off()

Created by Pretty R at inside-R.org