Project

General

Profile

« Previous | Next » 

Revision 8009a075

Added by Adam Wilson almost 11 years ago

Updating validation to use new cloud data

View differences:

climate/procedures/NDP-026D.R
6 6 
library(doMC)
7 7 
library(rasterVis)
8 8 
library(rgdal)
9 

  		




  
  9
  
    
library(reshape)


  		




  10
  10
  
    

  		




  11
  11
  
    

  		




  12
  12
  
    
## Data available here http://cdiac.ornl.gov/epubs/ndp/ndp026d/ndp026d.html


  		




  ......




  22
  22
  
    
colnames(st)=c("id","elev","lat","lon")


  		




  23
  23
  
    
write.csv(st,"stations.csv",row.names=F)


  		




  24
  24
  
    
coordinates(st)=c("lon","lat")


  		




  
  25
  
    
projection(st)="+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"


  		




  
  26
  
    

  		




  25
  27
  
    
## download data


  		




  26
  28
  
    
system("wget -N -nd ftp://cdiac.ornl.gov/pub/ndp026d/cat67_78/* -A '.tc.Z' -P data/")


  		




  27
  29
  
    

  		




  ......




  33
  35
  
    

  		




  34
  36
  
    
## use monthly timeseries


  		




  35
  37
  
    
cld=do.call(rbind.data.frame,mclapply(sprintf("%02d",1:12),function(m) {


  		




  36
  
  
    
  d=read.fwf(list.files("data",pattern=paste("MNYDC.",m,".tc",sep=""),full=T),skip=1,widths=f162)


  		




  
  38
  
    
  d=read.fwf(list.files("data",pattern=paste("MNYDC.",m,".tc$",sep=""),full=T),skip=1,widths=f162)


  		




  37
  39
  
    
  colnames(d)=c162


  		




  38
  40
  
    
  d$month=as.numeric(m)


  		




  39
  41
  
    
  print(m)


  		




  ......




  53
  55
  
    

  		




  54
  56
  
    
## add the MOD09 data to cld


  		




  55
  57
  
    
#### Evaluate MOD35 Cloud data


  		




  56
  
  
    
mod09=brick("~/acrobates/adamw/projects/cloud/data/mod09.nc")


  		




  
  58
  
    
mod09=brick("~/acrobates/adamw/projects/cloud/data/cloud_ymonmean.nc")


  		




  57
  59
  
    

  		




  58
  60
  
    
## overlay the data with 32km diameter (16km radius) buffer


  		




  59
  61
  
    
## buffer size from Dybbroe, et al. (2005) doi:10.1175/JAM-2189.1.


  		




  60
  62
  
    
buf=16000


  		




  61
  63
  
    
bins=cut(1:nrow(st),100)


  		




  
  64
  
    
if(file.exists("valid.csv")) file.remove("valid.csv")


  		




  62
  65
  
    
mod09sta=lapply(levels(bins),function(lb) {


  		




  63
  66
  
    
  l=which(bins==lb)


  		




  64
  67
  
    
  td=extract(mod09,st[l,],buffer=buf,fun=mean,na.rm=T,df=T)


  		












Also available in:  Unified diff