Project

General

Profile

« Previous | Next » 

Revision 07f76f37

Added by Adam Wilson almost 12 years ago

MOD06 processing running successfully using Pleiades environment set up by Andrew Michaelis at NASA.

View differences:

climate/procedures/MOD06_L2_process.r
1 1
###################################################################################
2 2
###  R code to aquire and process MOD06_L2 cloud data from the MODIS platform
3 3

  
4
## load command line arguments (mname)
5
args=(commandArgs(TRUE)) ##args is now a list of character vectors
4
## import commandline arguments
5
library(getopt)
6
## get options
7
opta <- getopt(matrix(c(
8
                        'date', 'd', 1, 'character',
9
                        'tile', 't', 1, 'character',
10
                        'help', 'h', 0, 'logical'
11
                        ), ncol=4, byrow=TRUE))
12
if ( !is.null(opta$help) )
13
  {
14
       prg <- commandArgs()[1];
15
          cat(paste("Usage: ", prg,  " --date | -d <file> :: The date to process\n", sep=""));
16
          q(status=1);
17
     }
18

  
19
#gridfile = paste("HDF4_EOS:EOS_GRID:\"", opta$data,"\":MODIS_Grid_Daily_1km_LST:Night_view_angl", sep="")
20
#cat(paste("Checking : ", gridfile, " ...\n"));
21

  
22

  
23
## Load command line arguments (mname)
24
#args=(commandArgs(TRUE)) ##args is now a list of character vectors
6 25
## Then cycle through each element of the list and evaluate the expressions.
7
eval(parse(text=args))
26
#eval(parse(text=args))
27
#print(args)
8 28

  
9 29
#system("module list")
10 30
#system("source ~/moduleload")
11 31
#system("module list")
12 32

  
13
print(args)
14 33

  
15
tile="h11v08"
34
date=opta$date  #date="20030301"
35
tile=opta$tile #tile="h11v08"
16 36
outdir="2_daily"  #directory for separate daily files
17
outdir2="3_summary" #directory for combined daily files and summarized files
37
#outdir2="3_summary" #directory for combined daily files and summarized files
18 38

  
19 39
print(paste("Processing tile",tile," for date",date))
20 40

  
......
133 153
  print(paste("Set up temporary grass session in",tf))
134 154

  
135 155
  ## set up temporary grass instance for this PID
136
  initGRASS(gisBase="/nobackupp1/awilso10/software/grass-6.4.3svn",gisDbase=tf,SG=td,override=T,location="mod06",mapset="PERMANENT",home=tf,pid=Sys.getpid())
156
#  initGRASS(gisBase="/nobackupp1/awilso10/software/grass-6.4.3svn",gisDbase=tf,SG=td,override=T,location="mod06",mapset="PERMANENT",home=tf,pid=Sys.getpid())
157
  initGRASS(gisBase="/u/armichae/pr/grass-6.4.2/",gisDbase=tf,SG=td,override=T,location="mod06",mapset="PERMANENT",home=tf,pid=Sys.getpid())
158

  
137 159
  system(paste("g.proj -c proj4=\"",projection(td),"\"",sep=""))
138 160

  
139 161
  ## Define region by importing one MOD11A1 raster.
......
231 253
    execGRASS("i.group",group="mod06",input=c("CER_daily","COT_daily","CLD_daily")) ; print("")
232 254
   
233 255
  ncfile=paste(outdir,"/MOD06_",date,".nc",sep="")
234
  execGRASS("r.out.gdal",input="mod06",output=ncfile,type="Int16",nodata=-32768,flags=c("quiet"),createopt=c("WRITE_GDAL_TAGS=YES","WRITE_LONLAT=NO"),format="netCDF")
256
  file.remove(ncfile)
257
  execGRASS("r.out.gdal",input="mod06",output=ncfile,type="Int16",nodata=-32768,flags=c("quiet"),
258
#      createopt=c("FORMAT=NC4","ZLEVEL=5","COMPRESS=DEFLATE","WRITE_GDAL_TAGS=YES","WRITE_LONLAT=NO"),format="netCDF")
259
      createopt=c("FORMAT=NC","WRITE_GDAL_TAGS=YES","WRITE_LONLAT=NO"),format="netCDF")
260

  
235 261
  ncopath="/nasa/sles11/nco/4.0.8/gcc/mpt/bin/"
236 262
  system(paste(ncopath,"ncecat -O -u time ",ncfile," ",ncfile,sep=""))
237
  system(paste(ncopath,"ncap2 -O -s 'time[time]=",as.integer(fs$date[fs$dateid==date]-as.Date("2000-01-01")),"'",ncfile," ",ncfile,sep=""))
238
  system(paste(ncopath,"ncatted -a calendar,time,c,c,\"standard\" -a long_name,time,c,c,\"time\" -a units,time,c,c,\"days since 2000-01-01 12:00:00\" ",ncfile,sep=""))
263
#  system(paste("cdo setdate,",as.integer(fs$date[fs$dateid==date]-as.Date("2000-01-01"))," -settaxis,2000-01-01,12:00:00,1day ",ncfile," ",ncfile,sep=""))
264
#  timdimf=paste(tempdir(),"/time.nc",sep="")
265
#  nc_dtime=ncdim_def("time","days since 2000-01-01 12:00:00",
266
#                                  vals=as.integer(fs$date[fs$dateid==date]-as.Date("2000-01-01")),
267
#                                  calendar="gregorian",longname="time")
268
#  nc_vtime=ncvar_def("time","
269
#  nc=nc_create(timdimf,nc_time)
270
## create temporary nc file with time information to append to MOD06 data
271
  cat(paste("
272
    netcdf time {
273
      dimensions:
274
        time = 1 ;
275
      variables:
276
        int time(time) ;
277
      time:units = \"days since 2000-01-01 12:00:00\" ;
278
      time:calendar = \"gregorian\";
279
      time:long_name = \"time of observation\"; 
280
    data:
281
      time=",as.integer(fs$date[fs$dateid==date][1]-as.Date("2000-01-01")),";
282
    }"),file=paste(tempdir(),"/time.cdl",sep=""))
283
system(paste("ncgen -o ",tempdir(),"/time.nc ",tempdir(),"/time.cdl",sep=""))
284
system(paste(ncopath,"ncks -A ",tempdir(),"/time.nc ",ncfile,sep=""))
285

  
286
# ncwa -O -h -a record out.nc out.nc
287
#  ncvar_put( nc, "time", vals=as.integer(fs$date[fs$dateid==date]-as.Date("2000-01-01")), start=NA, count=NA, verbose=FALSE )
288
#  system(paste(ncopath,"ncap2 -O -s 'time[time]=",as.integer(fs$date[fs$dateid==date]-as.Date("2000-01-01")),"'",ncfile," ",ncfile,sep=""))
289
#  system(paste(ncopath,"ncatted -O calendar,time,c,c,\"standard\" -a long_name,time,c,c,\"time\" -a units,time,c,c,\"days since 2000-01-01 12:00:00\" ",ncfile,sep=""))
239 290
  system(paste(ncopath,"ncrename -v Band1,CER -v Band2,COT -v Band3,CLD ",ncfile,sep=""))
240 291
  system(paste(ncopath,"ncatted -a scale_factor,CER,o,d,0.01 -a units,CER,o,c,\"micron\" -a missing_value,CER,o,d,-32768 -a long_name,CER,o,c,\"Cloud Particle Effective Radius\" ",ncfile,sep=""))
241 292
  system(paste(ncopath,"ncatted -a scale_factor,COT,o,d,0.01 -a units,COT,o,c,\"none\" -a missing_value,COT,o,d,-32768 -a long_name,COT,o,c,\"Cloud Optical Thickness\" ",ncfile,sep=""))
......
244 295
  
245 296
### delete the temporary files 
246 297
  unlink_.gislock()
247
  system("/nobackupp1/awilso10/software/grass-6.4.3svn/etc/clean_temp")
298
#  system("/nobackupp1/awilso10/software/grass-6.4.3svn/etc/clean_temp")
248 299
  system(paste("rm -rR ",tf,sep=""))
249 300
}
250 301

  
......
283 334
#foreach(i=notdone[1:3],.packages=(.packages())) %dopar% mod06(i,tile)
284 335
#foreach(i=1:20) %dopar% print(i)
285 336

  
286

  
287
q("no")
337
q("no",status=0)
climate/procedures/Pleiades.R
42 42
table(done)
43 43
notdone=alldates[!done]  #these are the dates that still need to be processed
44 44

  
45
if(exists("fdly")){
46
  notdone=notdone[notdone%in%fdly$dateid[is.na(fdly$npar)]]
47

  
48

  
45 49
tile="h11v08"   #can move this to submit script if needed
46 50
script="/u/awilso10/environmental-layers/climate/procedures/MOD06_L2_process.r"
47
#write.table(paste("--verbose ",script," date=",notdone," tile=\"",tile,"\"",sep=""),file="notdone.txt",row.names=F,col.names=F,quote=F)
51
write.table(paste("--verbose ",script," --date ",notdone," --tile \"",tile,"\"",sep=""),file="notdone.txt",row.names=F,col.names=F,quote=F)
48 52
#write.table(paste("--verbose ",script," date=",notdone[1:30],sep=""),file="notdone.txt",row.names=F,col.names=F,quote=F)
49
write.table(notdone[1:30],file="notdone.txt",row.names=F,col.names=F,quote=F)
53
#write.table(notdone[1:30],file="notdone.txt",row.names=F,col.names=F,quote=F)
50 54

  
51 55
save(fs,alldates,gridfile,td,file="allfiles.Rdata")
52 56

  
53 57
## run script
54
cat(paste("
55
#! /bin/bash
56
source ~/moduleload
57
source ~/.bashrc
58
Rscript --verbose --vanilla /u/awilso10/environmental-layers/climate/procedures/MOD06_L2_process.r date=$1
58
#cat(paste("
59
##! /bin/bash
60
#source ~/moduleload
61
#source ~/.bashrc
62
#Rscript --verbose --vanilla /u/awilso10/environmental-layers/climate/procedures/MOD06_L2_process.r date=$1
59 63
#Rscript --verbose --vanilla rtest
60
",sep=""),file="MOD06_process2")
61
system("chmod +x MOD06_process2")
64
#",sep=""),file="MOD06_process2")
65
#system("chmod +x MOD06_process2")
62 66

  
63
cat(paste("
64
library(rgdal)
65
GDALinfo
66
",sep=""),file="rtest")
67
#cat(paste("
68
#library(rgdal)
69
#GDALinfo
70
#",sep=""),file="rtest")
67 71

  
68 72

  
69 73
## Submission script
70 74

  
71
cat(paste("
72
#PBS -S /bin/bash
73
#PBS -l select=2:ncpus=16:model=san
74
###PBS -l select=4:ncpus=8:model=neh
75
##PBS -l select=1:ncpus=12:model=wes
76
####### old: select=48:ncpus=8:mpiprocs=8:model=neh
77
#PBS -l walltime=2:00:00
78
#PBS -j oe
79
#PBS -m e
80
#PBS -V
81
#PBS -q devel
82
#PBS -o log/log_^array_index^
83
#PBS -o log/log_DataCompile.log
84
#PBS -M adam.wilson@yale.edu
85
#PBS -N MOD06
75
#cat(paste("
76
##PBS -S /bin/bash
77
##PBS -l select=2:ncpus=16:model=san
78
####PBS -l select=4:ncpus=8:model=neh
79
###PBS -l select=1:ncpus=12:model=wes
80
######## old: select=48:ncpus=8:mpiprocs=8:model=neh
81
##PBS -l walltime=2:00:00
82
##PBS -j oe
83
##PBS -m e
84
##PBS -V
85
##PBS -q devel
86
##PBS -o log/log_^array_index^
87
##PBS -o log/log_DataCompile.log
88
##PBS -M adam.wilson@yale.edu
89
##PBS -N MOD06
86 90

  
87 91
## cd to working directory
88
cd /nobackupp1/awilso10/mod06
92
#cd /nobackupp1/awilso10/mod06
89 93

  
90 94
## set some memory limits
91 95
#  ulimit -d 1500000 -m 1500000 -v 1500000  #limit memory usage
92
  source /usr/local/lib/global.profile
93
  source /u/awilso10/.bashrc
94
  source /u/awilso10/moduleload
95
## export a few important variables
96
  export NNODES=32
97
  export R_LIBS=\"/u/awilso10/R/x86_64-unknown-linux-gnu-library/2.15/\"
98
## Run the script!
96
#  source /usr/local/lib/global.profile
97
#  source /u/awilso10/.bashrc
98
#  source /u/awilso10/moduleload
99
### export a few important variables
100
#  export NNODES=32
101
#  export R_LIBS=\"/u/awilso10/R/x86_64-unknown-linux-gnu-library/2.15/\"
102
### Run the script!
99 103
## current version not parallelizing across nodes!
100 104
#  TMPDIR=$TMPDIR Rscript --verbose --vanilla /u/awilso10/environmental-layers/climate/procedures/MOD06_L2_process.r date=20000403
101 105

  
102
WORKLIST=notdone.txt
103
#EXE=\"Rscript\"
104
EXE="./MOD06_process2"
105
LOG=log/log_DataCompile.log
106
MQUEUE=/nobackupp4/pvotava/software/share/mqueue-eg/mqueue/mqueue
106
#WORKLIST=notdone.txt
107
##EXE=\"Rscript\"
108
#EXE="./MOD06_process2"
109
#LOG=log/log_DataCompile.log
110
#MQUEUE=/nobackupp4/pvotava/software/share/mqueue-eg/mqueue/mqueue
111

  
112
#TMPDIR=$TMPDIR mpiexec -np $NNODES $MQUEUE -l $WORKLIST -p $EXE -v -v -v --random-starts 2-4 --work-analyze #> $LOG
113
#exit 0
114
#",sep=""),file="MOD06_process")
115

  
116

  
117
### simplified qsub script
118
cat(paste("
119
#PBS -S /bin/bash
120
#PBS -l select=48:ncpus=8:mpiprocs=8
121
##PBS -l select=2:ncpus=4:mpiprocs=4
122
#PBS -l walltime=02:00:00
123
#PBS -j n
124
#PBS -m be
125
#PBS -N mod06
126
#PBS -q devel
127
#PBS -V
128

  
129
CORES=384
130
HDIR=/u/armichae/pr/
131
  source $HDIR/etc/environ.sh
132
  source /u/awilso10/.bashrc
133
IDIR=/nobackupp1/awilso10/mod06/
134
##WORKLIST=$HDIR/var/run/pxrRgrs/work.txt
135
WORKLIST=$IDIR/notdone.txt
136
EXE=Rscript
137
LOGSTDOUT=$IDIR/log/log.stdout
138
LOGSTDERR=$IDIR/log/log.stderr
139
mpiexec -np $CORES pxargs -a $WORKLIST -p $EXE -v -v -v --work-analyze 1> $LOGSTDOUT 2> $LOGSTDERR
140
"),file="mod06_qsub")
107 141

  
108
TMPDIR=$TMPDIR mpiexec -np $NNODES $MQUEUE -l $WORKLIST -p $EXE -v -v -v --random-starts 2-4 --work-analyze #> $LOG
109
exit 0
110
",sep=""),file="MOD06_process")
111 142

  
112 143
### Check the file
113
system("cat MOD06_process")
144
system("cat mod06_qsub")
114 145
#system("cat ~/environmental-layers/climate/procedures/MOD06_L2_process.r")
115 146

  
116 147
## check queue status
......
118 149
system("/u/scicon/tools/bin/qtop.pl 492352")
119 150

  
120 151
## Submit it (and keep the pid)!
121
system("qsub MOD06_process")
122
system("/u/scicon/tools/bin/pdsh_gdb -j 493281 -d tmp -s -u awilso10")
152
system("qsub mod06_qsub")
153
system("/u/scicon/tools/bin/pdsh_gdb -j 568835 -d tmp -s -u awilso10")
123 154

  
124 155
## work in interactive mode
125 156
# system("qsub -I -l walltime=2:00:00 -lselect=2:ncpus=16:model=san -q devel")
......
127 158

  
128 159
## check progress
129 160
system("qstat -u awilso10")
130
system(paste("/u/scicon/tools/bin/qps ",pid))
161
system(paste("/u/scicon/tools/bin/qps ",568835))
131 162
system(paste("qstat -t -x",pid))
132 163

  
133 164
system("qstat devel ") 
......
141 172
fdly=data.frame(
142 173
  path=list.files(outdir,pattern="nc$",full=T),
143 174
  file=list.files(outdir,pattern="nc$"))
144
fdly$date=as.Date(substr(fdly$file,7,14),"%Y%m%d")
145
fdly$month=format(fdly$date,"%m")
146
fdly$year=format(fdly$date,"%Y")
175
  fdly$dateid=substr(fdly$file,7,14)
176
  fdly$date=as.Date(substr(fdly$file,7,14),"%Y%m%d")
177
  fdly$month=format(fdly$date,"%m")
178
  fdly$year=format(fdly$date,"%Y")
147 179

  
148 180
## check validity (via npar and ntime) of nc files
149 181
for(i in 1:nrow(fdly)){
150
  fdly$ntime[i]=as.numeric(system(paste("cdo  sinfo ",fdly$path[i]),intern=T))
151
  fdly$npar[i]=as.numeric(system(paste("cdo -s npar ",fdly$path[i]),intern=T))
182
 # fdly$ntime[i]<-as.numeric(system(paste("cdo  sinfo ",fdly$path[i]),intern=T))
183
  fdly$npar[i]<-as.numeric(system(paste("cdo -s npar ",fdly$path[i]),intern=T))
152 184
  print(i)
153 185
}
154 186

  
......
156 188
tsdir=paste(tempdir(),"/summary",sep="")
157 189
dir.create(tsdir)
158 190
lapply(unique(fdly$year),function(y){
159
  system(paste("cdo -O mergetime ",paste(fdly$path[fdly$year==y],collapse=" ")," ",tsdir,"/MOD09_",tile,"_",y,"_daily.nc",sep=""))
191
  system(paste("cdo -O mergetime ",paste(fdly$path[!is.na(fdly$npar)&fdly$year==y],collapse=" ")," ",tsdir,"/MOD09_",tile,"_",y,"_daily.nc",sep=""))
160 192
  print(paste("Finished merging daily files for year",y))
161 193
})
162 194
## Combine the year-by-year files into a single daily file
......
168 200

  
169 201
print("Finished!   %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
170 202
## quit R
171
q("no")
203
#q("no")
172 204
 
173 205

  
174 206
#################################################################
175 207

  
176 208
### copy the files back to Yale
177 209
list.files("2_daily")
178
system("scp 2_daily/* adamw@acrobates.eeb.yale.edu:/data/personal/adamw/projects/interp/data/modis/Venezuela")
210
system(paste("scp ",outdir2,"/*_ymonmean.nc adamw@acrobates.eeb.yale.edu:/data/personal/adamw/projects/interp/data/modis/Venezuela/summary",sep=""))
179 211

  
180 212
system("scp  /tmp/Rtmp6I6tFn/MOD06_L2.A2000061.1615.051.2010273184629.hdf adamw@acrobates.eeb.yale.edu:/data/personal/adamw/projects/interp/data/modis/Venezuela")
181 213
system("scp 2_daily/MOD06_20000410.nc adamw@acrobates.eeb.yale.edu:/data/personal/adamw/projects/interp/data/modis/Venezuela")
182 214

  
183 215

  
184
list.files(" /tmp/Rtmp6I6tFn")
185 216

  
186 217

  
187 218

  

Also available in: Unified diff