I downloaded a netcdf file that shows chlorophyll concentrations. I am trying to plot this information as a raster file over a shapefile map. When I plot the raster layer, I am having trouble with two things
I can't figure out how to display coordinates with the chlorophyll layer - it won't let me turn my lat and lon into a raster layer
I want to overlay the raster layer over a shapefile map ... but it seems like the NA values are also showing a a color, so it won't overlay over anything?
Here is the code I have so far:
library(sp)
library(maptools)
library(RNetCDF)
library(ncdf4)
library(raster)
# retrieving chlorophyll info
attributes(mycdf$var)
attributes(mycdf$var$chlor_a)
# get lat lon info
attributes(mycdf$dim)
nc_lat <- ncvar_get(mycdf, attributes(mycdf$dim)$names[1])
nc_lon <- ncvar_get(mycdf, attributes(mycdf$dim)$names[2])
mycdf <- nc_open(file.choose(), verbose = TRUE, write = FALSE)
rasternew <- raster(chlor)
projection(rasternew)<- CRS("+proj=longlat +ellps=WGS84 +datum=WGS84
+no_defs")
rasternew<- t(rasternew)
# get map
library(rworldmap)
library(rworldxtra)
newmap<- getMap(resolution = "high")
View(newmap)
projection(newmap)<- CRS("+proj=longlat +ellps=WGS84 +datum=WGS84
+no_defs")
# plot
brks <- c(0, 01, 0.2, 0.3, 0.4, 0.5)
col<- c("aliceblue", "lightblue", "green", "yellow", "orange", "red")
par(mfrow=c(1,1), font.axis = 4, mai = rep(.1,4), pin = c(5, 4))
plot(rasternew, breaks = brks, col=col, legend = F)
par(mfrow=c(1,1), new = TRUE, font.axis = 4, mai = rep(.1,4), pin =
c(5, 4))
plot(newmap, xlim=c(150,152), ylim=c(-12, -9), asp = 1, axes = F, col =
"lightgrey", yaxt = "n", xaxt = "n", add = T)
legend(x='bottomleft', legend = c("0", "0.1", "0.2", "0.3", "0.4",
"0.5"), fill = col)
