Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.
The Atmospheric Science Data Center (ASDC) at NASA Langley Research Center offers several data sources. For example, it is possible to download a text file with the 22-year (July 1983 – June 2005) monthly and annual average of global horizontal irradiation.
nasafile <- 'http://eosweb.larc.nasa.gov/sse/global/text/global_radiation' nasa <- read.table(file=nasafile, skip=13, header=TRUE)
With this data, R and the solaR package we can calculate (for example) the global effective irradiation incident on several surfaces: a fixed PV generator, a two-axis tracker and a N-S horizontal tracker. First, let’s plot the original data using some spatial packages:
library(lattice) library(latticeExtra) library(solaR) library(sp) library(maps) library(mapdata) library(maptools) library(gstat) library(hexbin) proj <- CRS('+proj=latlon +ellps=WGS84') coords=nasa[,2:1] datos=nasa[,3:15] nasaSP <- SpatialPixelsDataFrame(points=coords, data=datos, proj4string=proj) world <- map("world", plot=FALSE) llCRS <- CRS("+proj=latlong +ellps=WGS84") world_sp <- map2SpatialLines(world, proj4string=llCRS) paleta=colorRampPalette(rev(brewer.pal('YlOrRd', n=9))) paleta2=colorRampPalette((brewer.pal('RdBu', n=11))) paleta3=heat.ob ##From hexbin mapa <- list('sp.lines', world_sp) spplot(nasaSP["Ann"], col.regions=paleta3, sp.layout=mapa)
Now we can calculate the global effective irradiation with solaR (NOTE: this loop spends some time to finish. The results are available here.)
N=nrow(nasa) gefNasa <- matrix(nrow=N, ncol=3) for (i in 1:N){ prom=list(G0dm=as.numeric(nasa[i,3:14]*1000)) lat=nasa[i,1] gefFixed <- calcGef(lat=lat, prom=prom) gef2x <- calcGef(lat=lat, modeRad='prev', prev=gefFixed, modeTrk='two') gefHoriz <- calcGef(lat=lat, modeRad='prev', prev=gefFixed, modeTrk='horiz') gefNasa[i, ] <- sapply(list(gefFixed, gef2x, gefHoriz), function(x)as.data.frameY(x)$Gefd) print(i) } gefNasaDF <- as.data.frame(gefNasa) names(gefNasaDF) <- c('Fixed', 'Two', 'Horiz')
Ok, we got it. Let’s build an Spatial object with it:
gefNasaSP <- SpatialPixelsDataFrame(points=coords, data=gefNasaDF, proj4string=proj)
And now we can plot the results. First, the two-axis tracker:
ncuts <- 7 colContour <- paleta3(ncuts)[2] spplot(gefNasaSP['Two'], sp.layout=mapa, col.regions=paleta3, colorkey=list(space='bottom'), contour=TRUE, cuts=ncuts, col=colContour, lwd=0.5, scales=list(draw=TRUE) )
then, the N-S horizontal tracker:
spplot(gefNasaSP['Horiz'], sp.layout=mapa, col.regions=paleta3, colorkey=list(space='bottom'), contour=TRUE, cuts=ncuts, col=colContour, lwd=0.5, scales=list(draw=TRUE) )
and now the fixed generator:
spplot(gefNasaSP['Fixed'], sp.layout=mapa, col.regions=paleta3, colorkey=list(space='bottom'), contour=TRUE, cuts=ncuts, col=colContour, lwd=0.5, scales=list(draw=TRUE) )
Finally, let’s plot the irradiation incident on a NS horizontal tracker versus the irradiation incident on a fixed surface for the latitudes between -60º and 60º. I use the violin and box plots as described here.
gefNasaSP$HorizFixed <- gefNasaSP$Horiz/gefNasaSP$Fixed bwplot(HorizFixed~cut(Lat, pretty(Lat, 40)), xlab='Latitude', ylab=expression(G[ef]^{horiz}/G[ef]^{fixed}), data=as.data.frame(gefNasaSP), horizontal=FALSE, panel = function(..., box.ratio) { panel.violin(..., col = "lightblue", varwidth = FALSE, box.ratio = box.ratio) panel.bwplot(..., col='black', cex=0.8, pch='|', fill='gray', box.ratio = .1) }, par.settings = list(box.rectangle=list(col='black'), plot.symbol = list(pch='.', cex = 0.1)), scales=list(x=list(rot=45, cex=0.6)), subset=(abs(Lat)<60))
R-bloggers.com offers daily e-mail updates about R news and tutorials about learning R and many other topics. Click here if you're looking to post or find an R/data-science job.
Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.