############## RClimate Script: Global SST Anomaly - C   ###################################
## Script stored on http://chartsgraphs.wordpress.com account for Users to source()       ##                                                                                       ## 
## Download and process     NCDC - NOAA monthly global sst anomaly                        ##
## Developed by D Kelly O'Day to demonstrate use of source() function for climate data    ##
##                   http:chartsgraphs.wordpress.com    1/31/10                           ##
############################################################################################
 par(ps=10); par(las=1); par(oma=c(2.5,1,1,1)); par(mar=c(2,4,2,0))
## Read Source Data
 link <- "ftp://ftp.ncdc.noaa.gov/pub/data/anomalies/monthly.ocean.90S.90N.df_1901-2000mean.dat"
 data <- read.table(link)
 names(data) <- c("yr", "mo","ssta")
 yr_frac <- data$yr + (data$mo-0.5)/12
 data <- data.frame(yr_frac, data)

## Determine month - yr of last reading
    c <- nrow(data)                        # Find number of data rows
    l_mo <- data$mo[c]                     # Find last month of data
    l_yr <- data$yr[c]                     # Find last year of data
    l_ssta <- signif(data$ssta[c],4)       # Find last ssta reading

## Calc decade averages
  dec_mean<- as.numeric(14)
  dec_st <- as.numeric(14)
  dec_end <- as.numeric(14)
  yr_n <- as.numeric(data$yr)
  base_yr <- 1870
  dec_n <-  (as.numeric((yr_n - base_yr) %/% 10) * 10) + base_yr
  data <- data.frame(data, dec_n)
  for (i in 1:13) {dec_st[i] = base_yr+ i*10
     dec_sub <- subset(data, dec_n == dec_st[i], na.rm=T)
     dec_mean[i] <- mean(dec_sub$ssta)
     }
dec_st[14] <- 2010              # Need to have for last step line across decade
dec_mean[14] <- dec_mean[13]
dec<- data.frame(dec_st, dec_mean)

### Generate last point & Yr-month notes
  title <- paste("Monthly Global SST Anomaly Trends \nNODC - NOAA (1/1880 - ", l_mo, "/", l_yr,")")
  l_mo_yr <- paste(l_mo, "/", l_yr, " @ ",l_ssta, sep="")
  y_lab <- expression(paste("SST Anomaly - ",degree, "C (1901 - 2000 baseline)", sep=""))



## Plot Function ##################################################################################
plot_func<- function() { plot(data$yr_frac, data$ssta, type = "l" , col = "grey",
       xlab = "", ylab = y_lab, main = title, cex.main =0.85, cex.lab=0.85)
 points(yr_frac[c], l_ssta, type = "p",pch=16, col = "red")
# points(1959, -0.4, type = "p",pch=16, col = "red")
# text(1961,-0.4, l_mo_yr, col = "red",adj=0)
 points(dec_st, dec_mean, type = "s", col = "blue")

legend(1955,-0.3, c(l_mo_yr, "Decade Mean Anomaly", "Monthly Anomaly"), col = c("red","blue", "grey"),
       text.col = "black", lty = c(0, 1,1), pch=c(16,0,0), pt.cex=c(1.1,0,0),
       merge = T, bg = "white", bty="n", cex=0.85)

# Generate and add bottom footer KOD , source, System date notes
       source_note <- paste("Data Source: ", link)
       mtext("D Kelly O'Day - http://chartgraphs.wordpress.com", 1,1, adj = 0, cex = 0.8, outer=TRUE)
       mtext(format(Sys.time(), "%m/%d/ %Y"), 1, 1, adj = 1, cex = 0.8, outer=TRUE)
       mtext(source_note, 1,0,adj=0.5, cex=0.8, outer=T)
#############################################################################################
 }

f_to <- "C:\\process_trends_fp_Jan_2010\\images\\RClimate_SST_A_latest.png"  
## Step 3 - call png()
  png(file=f_to,width=450,height=500, res=72)
## Step 4 - call plot_func()
  plot_func()
# Step 5 - turn off png dev
  dev.off()
# ## Step 6 plot to console
  plot_func()