# 2014-03-03 
# Function that loads 
# the monthly mean anomalies
# into the R-Studio environment
#
# this function has the following
# parameter arguments must be provided

# (1) Station ID (allowed values are: 'USW00014735','USW00094728')
# (2) variable   (allowed values are: 'tavg',tmin','tmax','prcp')
# (3) month      (allowed values are: 'Jan','Feb','Mar','Apr',May','Jun','Jul',
#                                     'Aug','Sep','Oct','Nov','Dec',
#                                  or 'Ann' (for annual averages))
# (4) start year      (must be in range of the data e.g. 1800-2013)
# (5) end year        (must be in range of the data and > start year)
#

# 2014-03-27 updated to work with state subdirectories 
# 
loadano2<-function(station="USW00014735",state="NY",variable="tavg",month='Jan',start=1950,end=2012) {
  # if no argument is given during call, the above
  # default values will be used 
  #(Station Albany, monthly mean temperature anomaly, in month January 1950-2012)
  file<-paste("data/",state,"/",station,"_",variable,"_mon_mean_ano.csv",sep='')
  print(paste("Load monthly mean anomalies from file ",file,sep=''))
  print(paste("State ID: ",state,sep=''))
  print(paste("Station ID: ",station,sep=''))
  print(paste("Climate variable: ",sep=''))
  print(paste("selected month: ",month,sep=''))
  print(paste("selected time range: ",start,"-",end,sep=''))
  checkmon<-c('Jan','Feb','Mar','Apr','May','Jun','Jul',
              'Aug','Sep','Oct','Nov','Dec','Ann')
  mhelp<-1:12
  imon<-mhelp[checkmon==month]
  
  if(!any(month==checkmon)){
    print(paste("Error: wrong selection of month!",
                "Please check month =",month," is case sensitve!",sep=''))
  }
  b<-read.csv(file)
  # first case extract one month
  if (!month=='Ann') {
    imonths<-length(b$time)
    iyears<-imonths/12
    mhelp<-rep(1:12,iyears)
    selmon<-mhelp==imon
    print(imon)
    time<-b$time[selmon]
    ano<-b$ano[selmon]
    # subsample year range
    # to get Decembers I increase the end year value
    ihelp<-(start<=time & time <= end+0.99999)
  } else {
    print("calculate annual mean anomalies")
    # reshape the vector into a matrix
    imonths<-length(b$time)
    iyears<-imonths/12
    mhelp<-rep(1:12,iyears)
    b1<-matrix(b$time,12,iyears)
    test<-matrix(mhelp,12,iyears)
    b2<-matrix(b$ano,12,iyears)
    # now I can use one of R's elegant functions
    # to repeat the mean calculation over all months
    # in each year apply()
    ano<-apply(b2,2,mean,na.rm=TRUE)
    time<-apply(b1,2,mean,na.rm=TRUE)
    # subsample year range
    # to get Decembers I increase the end year value
    ihelp<-(start<=time & time <= end+0.99999)
  }
  
  
  # return only the selected years
  # strip decimals from time 
  res<-list(time=floor(time[ihelp]),ano=ano[ihelp])
  
  return(res)
}