其他
【数据分享】月度基础气候栅格数据转年度数据处理与成图代码
前几天从数读城事公众号看到一篇文章,关于气候数据,具体这个包含啥我就不再赘述了,大家去数读城事看哈,链接在下面:
【数据下载】获取全球19类生物气候数据集和每月基础气候数据集
我就讲讲我用它干了嘛。
目前全球年累积降水量数据已经发布于国家地球系统科学数据中心网站,免费下载,欢迎大家使用!点击文末阅读原文即可跳转到下载页面!
01
年累积降水量
—
02
R语言数据处理
—
降水量数据下载下来后解压如下:我下载的是prec_2000-2009和prec_2010-2018,这样解压出来后,就是19年×12个月=228期降水数据。如果用ArcGIS处理就会比较麻烦,得提前把每年的筛出来,然后使用模型构建器循环。前段时间我试了试R进行影像的处理,这次当然就用R啦。
library(raster)
library(parallel) # R语言并行计算
tifnames <- list.files(path = "E:/WebGeoData/ClimateTiff/wc2.1_2.5m_prec_2000-2009", pattern = '.tif$')
tifdir <- paste("E:/WebGeoData/ClimateTiff/wc2.1_2.5m_prec_2000-2009/",tifnames,sep="")
cl.cores <- detectCores() #并行计算开始
cl <- makeCluster(cl.cores)
Prec2001<- parLapply(cl, tifdir[grep(pattern = "2001", tifnames)], raster)
Prec2002<- parLapply(cl, tifdir[grep(pattern = "2002", tifnames)], raster)
Prec2003<- parLapply(cl, tifdir[grep(pattern = "2003", tifnames)], raster)
Prec2004<- parLapply(cl, tifdir[grep(pattern = "2004", tifnames)], raster)
Prec2005<- parLapply(cl, tifdir[grep(pattern = "2005", tifnames)], raster)
Prec2006<- parLapply(cl, tifdir[grep(pattern = "2006", tifnames)], raster)
Prec2007<- parLapply(cl, tifdir[grep(pattern = "2007", tifnames)], raster)
Prec2008<- parLapply(cl, tifdir[grep(pattern = "2008", tifnames)], raster)
Prec2009<- parLapply(cl, tifdir[grep(pattern = "2009", tifnames)], raster)
stopCluster(cl) #并行计算结束
PrecBr2001<- brick(Prec2001) #R语言栅格处理块
PrecSum2001<- calc(PrecBr2001,sum) #对栅格进行求和,计算累积降水
writeRaster(PrecSum2001, "PrecSum2001", format="GTiff") #输出栅格
PrecBr2002<- brick(Prec2002)
PrecSum2002<- calc(PrecBr2002,sum)
writeRaster(PrecSum2002, "PrecSum2002", format="GTiff")
PrecBr2003<- brick(Prec2003)
PrecSum2003<- calc(PrecBr2003,sum)
writeRaster(PrecSum2003, "PrecSum2003", format="GTiff")
PrecBr2004<- brick(Prec2004)
PrecSum2004<- calc(PrecBr2004,sum)
writeRaster(PrecSum2004, "PrecSum2004", format="GTiff")
PrecBr2005<- brick(Prec2005)
PrecSum2005<- calc(PrecBr2005,sum)
writeRaster(PrecSum2005, "PrecSum2005", format="GTiff")
PrecBr2006<- brick(Prec2006)
PrecSum2006<- calc(PrecBr2006,sum)
writeRaster(PrecSum2006, "PrecSum2006", format="GTiff")
PrecBr2007<- brick(Prec2007)
PrecSum2007<- calc(PrecBr2007,sum)
writeRaster(PrecSum2007, "PrecSum2007", format="GTiff")
PrecBr2008<- brick(Prec2008)
PrecSum2008<- calc(PrecBr2008,sum)
writeRaster(PrecSum2008, "PrecSum2008", format="GTiff")
PrecBr2009<- brick(Prec2009)
PrecSum2009<- calc(PrecBr2009,sum)
writeRaster(PrecSum2009, "PrecSum2009", format="GTiff")
上面由月度数据转换为年度数据后,接下来就是出图,出图使用ArcMap很简单,一个个替换TIFF文件即可,但是手动弄太麻烦,在这里我继续使用R来完成。library(raster)
library(rasterVis)
library(RColorBrewer)
library(lattice)
plottifnames <- list.files(path = "E:/WebGeoData/ClimateTiff/ClimateRaster", pattern = '.tif$')
plottifdir <- paste("E:/WebGeoData/ClimateTiff/ClimateRaster/",plottifnames,sep="")
allrasters <- lapply(plottifdir, raster)
#制图并输出Jpg
for (i in 1:length(allrasters)) {
i <- i #可以省略,我是废话
jpeg(filename = paste("plot",i, ".jpeg") ,width = 864, height = 548, units = "px")
mapTheme <- rasterTheme(region=rev(brewer.pal(8,"RdBu")))
plt <-levelplot(allrasters[[i]],
main=paste(substr(allrasters[[i]]@data@names, start = 7, stop = 10), "年最高气温分布图"),
par.settings=mapTheme, #调色板
colorkey=list(labels=list(cex=1, font=2, col="brown"),
height=1, width=1.4,
title="℃" ,vjust=2)) #设置图例标注颜色和标题
print(plt) #必须有,否则空白图,plt也不行,必须带print
dev.off()
print(i)
}
单幅成图效果如下:数据为1961-2018,58年的图,很快就完了,出图速度完爆ArcMap