############## RClimate Script: Polar Amplification - GISS Temp ############################
## Script stored on http://chartsgraphs.wordpress.com account for Users to source()       ##                                                                                       ## 
## Download and process Monthly CO2 Data File                                             ##
## Developed by D Kelly O'Day to demonstrate use of source() function for climate data    ##
##                   http:chartsgraphs.wordpress.com    1/27/10                           ##
###################################################################

 par(mar=c(4,4,3,1)); par(oma=c(3,1,1,1));
 par(ps=10); par(las=1)

## Read GISSDecadal Anomaly data by Lattitude zone
 link  <- "http://processtrends.com/Files/GISS_zonal_mean_anon_2000_2009.csv"
 df<- read.table(link, skip=0,sep=",", na.strings=c(9999.0000)) 
 names(df) <- c("lat", "anom")

## calc Overall  avg & max
 overall_avg = signif(mean(df$anom, na.rm=T), 3)
 overall_max <- signif(max(df$anom, na.rm=T), 3)

## Plot annotation
  title <- "Polar Amplification: 2000-2009\n GISS (Land & Ocean) Temperature Anomaly by Lattitude Zone"
  x_lab <- expression(paste("Zonal Mean Decadal Anomaly - ",degree,"C (1951-1980 baseline)"))
  y_lab <- expression(paste("Latitude - ",degree))
  avg_note <- paste("Global mean =",overall_avg)
  max_note <-  paste("Zonal max = ", overall_max)
  data_note <- "Data Source: NASA GISS - http://data.giss.nasa.gov/gistemp/maps/"
## Plot
 y_ax <- seq(-90,90,30)  
 plot(df$anom, df$lat, ylim = c(-90,90), type = "l", col = "black",xlim=c(0,2.5),
     main=title, xlab = x_lab,ylab = y_lab, axes=F, xaxs="i",yaxs="i",
     cex.main=0.85,cex.lab=0.85)
 axis(1, at=NULL,col = "black", cex.axis=0.85)   
 axis(2, at=y_ax, col = "black", cex.axis=0.85)
 abline(h=0, col="grey")
 # Add overall avg & max point markers & lines
    points(overall_avg,-90,pch=4, col = "blue")
    points(overall_max,-90, pch=4, col = "red")
    points(x=c(overall_avg, overall_avg), y=c(-90, 37),col = "blue", type="l")
    points(x=c(overall_max, overall_max), y=c(-90, 80),col = "red", type="l")
 # Add blank text rectaangles & text annotations for overall avg & max
    rect(0.5,6,1 , 15, col = "white", border = "white") 
    text(overall_avg-0.05, 11, avg_note,  cex=0.75, adj=0,col="blue")
    rect(1.6,-8,1.8 ,-15, col = "white", border = "white") 
    text(overall_max, -11, max_note, cex=0.75, adj=0.5, col="red")

 # Plot Outer margin annotation
    my_date <- format(Sys.time(), "%m/%d/%y")
    mtext("D Kelly O'Day - http://chartsgraphs.wordpress.com", side = 1, line = 1, cex=0.8, outer = T, adj = 0)
    mtext(my_date, side = 1, line =1, cex = 0.8, outer = T, adj = 1) 
    mtext(data_note, side = 1, line =0, cex = 0.8, outer = T, adj = 0.5)