From: john seers (IFR) <john.seers_at_bbsrc.ac.uk>

Date: Mon 20 Nov 2006 - 13:49:23 GMT

Date: Mon 20 Nov 2006 - 13:49:23 GMT

--- John Seers Institute of Food Research Norwich Research Park Colney Norwich NR4 7UA tel +44 (0)1603 251497 fax +44 (0)1603 507723 e-mail john.seers@bbsrc.ac.uk e-disclaimer at http://www.ifr.ac.uk/edisclaimer/ Web sites: www.ifr.ac.uk www.foodandhealthnetwork.com -----Original Message----- From: r-help-bounces@stat.math.ethz.ch [mailto:r-help-bounces@stat.math.ethz.ch] On Behalf Of Jenny Barnes Sent: 20 November 2006 13:03 To: r-help@stat.math.ethz.ch Subject: [R] problem with loop to put data into array with missing data forsome files Dear R-help community, My main goal of this message is to find a way of skipping a file of a month/year in a loop that does not exist (and making it's output into an data.out array would be NA) and moving onto the next year/month in the loop to carry on filling data.out with real precipitation data. The situation so far: I downloaded 50 years worth of GRIB data files from the NCEP data site http://nomad3.ncep.noaa.gov/pub/reanalysis-1/month/grb2d.gau/ I then created a loop in R to read each month of each of the 50 years worth of files and only extract the precipitation records using wgrib and grep as show in the code at the end of this message. I had to use grep to extract the PRATE records each time, as unfortunately with these grib files precipitation is not always the same record number each time. This worked fine until Feb 1972 (which is file /home/jmb/sst_precip/gribdata/197202.grib). This file did not download completely and as a result the file that saved to my computer did not include the record for precipitation - only the first 13 records. At this point in the loop it just stopped with error message: Error in seq.default(lat.beg, lat.end, length = ny) : length must be non-negative number Which is obviously because ny does not actually exist for this file as it cannot be extracted from source as this record does not exist for this file! My question is: Does anybody have any suggestions as to how I could code it so that my loop will allow this file to just output NA into the data.out array in data.out$data and maybe keep FALSE in data.out$date for the month&years where this has occured and then continue to loop around all the other months and years in succession? Things that don't work that I have tried: 1)I cannot just restart the loop after any anomalous years as this "resets" data.out and changes its dimensions depending on the number of years looped (nyrs) etc. 2)The first point of trouble is when I call upon ny in the lineReceived on Tue Nov 21 21:52:55 2006

>lats <- seq(lat.beg, lat.end, length=ny)

as ny is NA in this record - so I tried putting the line

>test <- dims[15]

>if(test != "NA"){

with the closing } at the end of the rest of the loop ...this gives the error message

>Error in if (test != "NA") { : missing value where TRUE/FALSE needed

3) My thought was that maybe I would need something like

>if (parameter == "PRATE" ){

then continue with loop; but this will obviously not work by itself and gives the error message

>Error in if (test != "NA") { : missing value where TRUE/FALSE needed

- would need to specify what to do if parameter is not =="PRATE". I encourage anybody with another idea to contact me with their suggestions! Thank you for taking the time to read my problem, Sincerely, Jennifer Barnes r-help@stat.math.ethz.ch Here is my code: #takes precip data for every year and month needed #from file saved on computer already #and sorts it into a list (array) #set some control parameters startyr <-1954 #note needs to be year before data wanted nyrs <- 50 months <- nyrs*12 k <- 0 library(chron) library(utils) data.out <- list(lats=seq(88.542, -88.542, length=94), lons=seq(0, 360-1.875, length=192), date=vector(length=months), data=array(NA, c(months, 94*192)) ) #the loop for(yr in 1:nyrs){ for(mon in 1:12){ i<- startyr+yr if(mon < 10) {j <- 0} if(mon < 10) {mon2 <- mon} if(mon == 10) {(j <- 1)} if(mon == 10) {mon2 <- 0} if(mon == 11) {(j <- 1)} if(mon == 11) {mon2 <- 1} if(mon == 12) {(j <- 1)} if(mon == 12) {mon2 <- 2} outfile <- paste("/home/jmb/sst_precip/gribdata/",i, j, mon2, ".grib", sep="") #decode the grib cmd <- paste("/usr/local/bin/wgrib ", outfile, "| grep \":PRATE:\" | /usr/local/bin/wgrib", outfile, "-i -V -text -o /tmp/wgrib.output") status <- system(cmd, intern=T) #extract metadata needed for working with julian dates: dt <- unlist(strsplit(status[1], " "))[3] dt.yr <- as.numeric(substr(dt, 1,4)) dt.mon <- as.numeric(substr(dt, 5,6)) dt.day <- as.numeric(substr(dt, 7,8)) dt.hr <- as.numeric(substr(dt, 9,10)) dt.jul <- julian(dt.mon, dt.day+(dt.hr/24), dt.yr) parameter <- unlist(strsplit(status[1], " "))[4] param <- status[2] dims <- unlist(strsplit(status[3], " ")) nx <- as.numeric(dims[13]) ny <- as.numeric(dims[15]) temp <- unlist(strsplit(status[5]," ")) lat.beg <- as.numeric(temp[6]) lat.end <- as.numeric(temp[8]) lats <- seq(lat.beg, lat.end, length=ny) temp <- unlist(strsplit(status[6]," ")) lon.beg <- as.numeric(temp[14]) lon.end <- as.numeric(temp[16]) if(lon.beg > lon.end){lon.beg <- (lon.beg -360)} #since sp prefers -180 to 180 encoding lons <- seq(lon.beg, lon.end, length=nx) lons[lons > 180] <- lons[lons > 180]-360 #since sp prefers -180 to 180 encoding #read flat file into R (header first) dims <- scan("/tmp/wgrib.output", nlines=1, quiet=T) nrec <- dims[1]*dims[2] indata <- paste("indata.", i, j, mon2, sep="") indata <- scan("/tmp/wgrib.output", 'numeric', nlines=nrec, skip=1, quiet=T) indata <- matrix(as.numeric(indata), nrow=dims[1], ncol=dims[2], byrow=F) #to get precip into mm/month, not mm/sec: if((dt.mon==02)&(dt.yrendnumber%%4==0)){days.mon <- 29} else {days.mon <-28} #to take into account leap-years!!! if(dt.mon==01){days.mon <- 31} if(dt.mon==03){days.mon <- 31} if(dt.mon==05){days.mon <- 31} if(dt.mon==07){days.mon <- 31} if(dt.mon==08){days.mon <- 31} if(dt.mon==10){days.mon <- 31} if(dt.mon==12){days.mon <- 31} if(dt.mon==04){days.mon <- 30} if(dt.mon==06){days.mon <- 30} if(dt.mon==09){days.mon <- 30} if(dt.mon==11){days.mon <- 30} sec.mon <- days.mon*24*60*60 newindata <- paste("newindata.",i, j, mon2, sep="") newindata <- sec.mon*indata k<-k+1 data.out$date[k] <- paste(i, j, mon2, sep="") data.out$data[k,] <- newindata } } ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Jennifer Barnes PhD student - long range drought prediction Climate Extremes Department of Space and Climate Physics University College London Holmbury St Mary, Dorking Surrey RH5 6NT Web: http://climate.mssl.ucl.ac.uk ______________________________________________ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. ______________________________________________ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.

Archive maintained by Robert King, hosted by
the discipline of
statistics at the
University of Newcastle,
Australia.

Archive generated by hypermail 2.1.8, at Tue 21 Nov 2006 - 11:30:21 GMT.

*
Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help.
Please read the posting
guide before posting to the list.
*