| 1 |
9a19743f
|
Adam M. Wilson
|
### Script to compile the monthly cloud data from earth engine into a netcdf file for further processing
|
| 2 |
|
|
|
| 3 |
|
|
setwd("~/acrobates/adamw/projects/cloud")
|
| 4 |
|
|
library(raster)
|
| 5 |
|
|
library(doMC)
|
| 6 |
|
|
|
| 7 |
|
|
registerDoMC(10)
|
| 8 |
|
|
|
| 9 |
|
|
tempdir="tmp"
|
| 10 |
|
|
if(!file.exists(tempdir)) dir.create(tempdir)
|
| 11 |
|
|
|
| 12 |
|
|
## Get list of available files
|
| 13 |
|
|
df=data.frame(path=list.files("data/mod09",pattern="*.tif$",full=T),stringsAsFactors=F)
|
| 14 |
|
|
df[,c("region","year","month")]=do.call(rbind,strsplit(basename(df$path),"_|[.]"))[,c(2,3,4)]
|
| 15 |
|
|
df$date=as.Date(paste(df$year,"_",df$month,"_15",sep=""),"%Y_%m_%d")
|
| 16 |
|
|
|
| 17 |
|
|
table(df$year,df$month)
|
| 18 |
|
|
## drop some if not complete
|
| 19 |
|
|
df=df[df$year<=2006,]
|
| 20 |
|
|
|
| 21 |
|
|
## Loop over existing months to build composite netcdf files
|
| 22 |
|
|
foreach(date=unique(df$date)) %dopar% {
|
| 23 |
|
|
## get date
|
| 24 |
|
|
print(date)
|
| 25 |
|
|
## Define output and check if it already exists
|
| 26 |
|
|
ncfile=paste(tempdir,"/mod09_",date,".nc",sep="")
|
| 27 |
|
|
if(file.exists(ncfile)) next
|
| 28 |
|
|
## merge regions to a new netcdf file
|
| 29 |
|
|
system(paste("gdal_merge.py -o ",ncfile," -of netCDF -ot Byte -n 0 ",paste(df$path[df$date==date],collapse=" ")))
|
| 30 |
|
|
system(paste("ncecat -O -u time ",ncfile," ",ncfile,sep=""))
|
| 31 |
|
|
## create temporary nc file with time information to append to MOD06 data
|
| 32 |
|
|
cat(paste("
|
| 33 |
|
|
netcdf time {
|
| 34 |
|
|
dimensions:
|
| 35 |
|
|
time = 1 ;
|
| 36 |
|
|
variables:
|
| 37 |
|
|
int time(time) ;
|
| 38 |
|
|
time:units = \"days since 2000-01-01 00:00:00\" ;
|
| 39 |
|
|
time:calendar = \"gregorian\";
|
| 40 |
|
|
time:long_name = \"time of observation\";
|
| 41 |
|
|
data:
|
| 42 |
|
|
time=",as.integer(date-as.Date("2000-01-01")),";
|
| 43 |
|
|
}"),file=paste(tempdir,"/",date,"_time.cdl",sep=""))
|
| 44 |
|
|
system(paste("ncgen -o ",tempdir,"/",date,"_time.nc ",tempdir,"/",date,"_time.cdl",sep=""))
|
| 45 |
|
|
system(paste("ncks -A ",tempdir,"/",date,"_time.nc ",ncfile,sep=""))
|
| 46 |
|
|
## add other attributes
|
| 47 |
|
|
system(paste("ncrename -v Band1,CF ",ncfile,sep=""))
|
| 48 |
|
|
system(paste("ncatted ",
|
| 49 |
|
|
" -a units,CF,o,c,\"Proportion Days Cloudy\" ",
|
| 50 |
|
|
" -a missing_value,CF,o,b,255 ",
|
| 51 |
|
|
" -a _FillValue,CF,d,, ",
|
| 52 |
|
|
" -a valid_range,CF,o,b,\"0,100\" ",
|
| 53 |
|
|
" -a long_name,CF,o,c,\"Proportion cloudy days (%)\" ",
|
| 54 |
|
|
ncfile,sep=""))
|
| 55 |
|
|
## add the fillvalue attribute back (without changing the actual values)
|
| 56 |
|
|
system(paste("ncatted -a _FillValue,CF,o,b,255 ",ncfile,sep=""))
|
| 57 |
|
|
|
| 58 |
|
|
if(as.numeric(system(paste("cdo -s ntime ",ncfile),intern=T))<1) {
|
| 59 |
|
|
print(paste(ncfile," has no time, deleting"))
|
| 60 |
|
|
file.remove(ntime)
|
| 61 |
|
|
}
|
| 62 |
|
|
print(paste(basename(ncfile)," Finished"))
|
| 63 |
|
|
|
| 64 |
|
|
}
|
| 65 |
|
|
|
| 66 |
|
|
### merge all the tiles to a single global composite
|
| 67 |
|
|
#system(paste("ncdump -h ",list.files(tempdir,pattern="mod09.*.nc$",full=T)[10]))
|
| 68 |
|
|
system(paste("cdo -O mergetime ",paste(list.files(tempdir,pattern="mod09.*.nc$",full=T),collapse=" ")," data/mod09.nc"))
|
| 69 |
|
|
|
| 70 |
|
|
|
| 71 |
|
|
### generate the monthly mean and sd
|
| 72 |
|
|
system(paste("cdo -O merge -ymonmean data/mod09.nc -chname,CF,CF_sd -ymonstd data/mod09.nc -chname,CF,CF_annual -timmean data/mod09.nc data/mod09_clim.nc"))
|
| 73 |
|
|
|
