Dowloaded the relevant packages which were
library(ncdf4);
library(ncdf4.helpers);
## Warning: package 'ncdf4.helpers' was built under R version 3.5.2
library(maps);
library(RColorBrewer)
## Warning: package 'RColorBrewer' was built under R version 3.5.2
library(lattice)
library(latticeExtra)
As well as this variables and open plots were cleared before the workshop begun
graphics.off()
rm(list=ls())
The work directory was set up as follows to allow the assignement to take place
setwd("C:/Users/Kennf/Documents")
the NetCDF data was loaded up as it is a widely used format for climate data.
nc <- nc_open(file.path(getwd(),"HadCRUT.4.6.0.0.median.nc")); # open netcdf file
tobs <- ncvar_get(nc,"temperature_anomaly"); # load temperature anomaly
ts <- ncvar_get(nc,"time");
lat <- ncvar_get(nc,"latitude")
lon <- ncvar_get(nc,"longitude")
tunits<-ncatt_get(nc,"time",attname="units")
tustr<-strsplit(tunits$value, " ")
date<-as.character(as.Date(ts,origin=unlist(tustr)[3]))
tobs[t<-90] = NA;
tobs[t>490] = NA;
next the GMIT I estimate was plotted using the following function.
gmst = read.csv(file.path(getwd(),'global_temps_monthly.csv'),skip = 1,header=TRUE)
gmst <- data.frame(Year = gmst$Year, t = rowMeans(gmst[, 2:13], na.rm = T))
plot(gmst$Year, gmst$t, type='l')
plot(gmst$Year, gmst$t, type='l', col="Red", main="GMST", xlab="Years", ylab="T Anomaly")
The second graph that was made was using the ‘filled.contour’ and ‘map’ using the selected year 1896-07-16
year_index <- which(date=="1896-07-16")
mapmat= tobs[,,year_index];
int=seq(-6,6,length.out = 11);
colors <- rev(brewer.pal(10, "RdYlBu"))
rgb.palette <- colorRampPalette(colors)
#filled.contour(lon, lat, mapmat, color.palette=rgb.palette, levels=int,
# plot.title=title(main=paste0('HadCRUT Anom. in ', date[year_index]),
# xlab="Latitude",ylab="Longitude", cex.lab=1.5),
# plot.axes={axis(1, cex.axis=1.5);
# axis(2, cex.axis=1.5);
# map('world', add=TRUE);
# grid()},
# key.title=title(main="[oC]"),
# key.axes={axis(4, cex.axis=1.5)})
Once the map was created, the next goal was to find and plot in the same way two months that illustrate inhomogenties of the dataset through time
dim(tann) #[1] 72 36 168
The following function was put in to create the Merdional meansand mean of zonal means
#
#merid_means<-apply(tann,c(1,3),mean,na.rm=TRUE)
#mean_merid_means <- apply(merid_means,2,mean,na.rm=TRUE)
#zonal_means<-apply(tann,c(2,3),mean,na.rm=TRUE)
#mean_zonal_means <- apply(zonal_means,2,mean,na.rm=TRUE)
#plot(gmst$Year, gmst$t, type='l', col="Red", main="GMST", xlab="Years", ylab="T Anomaly")
#lines(years, mean_merid_means, col="blue")
#lines(years, mean_zonal_means, col="green")
#legend("topleft", legend=c("gis temp", "Merdional means", "Mean of zonal means"),col=c("red", "blue", "green"), lty=1)
#w<-replicate(length(years),t(replicate(length(lon),cos(lat*pi/180))))
# gmst_weighted_merid <- apply(tann*w,c(1,3),sum,na.rm=TRUE)/sum(cos(lat*pi/180)
#gmst_weighted <- apply(gmst_weighted_merid,2,mean,na.rm=TRUE)
#lines(years, gmst_weighted, col="black")
#legend("topleft", legend=c("gis temp", "Merdional means", "Mean of zonal means", "Weighted means"),col=c("red", "blue", "green", "pink"), lty=1)
The final part of this assignment was using the following function to get the following result
which(gmst$Year=="1960"):which(gmst$Year==max(gmst$Year))
## [1] 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97
## [18] 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114
## [35] 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131
## [52] 132 133 134 135 136 137 138
# [1] 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103
# [24] 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126
# [47] 127 128 129 130 131 132 133 134 135 136 137 138