R/wt.phase.image.R

In WaveletComp: Computational Wavelet Analysis

wt.phase.image <-
function(WT, my.series = 1,
       plot.coi = TRUE, plot.contour = TRUE, siglvl = 0.1, col.contour = "white",
       plot.ridge = TRUE, col.ridge = "black", 
       n.levels = 100, color.palette = "rainbow(n.levels, start = 0, end = .7)", 
       useRaster = TRUE, max.contour.segments = 250000,
       plot.legend = TRUE,
       legend.params = list(width = 1.2, shrink = 0.9, mar = 5.1, 
                            n.ticks = 6, pi.style = TRUE,
                            label.digits = 1, label.format = "f", 
                            lab = NULL, lab.line = 3),
       label.time.axis = TRUE, 
       show.date = FALSE, date.format = NULL, date.tz = NULL, 
       timelab = NULL, timetck = 0.02, timetcl = 0.5,
       spec.time.axis = list(at = NULL, labels = TRUE, 
                             las = 1, hadj = NA, padj = NA),
       label.period.axis = TRUE, 
       periodlab = NULL, periodtck = 0.02, periodtcl = 0.5,
       spec.period.axis = list(at = NULL, labels = TRUE, 
                               las = 1, hadj = NA, padj = NA),
       main = NULL,
       lwd = 2, lwd.axis = 1, 
       graphics.reset = TRUE,
       verbose = FALSE) {       
         
  ##########################################       
 
  if(verbose == T){
       out <- function(...){ cat(...) }
    }
    else{
       out <- function(...) { }
    }  
    
  default.options = options()   
 
  options(max.contour.segments = as.integer(max.contour.segments))
  
  
  series.data = WT$series
    
  ####################################
  ## Identify the scenario
  ####################################
    
  if (class(WT) == 'analyze.wavelet') {
    
      out("Your input object class is 'analyze.wavelet'...\n") 
      
      my.series = ifelse(names(series.data)[1] == 'date', names(series.data)[2], names(series.data)[1]) 
        
      Phase = WT$Phase  
      Power = WT$Power
      Power.pval = WT$Power.pval
      Ridge = WT$Ridge     
        
  } 
  if (class(WT) == 'analyze.coherency') {   
    
      out("Your input object class is 'analyze.coherency'...\n") 
            
      if (is.numeric(my.series)) { 
          if (!is.element(my.series,c(1,2))) { stop("Please choose either series number 1 or 2!") }
          my.series = ifelse(names(series.data)[1] == 'date', names(series.data)[my.series+1], names(series.data)[my.series])  
      }
        
      ind = which( names(series.data) == my.series ) 
      which.series.num = ifelse(names(series.data)[1] == 'date', ind-1, ind)
      if (!is.element(which.series.num, c(1,2))) { stop("Your series name is not available, please check!") }
      
      if (which.series.num == 1) {
          Phase = WT$Phase.x
          Power = WT$Power.x
          Power.pval = WT$Power.x.pval
          Ridge = WT$Ridge.x
      }
      if (which.series.num == 2) {
          Phase = WT$Phase.y
          Power = WT$Power.y
          Power.pval = WT$Power.y.pval
          Ridge = WT$Ridge.y
     }      
        
  }   
           
  out(paste("A phase image of your time series '", my.series, "' will be plotted...", sep=''), '\n')  
  
    
  ##########################################
                          
  axis.1 <- WT$axis.1
  axis.2 <- WT$axis.2 
  
  phase.levels = seq(-pi,pi,length.out=n.levels+1)
  key.cols = rev(eval(parse(text=color.palette)))
  
  # legend parameters  
  
  if (is.null(legend.params$width))         legend.params$width = 1.2
  if (is.null(legend.params$shrink))        legend.params$shrink = 0.9
  if (is.null(legend.params$mar)) legend.params$mar = ifelse(is.null(legend.params$lab), 5.1, 6.1)   
  if (is.null(legend.params$n.ticks))       legend.params$n.ticks = 6
  if (is.null(legend.params$pi.style))      legend.params$pi.style = TRUE
  if (is.null(legend.params$label.digits)) {
      if (legend.params$pi.style == TRUE) {legend.params$label.digits = 1}
      else {legend.params$label.digits = 2}
  }    
  if (is.null(legend.params$label.format))  legend.params$label.format = "f"
  if (is.null(legend.params$lab.line))      legend.params$lab.line = 3
  
  # date parameters
  
  if (is.null(date.format)) { date.format = WT$date.format }
  if (is.null(date.tz)) { date.tz = ifelse(is.null(WT$date.tz),"",WT$date.tz) }
  
  # individual time axis parameters
  
  if (!is.list(spec.time.axis))       spec.time.axis = list()
  if (is.null(spec.time.axis$at))     spec.time.axis$at = NULL
  if (is.null(spec.time.axis$labels)) spec.time.axis$labels = T
  if (is.null(spec.time.axis$las))    spec.time.axis$las = 1
  if (is.null(spec.time.axis$hadj))   spec.time.axis$hadj = NA
  if (is.null(spec.time.axis$padj))   spec.time.axis$padj = NA 
  
  # initialize warning indicators for the case of individual time axis specification
  # warning: reset to time axis default?
  time.axis.warning = F 
  # warning: NAs among time axis tick marks?
  time.axis.warning.na = F
  # warning: calendar dates not chronological and/or format not standard unambiguous 
  chronology.warning = F
    
  if ( (!is.null(spec.time.axis$at)) & (label.time.axis==F) )  warning("\nPlease set label.time.axis = TRUE to make time axis specification effective.", immediate. = TRUE)
  
  # individual period axis parameters
  
  if (!is.list(spec.period.axis))       spec.period.axis = list()
  if (is.null(spec.period.axis$at))     spec.period.axis$at = NULL
  if (is.null(spec.period.axis$labels)) spec.period.axis$labels = T
  if (is.null(spec.period.axis$las))    spec.period.axis$las = 1
  if (is.null(spec.period.axis$hadj))   spec.period.axis$hadj = NA
  if (is.null(spec.period.axis$padj))   spec.period.axis$padj = NA
  
  # initialize warning indicators for the case of individual period axis specification
  # warning: reset to period axis default?
  period.axis.warning = F  
  # warning: NAs among period axis tick marks?
  period.axis.warning.na = F
  
  if ( (!is.null(spec.period.axis$at)) & (label.period.axis==F) )  warning("\nPlease set label.period.axis = TRUE to make period axis specification effective.", immediate. = TRUE)
  
  #######################################################################################
  ## start plotting
  #######################################################################################
  
  op = par(no.readonly = TRUE)
    
  image.plt  = par()$plt
  legend.plt = NULL
  
  if (plot.legend == T) {
  
      # construct plot regions for image and legend
   
      legend.plt = par()$plt
      
      char.size = par()$cin[1]/par()$din[1]
    
      hoffset       = char.size * par()$mar[4]
      legend.width  = char.size * legend.params$width
      legend.mar    = char.size * legend.params$mar
      
      legend.plt[2] = 1 - legend.mar
      legend.plt[1] = legend.plt[2] - legend.width
  
      vmar = (legend.plt[4] - legend.plt[3]) * ((1 - legend.params$shrink)/2)
  
      legend.plt[4] = legend.plt[4] - vmar
      legend.plt[3] = legend.plt[3] + vmar

      image.plt[2] = min(image.plt[2], legend.plt[1] - hoffset)
      
      # plot legend first
               
      par(plt = legend.plt)

      key.marks  = round(seq(from = 0, to = 1, length.out=legend.params$n.ticks)*n.levels)
      if (legend.params$pi.style == FALSE) {
          line.add = 0
          key.labels = formatC(as.numeric(phase.levels), digits = legend.params$label.digits, format = legend.params$label.format)[key.marks+1]
      }
      if (legend.params$pi.style == TRUE) {
          line.add = 0.6
          pi.fac.seq = seq(-1,1,length.out=n.levels+1)
          neg.ind = (pi.fac.seq<0)[key.marks+1] 
          pos.ind = (pi.fac.seq>0)[key.marks+1] 
          pi.fac = formatC(abs(pi.fac.seq), digits = legend.params$label.digits, format = legend.params$label.format)[key.marks+1] 
          inner.key.labels = NULL
          for (label.i in (2:(legend.params$n.ticks-1))) { 
               pi.fac.i = pi.fac[label.i] 
               if (neg.ind[label.i] == TRUE) { inner.key.labels = c( inner.key.labels, bquote(-.(pi.fac.i)*pi)) }
               if ((neg.ind[label.i] == FALSE) & (pos.ind[label.i] == FALSE)) { inner.key.labels = c( inner.key.labels, 0) }
               if ((neg.ind[label.i] == FALSE) & (pos.ind[label.i] == TRUE)) { inner.key.labels = c( inner.key.labels, bquote(.(pi.fac.i)*pi)) }
               }    
          key.labels = c(expression(-pi), inner.key.labels, expression(pi))
      }
      
#       key.phase.levels = quantile(phase.levels, seq(from = 0, to = 1, length.out=legend.params$n.ticks))
#       key.labels = formatC(as.numeric(key.phase.levels), digits = legend.params$label.digits, format = legend.params$label.format)

      image(1, seq(from = 0, to = n.levels), matrix(phase.levels, nrow=1), col = key.cols, breaks = phase.levels, useRaster=T, xaxt='n', yaxt='n', xlab='', ylab='')
      axis(4, lwd=lwd.axis, at=key.marks, labels=NA, tck=0.02, tcl=(par()$usr[2]-par()$usr[1])*legend.params$width-0.04)
      mtext(key.labels, side = 4, at = key.marks, adj=1, line = 0.5 + (1 + 0.5*legend.params$label.digits + line.add)*par()$cex.axis, las=2, font = par()$font.axis, cex=par()$cex.axis)
      text(x = par()$usr[2] + (1.5+legend.params$lab.line)*par()$cxy[1], y=n.levels/2, labels=legend.params$lab, xpd=NA, srt = 270, font = par()$font.lab, cex=par()$cex.lab)
  
      box(lwd = lwd.axis)    
    
      par(new=TRUE, plt = image.plt)  
      
  }    
  
  #######################################################################################
  ## plot phase image
  #######################################################################################
  
  image(axis.1, axis.2, t(Phase), col = key.cols, breaks = phase.levels, useRaster = useRaster,
        ylab = "", xlab = '', axes = FALSE, main = main)         
       
  # plot contour lines?     
  if  ((plot.contour == T) & (is.null(Power.pval) == F)) {      
      contour(axis.1, axis.2, t(Power.pval) < siglvl, levels = 1, lwd = lwd, 
              add = TRUE, col = col.contour, drawlabels = FALSE)
  }
  
  # plot ridge?
  if  (plot.ridge == T) {
       contour(axis.1, axis.2, t(Ridge), levels = 1, lwd = lwd,
              add = TRUE, col = col.ridge, drawlabels = FALSE)
  }     

  # plot cone of influence?
  if (plot.coi == T) {
      polygon(WT$coi.1, WT$coi.2, border = NA, col = rgb(1, 1, 1, 0.5))
  }  
  
  box(lwd = lwd.axis)
  
  ##############################################
  
  # label period axis ?  
  if (label.period.axis == T) {
  
      # check if there is a user-defined period axis specification
      period.axis.default = ( is.null(spec.period.axis$at) )
      
      # first: check conditions of period axis specification
      # period.axis.warning = TRUE leads to default reset
      
      if ( !period.axis.default ) {
      
           # 1. Are tick marks (spec.period.axis$at) interpretable as numbers?
           # is.numeric returns TRUE iff there is no character value but at least one non-NA value which is interpretable as number
           # (further logical entries T and F are interpreted as 1 and 0 which could be accepted here... entry 0 will be removed later on)
           
           # if numeric: are period values all positive? all finite? (NAs, being also interpreted as not finite, will be removed during plotting)
           if ( is.numeric(spec.period.axis$at) ) { 
                period.axis.warning = ( ( sum(spec.period.axis$at<=0, na.rm=T)>0 ) | ( sum(!is.na(spec.period.axis$at)) != sum(is.finite(spec.period.axis$at)) ) )
           }
           else { period.axis.warning = TRUE }
           
           # 2. Are labels (spec.period.axis$labels) appropriate? Either logical and single-valued or non-logical and of the same length as tick marks?
           # is.logical returns TRUE iff there is not any character or numeric value, but there could be NA, even if it is the only value.
           
           # if logical:  single-valued and non-NA? 
           if ( is.logical(spec.period.axis$labels) ) { 
                period.axis.warning = ( period.axis.warning | ifelse(length(spec.period.axis$labels)!=1,TRUE, is.na(spec.period.axis$labels)) )
           }
           
           # if non-logical: do vectors of tick marks and labels have equal length?
           if ( !is.logical(spec.period.axis$labels) ) {
               period.axis.warning = ( period.axis.warning | (length(spec.period.axis$labels) != length(spec.period.axis$at)) ) 
           }
           
      } # end if ( !period.axis.default )
      
      # default reset in case of any warning
      period.axis.default = (period.axis.default | period.axis.warning)
  
      if ( (is.null(periodlab)) | (!is.null(periodlab) & period.axis.warning) ) {periodlab='period'}
    
      if (period.axis.default) {
          period.tick = unique(trunc(axis.2))
          period.tick[period.tick<log2(WT$Period[1])] = NA
          period.tick = na.omit(period.tick)
          period.tick.label = 2^(period.tick)   
          axis(2, lwd = lwd.axis, at = period.tick, labels = NA, tck=periodtck, tcl=periodtcl)
          axis(4, lwd = lwd.axis, at = period.tick, labels = NA, tck=periodtck, tcl=periodtcl)
          mtext(period.tick.label, side = 2, at = period.tick, las = 1, line = par()$mgp[2]-0.5, font = par()$font.axis, cex=par()$cex.axis)
      }
      
      if (!period.axis.default) {
      
          period.tick = log2(spec.period.axis$at)
          period.tick[(period.tick<log2(WT$Period[1]))] = NA
          # (works even if any period.tick value is NA beforehand)   
          period.tick.label = spec.period.axis$labels
          
          # NAs among tick marks?
          which.na = which(is.na(period.tick))
          # NA warning
          if ( length(which.na)>0 ) { period.axis.warning.na = T } 
                    
          # in case of NAs among tick marks:
          # NA tick marks (and corresponding labels) are omitted automatically by the axis command
          
          if ( is.logical(period.tick.label) ) {
               # then period.tick.label has length 1 by criterion
               # The following works even if any period.tick value is NA   
               if (period.tick.label == T) { period.tick.label = 2^(period.tick) }
          }     
          
          axis(2, lwd = lwd.axis, at = period.tick, labels = period.tick.label, tck=periodtck, tcl=periodtcl,
               las = spec.period.axis$las, hadj = spec.period.axis$hadj, padj=spec.period.axis$padj, 
               mgp = par()$mgp - c(0,0.5,0), font = par()$font.axis, cex.axis=par()$cex.axis)
          axis(4, lwd = lwd.axis, at = period.tick, labels = NA, tck=periodtck, tcl=periodtcl)  
          
      } # end if (!period.axis.default)
      
      mtext(periodlab, side = 2, line = par()$mgp[1]-0.5, font = par()$font.lab, cex=par()$cex.lab)
          
  }
  
  #################################

  # label time axis ? --- checks whether time.axis.default or not is applicable (including format and chronology tests of date in case show.date=T
  if (label.time.axis == T) {
  
      # check if there is a user-defined time axis specification
      time.axis.default = ( is.null(spec.time.axis$at) )
  
      # check chronology and format of calendar date in case show.date = T
      if (show.date==T) { 
      
           if (is.element('date',names(series.data))) { my.date = series.data$date } else { my.date = rownames(series.data) }
           
           if (is.null(date.format)) { 
               
                   chronology.warning = inherits(try(as.Date(my.date, tz=date.tz), silent=T),'try-error')
                   
                   if (!chronology.warning) {
                       my.date = as.Date(my.date, tz=date.tz) 
                       chronology.warning = ifelse( sum(is.na(my.date))> 0, TRUE, sum(diff(my.date, tz=date.tz)<0)>0 )
                   }
                   
           }
           if (!is.null(date.format)) {
               
                   chronology.warning = inherits(try(as.POSIXct(my.date, format=date.format, tz=date.tz), silent=T),'try-error')
               
                   if (!chronology.warning) {
                       my.date = as.POSIXct(my.date, format=date.format, tz=date.tz) 
                       chronology.warning = ifelse( sum(is.na(my.date))> 0, TRUE, sum(diff(my.date, tz=date.tz)<0)>0 )
                   }
            }
            if (chronology.warning) { 
                show.date = F
                time.axis.default = T
                timelab = 'index' 
            }
            
      } # end if (show.date==T)      
       
      # first: check conditions of time axis specification
      # time.axis.warning = TRUE leads to default reset
      
      if ( !time.axis.default ) { 
      
           # 1. According to show.date: Are tick marks (spec.time.axis$at) interpretable as numbers or dates?
           # show.date = F requires numeric values
           # show.date = T requires dates
           # mismatch results in default reset
           
           # case show.date==FALSE
           # Are tick marks (spec.time.axis$at) interpretable as numbers?
                      
           # is.numeric returns TRUE iff there is no character value but at least one non-NA value which is interpretable as number
           # (further logical entries T and F are interpreted as 1 and 0 which could be accepted here, non-positive values as well,
           # since values are matched with the true range of values in the axis command)
           
           if (show.date==F) { time.axis.warning = (!is.numeric(spec.time.axis$at)) }
           
           # case show.date==TRUE
           # Are tick marks (spec.time.axis$at) interpretable as dates using given date format and time zone specifications?
           
           # checking the first value given in spec.time.axis$at:
           # inherits is TRUE (and thus gives a time axis warning) iff the first value cannot be formatted, regardless of subsequent values
           # date formatting may produce NAs, however, therefore checking is necessary whether there are NAs ONLY
           # (NA tick marks are omitted automatically by the axis command as long as there are non-NAs)
           
           if (show.date==T) { 
           
               if (is.null(date.format)) { 
                   time.axis.warning = inherits(try(as.Date(spec.time.axis$at, tz=date.tz), silent=T), 'try-error') 
                   if (!time.axis.warning) { time.axis.warning = (sum(!is.na(as.Date(spec.time.axis$at, tz=date.tz)))==0) }
               }
            
               if (!is.null(date.format)) { 
                   time.axis.warning = inherits(try(as.POSIXct(spec.time.axis$at, format=date.format, tz=date.tz), silent=T),'try-error') 
                   if (!time.axis.warning) { time.axis.warning =  (sum(!is.na(as.POSIXct(spec.time.axis$at, format=date.format, tz=date.tz)))==0) }
               }            
           
           } # end if (show.date==T)
           
           if (!time.axis.warning) { 
           # 2. Are labels (spec.time.axis$labels) appropriate? Either logical and single-valued or non-logical and of the same length as tick marks?
           # is.logical returns TRUE iff there is not any character or numeric value, but there could be NA, even if it is the only value.
           
             # if logical:  single-valued and non-NA? 
             if ( is.logical(spec.time.axis$labels) ) { 
                  time.axis.warning = ( ifelse(length(spec.time.axis$labels)!=1,TRUE, is.na(spec.time.axis$labels)) )
             }
           
             # if non-logical: do vectors of tick marks and labels have equal length?
             if (!is.logical(spec.time.axis$labels)) {
                 time.axis.warning = ( length(spec.time.axis$labels) != length(spec.time.axis$at) ) 
             }
           } # end if (!time.axis.warning)
           
      } # end if ( !time.axis.default )
       
      # default reset in case of a warning
      time.axis.default = (time.axis.default | time.axis.warning)
  
  }  # end if (label.time.axis == T)   
   
  ##############################################
  
  half.increment.2 = (axis.2[2]-axis.2[1])/2
  my.ylim = c(min(axis.2)-half.increment.2, max(axis.2)+half.increment.2)
  
  index.condition = ( (show.date==F) | (label.time.axis==F) )
  
  if ((WT$dt != 1) & (index.condition)) {
  
  par(new=TRUE)
  
  # empty plot, but index axis (w/o plotting)
  plot(1:WT$nc, seq(min(axis.2),max(axis.2), length.out=WT$nc),
       xlim = c(0.5,WT$nc+0.5), 
       ylim = my.ylim, 
       type="n", 
       xaxs = "i", yaxs ='i', 
       xaxt='n', yaxt='n', 
       xlab="", ylab="")
       
  }     
       
  ##############################################
  
  # finally: label time axis?  
  if (label.time.axis == T) {     
      
    
      if ( (is.null(timelab) & time.axis.default)  | (!is.null(timelab) & time.axis.warning) ) {timelab= ifelse(show.date, 'calendar date','index')}
            
      #####################
      
      if (show.date == F) {  
       
               if (time.axis.default) {
               
                   A.1 = axis(1, lwd = lwd.axis, labels=NA, tck = timetck, tcl = timetcl)
                   mtext(A.1, side = 1, at = A.1, line = par()$mgp[2]-0.5, font = par()$font.axis, cex=par()$cex.axis)
                   
               } # end if (time.axis.default)
               
               if (!time.axis.default) {
               
                   time.tick = spec.time.axis$at
                   time.tick.label = spec.time.axis$labels
                   
                   # NAs among tick marks?
                   which.na = which(is.na(time.tick))
                   # NA warning
                   if ( length(which.na)>0 ) { time.axis.warning.na = T } 
                   
                   # in case of NAs among tick marks:
                   # NA tick marks (and corresponding labels) are omitted automatically by the axis command
          
                   axis(1, lwd = lwd.axis, at=time.tick, labels=time.tick.label, tck = timetck, tcl = timetcl,
                        las = spec.time.axis$las, hadj = spec.time.axis$hadj, padj=spec.time.axis$padj,
                        mgp = par()$mgp - c(0,0.5,0), font = par()$font.axis, cex.axis=par()$cex.axis)
                       
               } # end if (!time.axis.default)
               
               mtext(timelab, side = 1, line = par()$mgp[1]-1, font = par()$font.lab, cex=par()$cex.lab)
               
      } # end if (show.date == F)
      
      #####################
      
      if (show.date == T) {  
      
               half.increment = difftime(my.date[2],my.date[1], tz=date.tz)/2
                     
               if (is.null(date.format)) { my.xlim = as.Date(range(my.date) - c(half.increment,-half.increment),tz=date.tz) } 
               else { my.xlim = as.POSIXct(range(my.date) - c(half.increment,-half.increment),format=date.format, tz=date.tz) }
                                           
               par(new=TRUE)
               
               if (time.axis.default) {
               
                   # empty plot, but calendar axis (w/ plotting)
                   plot(my.date, seq(min(axis.2),max(axis.2), length.out=WT$nc), 
                        xlim = my.xlim, 
                        ylim = my.ylim, 
                        type="n", 
                        xaxs = "i", yaxs ='i', 
                        yaxt='n', 
                        xlab="", ylab="",
                        lwd = lwd.axis, tck=timetck, tcl=timetcl, 
                        mgp = par()$mgp - c(0,0.5,0), font = par()$font.axis, cex.axis=par()$cex.axis)
                        
               } # end if (time.axis.default)
               
               if (!time.axis.default) {
               
                   # empty plot, but calendar axis (w/o plotting)
                   plot(my.date, seq(min(axis.2),max(axis.2), length.out=WT$nc),
                        xlim = my.xlim, 
                        ylim = my.ylim, 
                        type="n", 
                        xaxs = "i", yaxs ='i', 
                        xaxt='n', yaxt='n', 
                        xlab="", ylab="")
                   
                   # user-defined calendar axis
                   if (is.null(date.format)) { time.tick = as.Date(spec.time.axis$at, tz=date.tz) }
                   if (!is.null(date.format)) { time.tick = as.POSIXct(spec.time.axis$at, format=date.format, tz=date.tz) }
                   
                   time.tick.label = spec.time.axis$labels
                   
                   # NAs among tick marks?
                   which.na = which(is.na(time.tick))
                   # NA warning
                   if ( length(which.na)>0 ) { time.axis.warning.na = T } 
                   
                   # in case of NAs among tick marks:
                   # NA tick marks (and corresponding labels) are omitted automatically by the axis command
                   
                   if ( is.logical(time.tick.label) ) {
                        # then time.tick.label has length 1 by criterion and IS TRUE OR FALSE
                        if (time.tick.label == T) {time.tick.label = time.tick}
                   }   
          
                   axis(1, lwd = lwd.axis, at=time.tick, labels=time.tick.label, tck = timetck, tcl = timetcl,
                        las = spec.time.axis$las, hadj = spec.time.axis$hadj, padj=spec.time.axis$padj,
                        mgp = par()$mgp - c(0,0.5,0), font = par()$font.axis, cex.axis=par()$cex.axis)
                   
               } # end if (!time.axis.default)             
               
               mtext(timelab, side = 1, line = par()$mgp[1]-1, font = par()$font.lab, cex=par()$cex.lab)      
               
      }  # end if (show.date == T)  
      
  }  # end finally: if (label.time.axis == T)  
   
  #######################################################################################
  ## apropos graphical parameters
  #######################################################################################
  
  # reset contour line options
  options(default.options)
  
  # reset graphical parameters?
  if (graphics.reset == T) {
      par(op)
  }  
  
  if (chronology.warning) { 
      warning("\nPlease check your calendar dates, format and time zone: dates may not be in an unambiguous format or chronological. The default numerical axis was used instead.") 
  }  
  
  if (period.axis.warning == T) {
      warning("\nPlease check your period axis specifications. Default settings were used.")
  }
  if (period.axis.warning.na == T) {
      warning("\nNAs were produced with your period axis specifications.")
  }
  
  if (time.axis.warning == T) {
      warning("\nPlease check your time axis specifications. Default settings were used.")
  }
  if (time.axis.warning.na == T) {
      warning("\nNAs were produced with your time axis specifications.")
  }
  
  # output of graphical parameters
  
  output = list(op = op, image.plt = image.plt, legend.plt=legend.plt)
  class(output) = "graphical parameters"

  out("Class attributes are accessible through following names:\n")
  out(names(output), "\n")    
    
  return(invisible(output)) 
  
}