Revision d6ddb3e2
Added by Adam Wilson almost 11 years ago
| climate/procedures/GetCloudProducts.R | ||
|---|---|---|
| 1 |
## download cloud products from GEWEX for comparison |
|
| 2 |
setwd("~/acrobates/adamw/projects/cloud/")
|
|
| 3 |
|
|
| 4 |
## get PATMOS-X data |
|
| 5 |
for(y in 2000:2009) system(paste("wget -nc -P data/gewex http://climserv.ipsl.polytechnique.fr/gewexca/DATA/instruments/PATMOSX/variables/CA/CA_PATMOSX_NOAA_0130PM_",y,".nc.gz",sep=""))
|
|
| climate/procedures/MOD09_Visualize.R | ||
|---|---|---|
| 15 | 15 |
|
| 16 | 16 |
#### Evaluate MOD35 Cloud data |
| 17 | 17 |
mc=brick("~/acrobates/adamw/projects/cloud/data/mod09.nc",varname="CF")
|
| 18 |
NAvalue(mod09)=-1
|
|
| 18 |
NAvalue(mc)=-1
|
|
| 19 | 19 |
|
| 20 | 20 |
cols=colorRampPalette(c("#000000","#00FF00","#FF0000"))#"black","blue","red"))
|
| 21 | 21 |
for(i in 1:156){
|
| ... | ... | |
| 60 | 60 |
|
| 61 | 61 |
|
| 62 | 62 |
pdf("output/mod09_worldclim.pdf",width=11,height=8.5)
|
| 63 |
p1=levelplot(wc_map,col.regions=grey(seq(0,1,len=100)),cuts=99,margin=F,maxpixels=1e5,colorkey=list("top"))
|
|
| 64 |
p2=levelplot(mac,col.regions=grey(seq(0,1,len=100)),cuts=99,margin=F,maxpixels=1e5,colorkey=list("bottom"))
|
|
| 65 |
p3=c(p1,p2,x.same=T,y.same=T,merge.legends=T) |
|
| 66 |
print(p3) |
|
| 67 |
|
|
| 68 |
|
|
| 69 | 63 |
regs=list( |
| 64 |
Cascades=extent(c(-122.8,-118,44.9,47)), |
|
| 70 | 65 |
Hawaii=extent(c(-156.5,-154,18.75,20.5)), |
| 71 | 66 |
Boliva=extent(c(-71,-63,-20,-15)), |
| 72 | 67 |
Venezuela=extent(c(-69,-59,0,7)), |
| 73 |
reg=extent(c(-83,-45,-5,13)), |
|
| 74 |
reg2=extent(c(-81,-70,-4,10)) |
|
| 68 |
CFR=extent(c(17.75,22.5,-34.8,-32.6)), |
|
| 69 |
Madagascar=extent(c(46,52,-17,-12)) |
|
| 70 |
#reg2=extent(c(-81,-70,-4,10)) |
|
| 75 | 71 |
) |
| 76 |
|
|
| 77 |
r=3 #set region |
|
| 72 |
for(r in 1:length(regs)){
|
|
| 78 | 73 |
p_map=levelplot(crop(wc_map,regs[[r]]),col.regions=grey(seq(0,1,len=100)),cuts=99,margin=F,maxpixels=1e5,colorkey=list(space="bottom",height=.75),xlab="",ylab="",main=names(regs)[r],useRaster=T) |
| 79 | 74 |
p_mac=levelplot(crop(mac,regs[[r]]),col.regions=grey(seq(0,1,len=100)),cuts=99,margin=F,maxpixels=1e5,colorkey=list(space="bottom",height=.75),useRaster=T) |
| 80 |
p_npp=levelplot(crop(npp,regs[[r]]),col.regions=grey(seq(0,1,len=100)),cuts=99,margin=F,maxpixels=1e5,colorkey=list(space="bottom",height=.5),zscaleLog=T,useRaster=T)
|
|
| 81 |
p3=c("WorldClim Mean Annual Precip (mm)"=p_map,"MOD09 Cloud Frequency (%)"=p_mac,"NPP"=p_npp,x.same=T,y.same=T,merge.legends=T,layout=c(3,1))
|
|
| 75 |
#p_npp=levelplot(crop(npp,regs[[r]]),col.regions=grey(seq(0,1,len=100)),cuts=99,margin=F,maxpixels=1e5,colorkey=list(space="bottom",height=.5),zscaleLog=T,useRaster=T) #"NPP"=p_npp,
|
|
| 76 |
p3=c("WorldClim Mean Annual Precip (mm)"=p_map,"MOD09 Cloud Frequency (%)"=p_mac,x.same=T,y.same=T,merge.legends=T,layout=c(2,1))
|
|
| 82 | 77 |
print(p3) |
| 83 |
|
|
| 78 |
} |
|
| 84 | 79 |
dev.off() |
| 85 | 80 |
|
| 86 | 81 |
|
| climate/procedures/NDP-026D.R | ||
|---|---|---|
| 62 | 62 |
## overlay the data with 32km diameter (16km radius) buffer |
| 63 | 63 |
## buffer size from Dybbroe, et al. (2005) doi:10.1175/JAM-2189.1. |
| 64 | 64 |
buf=16000 |
| 65 |
mod09sta=lapply(1:nlayers(mod09),function(l) {print(l); extract(mod09[[l]],st,buffer=buf,fun=mean,na.rm=T,df=T)[,2]})
|
|
| 66 |
mod09st=do.call(cbind.data.frame,mod09sta) |
|
| 65 |
#mod09sta=lapply(1:nlayers(mod09),function(l) {print(l); extract(mod09[[l]],st,buffer=buf,fun=mean,na.rm=T,df=T)[,2]})
|
|
| 66 |
bins=cut(1:nrow(st),100) |
|
| 67 |
mod09sta=lapply(levels(bins),function(lb) {
|
|
| 68 |
l=which(bins==lb) |
|
| 69 |
td=extract(mod09,st[l,],buffer=buf,fun=mean,na.rm=T,df=T) |
|
| 70 |
td$id=st$id[l] |
|
| 71 |
print(lb)#as.vector(c(l,td[,1:4]))) |
|
| 72 |
write.table(td,"valid.csv",append=T,col.names=F,quote=F,sep=",",row.names=F) |
|
| 73 |
td |
|
| 74 |
})#,mc.cores=3) |
|
| 75 |
|
|
| 76 |
#mod09sta=extract(mod09,st,buffer=buf,fun=mean,na.rm=T,df=T) |
|
| 77 |
mod09st=read.csv("valid.csv",header=F)[,-c(1,2)]
|
|
| 78 |
|
|
| 79 |
#mod09st=do.call(rbind.data.frame,mod09sta) |
|
| 67 | 80 |
#mod09st=mod09st[,!is.na(colnames(mod09st))] |
| 68 |
colnames(mod09st)=names(mod09)
|
|
| 69 |
mod09st$id=st$id |
|
| 81 |
colnames(mod09st)=c(names(mod09),"id")
|
|
| 82 |
#mod09st$id=st$id
|
|
| 70 | 83 |
mod09stl=melt(mod09st,id.vars="id") |
| 71 | 84 |
mod09stl[,c("year","month")]=do.call(rbind,strsplit(sub("X","",mod09stl$variable),"[.]"))[,1:2]
|
| 72 | 85 |
## add it to cld |
| ... | ... | |
| 85 | 98 |
#lulc=crop(lulc,mod09) |
| 86 | 99 |
Mode <- function(x) {
|
| 87 | 100 |
ux <- na.omit(unique(x)) |
| 88 |
ux[which.max(tabulate(match(x, ux)),na.rm=T)]
|
|
| 101 |
ux[which.max(tabulate(match(x, ux)))] |
|
| 89 | 102 |
} |
| 90 |
lulcst=extract(lulc,st,fun=Mode,na.rm=T,buffer=buf,df=T)
|
|
| 103 |
lulcst=extract(lulc,st,fun=Mode,buffer=buf,df=T) |
|
| 91 | 104 |
colnames(lulcst)=c("id","lulc")
|
| 92 | 105 |
## add it to cld |
| 93 | 106 |
cld$lulc=lulcst$lulc[match(cld$StaID,lulcst$id)] |
| 94 |
cld$lulc=factor(as.integer(cld$lulc),labels=IGBP$class)
|
|
| 107 |
#cld$lulc=factor(as.integer(cld$lulc),labels=IGBP$class[sort(unique(cld$lulc))])
|
|
| 95 | 108 |
|
| 96 | 109 |
## update cld column names |
| 97 | 110 |
colnames(cld)[grep("Amt",colnames(cld))]="cld"
|
Also available in: Unified diff
udpating validation to full globe