11 |
11 |
wd="~/acrobates/adamw/projects/cloud"
|
12 |
12 |
setwd(wd)
|
13 |
13 |
|
|
14 |
tempdir="tmp"
|
|
15 |
if(!file.exists(tempdir)) dir.create(tempdir)
|
14 |
16 |
|
15 |
17 |
## Get list of available files
|
16 |
18 |
df=data.frame(path=list.files("/mnt/data2/projects/cloud/mod09",pattern="*.tif$",full=T,recur=T),stringsAsFactors=F)
|
... | ... | |
29 |
31 |
#df=df[df$month==1,]
|
30 |
32 |
table(df$year,df$month)
|
31 |
33 |
|
|
34 |
writeLines(paste("Tiling options will produce",nrow(tiles),"tiles and ",nrow(jobs),"tile-months. Current todo list is ",length(todo)))
|
|
35 |
|
|
36 |
|
32 |
37 |
## drop some if not complete
|
33 |
38 |
#df=df[df$month%in%1:9&df$year%in%c(2001:2012),]
|
34 |
|
rerun=F # set to true to recalculate all dates even if file already exists
|
|
39 |
rerun=T # set to true to recalculate all dates even if file already exists
|
35 |
40 |
|
36 |
41 |
## Loop over existing months to build composite netcdf files
|
37 |
42 |
foreach(date=unique(df$date)) %dopar% {
|
... | ... | |
45 |
50 |
if(!rerun&file.exists(ncfile)) return(NA)
|
46 |
51 |
## merge regions to a new netcdf file
|
47 |
52 |
# system(paste("gdal_merge.py -o ",tffile," -init -32768 -n -32768.000 -ot Int16 ",paste(df$path[df$date==date],collapse=" ")))
|
48 |
|
system(paste("gdalbuildvrt -overwrite -srcnodata -32768 ",vrtfile," ",paste(df$path[df$date==date],collapse=" ")))
|
|
53 |
system(paste("gdalbuildvrt -overwrite -srcnodata -32768 -vrtnodata -32768 ",vrtfile," ",paste(df$path[df$date==date],collapse=" ")))
|
49 |
54 |
## Warp to WGS84 grid and convert to netcdf
|
50 |
|
ops="-t_srs 'EPSG:4326' -multi -r cubic -te -90 -90 0 90 -tr 0.008333333333333 -0.008333333333333"
|
51 |
|
ops="-t_srs 'EPSG:4326' -multi -r cubic -te -180 -90 180 90 -tr 0.008333333333333 -0.008333333333333"
|
52 |
|
|
53 |
|
system(paste("gdalwarp -overwrite ",ops," -srcnodata -32768 -dstnodata -32768 -of netCDF ",vrtfile," ",ncfile," -ot Int16"))
|
54 |
|
# system(paste("gdalwarp -overwrite ",ops," -srcnodata -32768 -dstnodata -32768 -of netCDF ",vrtfile," ",tffile," -ot Int16"))
|
|
55 |
##ops="-t_srs 'EPSG:4326' -multi -r cubic -te -180 0 180 10 -tr 0.008333333333333 -0.008333333333333 -wo SOURCE_EXTRA=50" #-wo SAMPLE_GRID=YES -wo SAMPLE_STEPS=100
|
|
56 |
ops="-t_srs 'EPSG:4326' -multi -r cubic -te -180 -90 180 90 -tr 0.008333333333333 -0.008333333333333 -wo SOURCE_EXTRA=50"
|
55 |
57 |
|
|
58 |
system(paste("gdalwarp -overwrite ",ops," -et 0 -srcnodata -32768 -dstnodata -32768 ",vrtfile," ",tffile," -ot Int16"))
|
|
59 |
system(paste("gdal_translate -of netCDF ",tffile," ",ncfile))
|
|
60 |
# system(paste("gdalwarp -overwrite ",ops," -srcnodata -32768 -dstnodata -32768 -of netCDF ",vrtfile," ",tffile," -ot Int16"))
|
|
61 |
file.remove(tffile)
|
56 |
62 |
setwd(wd)
|
57 |
63 |
system(paste("ncecat -O -u time ",ncfile," ",ncfile,sep=""))
|
58 |
64 |
## create temporary nc file with time information to append to MOD06 data
|
Finished second complete ee download of mod09 cloud mask using python tiler