Source Code for Renewable Grid Integration Series: Part 1

Source Code – Programmed in R 3.2.2

Data is accessed from http://www.aemo.com.au:

http://www.aemo.com.au/Electricity/Data/Price-and-Demand/Aggregated-Price-and-Demand-Data-Files/Aggregated-Price-and-Demand-2011-to-2016

http://www.nemweb.com.au/REPORTS/ARCHIVE/Next_Day_Actual_Gen/

http://www.nemweb.com.au/REPORTS/ARCHIVE/Next_Day_Dispatch/

http://www.aemo.com.au/About-the-Industry/Registration/Current-Registration-and-Exemption-lists

Data must be loaded into working directory.

File: load_dataset.R

# This file sources data for system demand and total generation. It then filters
# the generation data for wind and solar generation. It also loads relevant
# libraries and files

library(XML)
library(plyr)

# Create a listing of files in the directory

listing <- list.files()

# Create a datafile with all of the load data

loadfiles <- grepl(“DATA20″,listing)
loadlist <- listing[loadfiles]
l <- length(loadlist)
load_data <- read.csv(loadlist[1])

for (i in loadlist[2:l]) {
curr_file <- read.csv(i)
load_data <- rbind(load_data,curr_file)
}

load_data$SETTLEMENTDATE<-as.POSIXct(strptime(load_data$SETTLEMENTDATE,”%Y/%m/%d %H:%M:%S”))
write.csv(load_data,file=”usage_info.csv”)

# Create a datafile for Scheduled generation data

genfiles <- grepl(“PUBLIC_NEXT_DAY_DISPATCH”,listing)
genlist <- listing[genfiles]
l <- length(genlist)
inter <- read.csv(“Windsolar.csv”,header=FALSE)
gen_data <- list()

# Initialisation and testing of first gen_data input

for (i in 1:l) {
curr_file <- read.csv(genlist[i],skip=1)
curr_file <-subset(curr_file,UNIT_SOLUTION==”UNIT_SOLUTION”,select=c(SETTLEMENTDATE,DUID,INITIALMW))
curr_file <- subset(curr_file,curr_file$DUID %in% inter$V1)
gen_data[[i]] <- curr_file
}

# Create a datafile for Non-Scheduled generation data

genfiles <- grepl(“PUBLIC_NEXT_DAY_ACTUAL”,listing)
genlist <- listing[genfiles]
n <- length(genlist)

for (i in 1:n) {
curr_file <- read.csv(genlist[i],skip=2,header=FALSE)
curr_file <- subset(curr_file,select=V5:V7)
names(curr_file)<-c(“SETTLEMENTDATE”,”DUID”,”INITIALMW”)
curr_file <- subset(curr_file,curr_file$DUID %in% inter$V1)
gen_data[[i+l]] <- curr_file
}

gen_data<-rbind.fill(gen_data)
gen_data$SETTLEMENTDATE <- as.POSIXct(strptime(gen_data$SETTLEMENTDATE,”%Y/%m/%d %H:%M:%S”))
write.csv(gen_data,file=”gener_info_mult.csv”)

# Load data is stored in “usage_info.csv”.
# Generation data is stored in “gener_info_mult.csv”
# We have incorporated an excel shortcut workaround that assigns the relevant NEM
# region to each generator in the list. This can also be undertaken via R functionality
# as required.

# Additional (currently unused) data is collected for further work.

File: analyse_dataset.R

# This file aggregates all the wind data per region and offsets it against demand,
# producing an Rdata file with all Net Load information.

library(XML)
library(plyr)
library(ggplot2)
library(lattice)
library(gtools)

# Loads relevant files

loadinfo <- read.csv(“usage_info.csv”)
loadinfo$SETTLEMENTDATE <- as.POSIXct(strptime(loadinfo$SETTLEMENTDATE,”%Y-%m-%d %H:%M:%S”))
geninfo <- read.csv(“gener_info_mult.csv”)
geninfo$SETTLEMENTDATE <- as.POSIXct(strptime(geninfo$SETTLEMENTDATE,”%m/%d/%y %H:%M”))

# Aggregation of wind and solar data

intermit <- aggregate(INITIALMW ~ Region + SETTLEMENTDATE,data=geninfo,sum)

jdata <- merge(loadinfo, intermit, by.x=c(“SETTLEMENTDATE”,”REGION”), by.y=c(“SETTLEMENTDATE”,”Region”),all.x = TRUE, all.y=FALSE)

# NAs are assigned 0.

nas <- is.na(jdata$INITIALMW)
jdata$INITIALMW[nas]<-0

rm(geninfo)
rm(loadinfo)

jdata$Mthday <- format(jdata$SETTLEMENTDATE,’%m-%d’)
jdata$NetLoad <- jdata$TOTALDEMAND-jdata$INITIALMW
jdata$Yr <- format(jdata$SETTLEMENTDATE,’%Y’)
jdata$tim <- format(jdata$SETTLEMENTDATE,’%H:%M’)
jdata$dateun <- format(jdata$SETTLEMENTDATE,’%Y %H:%M’)

# calculate the mean and standard deviation for each 30 min period, aggregating
# across the month to smooth for weather variations.

meandata <- aggregate(NetLoad ~ dateun + Yr + tim + REGION,data=jdata,mean)
sddata <- aggregate(NetLoad ~ dateun + Yr + tim + REGION,data=jdata,sd)
statedata <- merge(meandata,sddata,by=c(“dateun”,”REGION”),all=FALSE)

# total data for the NEM is aggregated below.

tmeandata <- aggregate(NetLoad ~ dateun+Yr + tim,data=meandata,sum)
tsddata <- aggregate(NetLoad ~ dateun + Yr + tim,data=jdata,sd)
tdata <- merge(tmeandata,tsddata,by=”dateun”,all=FALSE)
tdata$REGION <- “NEM”
alldata <- rbind(statedata,tdata)

# final data is saved to an image.

save.image(file=”alldata.RData”)

File: graph_dataset.R

# This file graphs the datasets for the relevant data and can be broken down by
# NEM region and subregion.

library(ggplot2)
library(lattice)
library(gtools)
library(ggthemes)

load(file=”alldata.RData”)

reginput <- readline(prompt =”Enter Region (NEM, SA1, NSW1, QLD1, VIC1, TAS1): “)
subdata <- subset(alldata,REGION==reginput)

# optional subset for a specific day if required (Un-comment line 15):
# jdata <- subset(jdata,Mthday==”09-30”)

# Graph of Intraday Data over Years 2013-2015

p <- ggplot(subdata,aes(tim.x,NetLoad.x,colour=as.factor(Yr.x)))
# Choice of theme: p <- p + theme_classic()
p <- p+ geom_line(aes(colour=Yr.x,group=Yr.x))
p <- p + theme(legend.title=element_blank())
p <- p + labs(title = paste(“Intraday Net Load (September) for”,reginput), x= “Time”, y = “Net Load (MW)”)
p <- p + scale_x_discrete(breaks=c(“00:00″,”06:00″,”12:00″,”18:00″,”23:30″))
print(p)

# # Graph of Intraday Data with 95% Confidence Envelope

subdata15 <- subset(subdata,Yr.x==”2015”)
subdata15$timeo <- format(subdata15$tim.x, format = “%H:%M”)
subdata15$timeo <- as.POSIXct(subdata15$timeo, format = “%H:%M”, tz = “Australia”)

a <- ggplot(subdata15,aes(timeo,NetLoad.x,colour=as.factor(Yr.x)))
a <- a + geom_line(aes(colour=Yr.x,group=Yr.x))
minl <- subdata15$NetLoad.x-(2*subdata15$NetLoad.y)
maxl <- subdata15$NetLoad.x+(2*subdata15$NetLoad.y)
a <- a + geom_ribbon(aes(ymin=minl,ymax=maxl),alpha=0.3,colour=NA)
a <- a + theme(legend.title=element_blank(),legend.position=”none”)
a <- a + labs(title = paste(“Intraday Net Load (September 2015) for”,reginput), x= “Time”, y = “Net Load (MW)”)
a <- a + scale_x_datetime(date_labels = “%H:%M”)
print(a)

 

 

 

Advertisements

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s