Skip to content
Snippets Groups Projects
Commit 4230b011 authored by pbac's avatar pbac
Browse files

more

parent 6330692b
No related branches found
No related tags found
No related merge requests found
......@@ -3,8 +3,8 @@
library(influxdbr)
## connection
con <- influx_connection(host="localhost", port=8086)
##con <- influx_connection(host="vmpbac1.compute.dtu.dk", port=8086)
##con <- influx_connection(host="localhost", port=8086)
con <- influx_connection(host="vmpbac1.compute.dtu.dk", port=8086)
##influx_ping(con)
## ## Source the files in the "functions" folder
......
## ----init, echo=FALSE, cache=FALSE--------------------------------------------
##----------------------------------------------------------------
## Set the working directory
setwd(".")
## Source the files in the "functions" folder
files <- dir("functions",full.names=TRUE)
for(i in 1:length(files)) source(files[i])
##----------------------------------------------------------------
## ----prbsacf------------------------------------------------------------------
##----------------------------------------------------------------
## Generate a PRBS signal
## Use the function defined in the file "functions/prbs.R", which generates a PRBS signal
## - n is the length of the register
## - initReg just needs to be some initial value of 1,2,...
## it is the initial value of the registers and therefore only
## determines the start of the cycle
## - lambda is the length of the smallest period in
## which the signal can change, given in samples
x <- prbs(n=6)
acf(x, lag.max=length(x))
## prbs
plot(x,type="s")
acf(x, lag.max=length(x))
## See to PRBS
x <- prbs(n=6)
## LagWithCycling to the right
y <- lagWithCycling(x, lag=10)
## cross cor.
ccf(x,y, lag.max=length(x))
## ----prbshist-----------------------------------------------------------------
##----------------------------------------------------------------
## Generate the PRBS signals for the PRBS1 experiment, where a single signal controls all heaters in the building
## PRBS med n=6, lambda=4
## - Smallest period in one state for 5 minute sample period is then 4*5min=20min
## - Settling time of the system (T_s in the Godfrey 1980 paper) below the period (T_0 in the paper) in: lambda * (2^n-1) * 5 / 60 = 21 hours
n <- 5
lambda <- 4
(x <- prbs(n, 37, lambda))
length(x) / 24
(xper <- c(0,cumsum(abs(diff(x)))))
perlens <- unlist(lapply(split(xper,xper), length))
table(perlens)
hist(perlens, breaks=(0:n)*lambda+lambda*0.5)
## ----prbsplotts---------------------------------------------------------------
plot(x, type="s")
## ----prbsplotts2--------------------------------------------------------------
## Create sequence
testseq <- function(n, lambda, nhours, starthour, Tmin, Tmax, plotit=TRUE){
x <- prbs(n, 42, lambda)
t1 <- starthour + ((1:length(x))[abs(c(0,diff(x)))==1] / length(x)) * nhours
## Round to 5 min steps
t1 <- round(t1*12) / 12
val1 <- rep(c(Tmin,Tmax), length(t1))[1:length(t1)]
##
if(plotit){ plot(t1, val1, type="s") }
##
return(data.frame(t=t1, Tset=val1))
}
## val <- testseq(n=5, lambda=1, nhours=96, Tmin=15, Tmax=24)
## diff(val$t)
## Four days between 10 and 25 degrees C
val1 <- testseq(n=5, lambda=1, nhours=96, starthour=0, Tmin=10, Tmax=25)
## Two days between 17 and 20 degrees C
val2 <- testseq(n=5, lambda=1, nhours=48, starthour=96, Tmin=17, Tmax=20)
## Two days between 20 and 22 degrees C
val3 <- testseq(n=5, lambda=2, nhours=48, starthour=96+48, Tmin=20, Tmax=22)
## One day between 21 and 22 degrees C
val4 <- testseq(n=5, lambda=1, nhours=24, starthour=96+48+48, Tmin=21, Tmax=22)
## Combined
X <- rbind(val1, val2, val3, val4)
plot(X$t/24, X$Tset, type="s")
## ----outputs------------------------------------------------------------------
nms <- read.table("../data/outchannels.csv", header=TRUE, sep=",", as.is=TRUE)$ResampleVarname
nms <- nms[grep("^Room_C", nms)]
rooms <- sapply(strsplit(nms, "_"), function(x){x[2]})
rooms <- rooms[order(as.numeric(gsub("[[:alpha:]]","",rooms)))]
lagval <- numeric(length(rooms))
lagval[rooms=="C0.00"] <- 0
lagval[rooms=="C0.00B"] <- 0
lagval[rooms=="C0.01"] <- 1
lagval[rooms=="C0.02"] <- 1
lagval[rooms=="C0.03"] <- 1
lagval[rooms=="C0.04"] <- 1
lagval[rooms=="C0.05"] <- 0#?
lagval[rooms=="C0.06"] <- 2
lagval[rooms=="C0.07"] <- 3
lagval[rooms=="C0.07A"] <- 3
lagval[rooms=="C0.08"] <- 4
lagval[rooms=="C0.10"] <- 0#?
lagval[rooms=="C0.11"] <- 5
lagval[rooms=="C0.12"] <- 7
lagval[rooms=="C0.13"] <- 6
lagval[rooms=="C0.14"] <- 7
lagval[rooms=="C0.15"] <- 8
lagval[rooms=="C1.00"] <- 3
lagval[rooms=="C1.02"] <- 4
lagval[rooms=="C1.03"] <- 5
lagval[rooms=="C1.04"] <- 6
lagval[rooms=="C1.07"] <- 7
lagval[rooms=="C1.08"] <- 0
lagval[rooms=="C1.09"] <- 1
lagval[rooms=="C1.12"] <- 2
lagval[rooms=="C1.13"] <- 3
lagval[rooms=="C2.01"] <- 6
lagval[rooms=="C2.03"] <- 7
lagval[rooms=="C2.04"] <- 8
lagit <- function(val, nhours){
maxk <- max(lagval)
starthour <- floor(val$t[1] / 24) * 24
lapply(1:length(rooms), function(i){
k <- lagval[i]
val$t <- (val$t + k/maxk * nhours) %% nhours
val$t <- val$t + starthour
return(val[order(val$t), ])
})
}
L1 <- lagit(val1, 96)
L2 <- lagit(val2, 48)
L3 <- lagit(val3, 48)
L4 <- lagit(val4, 24)
L <- lapply(1:length(rooms), function(i){
rbind(L1[[i]], L2[[i]], L3[[i]], L4[[i]])
})
## ----plotrooms----------------------------------------------------------------
##setpar("ts", mfrow=c(5,5))
lapply(1:length(L), function(i){
x <- L[[i]]
plot(x$t, x$Tset, type="s", ylab="")
title(main=rooms[i], line=-1)
})
##dev.copy2pdf(file="setpoints/setpoints.pdf")
## ----write-to-files-----------------------------------------------------------
lapply(1:length(L), function(i){
x <- L[[i]]
## Set starting date and time in hours
x$t <- asp("2019-04-13 00:00") + x$t * 3600
## Write to separate .csv files
write.table(x, file=pst("outchannels/",rooms[i],".csv"), sep=",", row.names=FALSE)
})
## Write the experiment to the incluxdb
# Start of the averaging period
tstartAvg <- ct("2021-02-26 22:00")
# End of the averaging period
tendAvg <- ct("2021-02-27 06:00")
Tset <- 22
# Start of the experiment
tstart <- trunc(Sys.time(), units="hour") + 3600 #ct("2021-02-26 21:00")
# End of the experiment
(tend <- tstart + 43 * 3600)
library(influxdbr)
## connection
##con <- influx_connection(host="localhost", port=8086)
con <- influx_connection(host="vmpbac1.compute.dtu.dk", port=8086)
##influx_ping(con)
## ## Source the files in the "functions" folder
## files <- dir("functions",full.names=TRUE)
## for(i in 1:length(files)) source(files[i])
##library(onlineforecast)
library(devtools)
load_all(as.package("~/g/onlineforecast"))
##
wd <- "../buildingctrl/main/"
source(paste0(wd,"lib/r/read_idtags.R"))
source(paste0(wd,"lib/r/write_values.R"))
## Find the output channels, to be able to check that the channelids are correct
idtags <- read_idtags(wd)
##
i <- grep("^Room_.*spTemperature__northq", idtags$VarName)
##
(nmRooms <- idtags[i, "VarName"])
(nmTforwardRC <- idtags[grep("^HX_RC_spTforward", idtags$VarName), "VarName"])
(nmTforwardVE <- idtags[grep("^HX_VE_spTforward", idtags$VarName), "VarName"])
##
## ------------------------------------------------------------------------
## Generate constant experiment
X <- data.frame(t=seq(tstart, tend, by=0.5*60*60))
range(X$t)
## Set RegI and RegP
##grep("^[HX_].*_spReg", idtags$VarName, value=TRUE)
## Larger I means less integration! thus its probably reciprocal: see 16. Nov. 2019 stepping.
X$HX_RC_spRegP_building_ltech <- 0.5 #rep(seq(0.3,by=0.1,len=6), each=2)
X$HX_RC_spRegI_building_ltech <- 80 #rep(seq(30,by=10,len=6), each=2)
X$HX_VE_spRegP_building_ltech <- 0.5
X$HX_VE_spRegI_building_ltech <- 80 #rep(seq(30,by=10,len=6), each=2)
## The set points
set.seed(9387)
val <- #c(t(matrix(c(sample(seq(20,42,len=8)), sample(seq(48,65,len=8))), ncol=2)))
X[ ,nmTforwardRC] <- 60#rep(c(50,20), 6)
X[ ,nmTforwardVE] <- 45#c(val)#sample(seq(20,65,len=17))#c(65,60,55,50,45,40,35,30,25,30,35,40,45,50,55,60,65))#rep(c(50,20), 6)
## ------------------------------------------------------------------------
## Room set points
head(idtags)
L <- list()
db <- "raw"
rooms <- unlist(getse(strsplit(nmRooms, "_"), 2))
for(i in 1:length(rooms)){
val <- idtags[which(idtags$Node == rooms[i] & idtags$Var == "Temperature" & idtags$VarPrefix == "" & idtags$direction == "in" & idtags$Connection == "elsys"), ]
if( nrow(val) > 0 ){
x <- influx_select(con, db, "value", measurement = "observations", where = pst("channelid = '", val$channelid[1] ,"'"))
res <- x[[1]][[1]]
## if(nrow(val) >= 2){
## for(ii in 2:length(val)){
## x <- influx_select(con, db, "value", measurement = "observations", where = pst("channelid = '", val$channelid[ii] ,"'"))
## res <- merge(res, x[[1]][[1]])
## }
## }
L[[i]] <- data.frame(t = xts:::index.xts(res), val = unlist(res$value))
}
}
# SET IT FOR EACH ROOM AND WRITE IT!!
# First set to 22 for all rooms
plot(0, type="n", xlim=c(tstartAvg,tendAvg), ylim=c(18,24))
for(i in 1:length(nmRooms)){
X[ ,nmRooms[i]] <- Tset
x <- L[[i]][in_range(tstartAvg, L[[i]]$t, tendAvg), ]
if(!is.null(x)){
X[ ,nmRooms[i]] <- Tset + (Tset - mean(x$value, na.rm=TRUE) )
lines(x$t, x$value, col=i)
}
}
## ------------------------------------------------------------------------
plot_ts(X, c("RegP","RegI","spTforward"), plotfun=function(x,y,col){lines(x,y,col=col,type="s")})
plot_ts(X, nmRooms, plotfun=function(x,y,col){lines(x,y,col=col,type="s")})
Xres <- X
#Xres$t[nrow(Xres)] <- ct("2019-12-24")
Xres
##plot.ts(Xres, type="s")
## Replace the column names with channelid
it <- which(names(Xres) == "t")
names(Xres)[-it] <- unlist(sapply(names(Xres)[-it], function(nm){ idtags$channelid[idtags$VarName == nm] }))
#write.table(Xres, "../buildingctrl/sender/senddata.csv", sep=",", row.names=FALSE)
## Copy to the server
#system("scp ../buildingctrl/sender/senddata.csv pbac@vmpbac1.compute.dtu.dk:~/buildingctrl/sender/")
## Write directly in influxdb
tout <- ct(approx(1:nrow(X), X$t, seq(1,nrow(X),by=1/4))$y)
iseq <- 2:ncol(Xres)
Xres2 <- data.frame(t=tout, sapply(iseq, function(i){
approx(Xres$t, Xres[ ,i], tout, "constant")$y
}))
names(Xres2) <- names(Xres)
# Delete all set points ahead in time
influx_query(con=con, db=db, query=pst("DELETE FROM observations WHERE Connection = 'northq' AND direction = 'out' AND time > '", as.character(tstart), "'"), return_xts = FALSE)
#
write_values(Xres2, idtags, con, colnamesTag = "channelid")
## ## ------------------------------------------------------------------------
## ## Generate the experiment
## Xsetup <- data.frame(nhours = c(48,24,24),
## Tmin = c(10,30,40),
## Tmax = c(65,57,50),
## registerLen = c(4,4,4),
## seeds = c(10,24,738),
## shortPerLen = c(1,1,1))
## X <- experiment_prbs_Tforward(Xsetup,
## nameTforwardRC=nmTforwardRC,
## nameTforwardVE=nmTforwardVE,
## TforwardVE=10,
## nameRooms=nmRooms,
## Troom=22)
## ## plot(X$t, X$HX_RC_spTforward_building_ltech, type="s", xlab="Hours")
## ## plot_ts.data.frame(X, c("Tforward","^Room"))
## ## Remember that it is GMT start time
## Xres <- gen_equidistant_ts(X, tstart="2019-10-24 14:00", tdelta=3600, t_round_fun=trunc, units="hours")
## plot_ts.data.frame(Xres, c("Tforward","^Room"))
## ## write_values(Xres, idtags, con, colnamesTag="VarName")
## source("~/r/sourceFunctions.R")
## pdf("~/tmp/plot.pdf", height=20, width=20)
## setpar("ts", mfrow=c(length(L),1))
## for(i in 2:ncol(Xres)){
## plot(Xres$t, Xres[ ,i], type="b")
## title(main=names(X)[i], line=-1)
## }
## plotTSXAxis(Xres$t)
## dev.off()
## system("evince ~/tmp/plot.pdf &")
## Write the experiment to the incluxdb
library(onlineforecast)
set.seed(389)
# Start of the experiment
tstart <- ct("2021-03-05 01:00")
tend <- ct("2021-03-15 16:00")
library(influxdbr)
## Source the files in the "functions" folder
files <- dir("functions",full.names=TRUE)
for(i in 1:length(files)) source(files[i])
## connection
##con <- influx_connection(host="localhost", port=8086)
con <- influx_connection(host="borgerskolen.centerdenmark.com", port=8086)
##influx_ping(con)
## ## Source the files in the "functions" folder
## files <- dir("functions",full.names=TRUE)
## for(i in 1:length(files)) source(files[i])
##
wd <- "../buildingctrl/main/"
source(paste0(wd,"lib/r/read_idtags.R"))
source(paste0(wd,"lib/r/write_values.R"))
## Find the output channels, to be able to check that the channelids are correct
idtags <- read_idtags(wd)
##
i <- grep("^Room_.*spTemperature__northq", idtags$VarName)
##
(nmRooms <- idtags[i, "VarName"])
(nmTforwardRC <- idtags[grep("^HX_RC_spTforward", idtags$VarName), "VarName"])
(nmTforwardVE <- idtags[grep("^HX_VE_spTforward", idtags$VarName), "VarName"])
##
## ------------------------------------------------------------------------
## Night values
## Generate room experiment, where all rooms are following same Tset, with \pm 1 PRBS making differences between the rooms
nmRooms
testseq <- function(n, nhours, starthour, Tmin, Tmax, lambda=1, plotit=TRUE, seed=42){
x <- prbs(n, seed, lambda)
t1 <- starthour + ((1:length(x))[abs(c(0,diff(x)))==1] / length(x)) * nhours
## Round to 5 min steps
t1 <- round(t1*12) / 12
val1 <- rep(c(Tmin,Tmax), length(t1))[1:length(t1)]
##
if(plotit){ plot(t1, val1, type="s") }
##
return(data.frame(t=t1, Tset=val1))
}
testseqPM <- function(n, nhours, starthour, seed, lambda=1, plotit=TRUE){
x <- prbs(n, 42, lambda)
ishift <- (1:length(x))[abs(c(0,diff(x)))==1]
X <- rbind(data.frame(ishift=ishift, vals=x[ishift]), data.frame(ishift=ishift, vals=-x[ishift]))
set.seed(seed)
X <- X[sample(1:nrow(X)), ]
t1 <- c(0, cumsum(X$ishift))
val1 <- X$vals
#
t1 <- starthour + t1 / max(t1) * nhours
t1 <- t1[-length(t1)]
##
if(plotit){
plot(t1, val1, type="s")
}
##
return(data.frame(t=t1, Tset=val1))
}
# The collective night changes
# Use n and nhours change the properties (lambda doesn't do anything)
nhours <- c(96,72,72)
val1 <- testseq(n=4, nhours[1], starthour=0, Tmin=19, Tmax=21, seed=87)
## Two days between 17 and 20 degrees C
val2 <- testseq(n=4, nhours[2], starthour=sum(nhours[1]), Tmin=19, Tmax=21, seed=878)
## Two days between 20 and 22 degrees C
val3 <- testseq(n=4, nhours[3], starthour=sum(nhours[1:2]), Tmin=19, Tmax=22, seed=278)
## Combined
X <- rbind(val1, val2, val3)
plot(X$t/24, X$Tset, type="s")
# left interval time points
X$t <- X$t - X$t[1]
# The individual changes
val1 <- testseqPM(n=6, sum(nhours), starthour=0, seed=279)
# Round time to hours
X$t <- round(X$t)
# Round time to 30 min
val1$t <- round(val1$t*2)/2
# remove double time points
val1 <- val1[!duplicated(val1$t), ]
to_equidist_and_lag <- function(val1, i, nlags, X, hoursmax){
tseq <- seq(0, hoursmax, by=0.5)
val1$t <- val1$t - val1$t[1]
val1[nrow(val1)+1, ] <- c(hoursmax, val1$Tset[nrow(val1)])
val1$t <- (val1$t + round(max(val1$t) * (i-1)/nlags)) %% max(val1$t)
# Insert end point
X[nrow(X)+1, ] <- c(hoursmax, X$Tset[nrow(X)])
# Make random delay of start of transient
iseq <- c(-1,-nrow(X))
X$t[iseq] <- X$t[iseq] + round(rexp(nrow(X)-2, rate=2)*2) / 2
# equidistant
X <- approx(X$t, X$Tset, tseq, method="constant")
val1[nrow(val1)+1, ] <- c(hoursmax, val1$Tset[nrow(val1)])
val1 <- approx(val1$t, val1$Tset, tseq, method="constant", rule=2)
return(data.frame(t=X$x, Tset=X$y+val1$y))
}
# How many different room lags
nlags <- 8
setpar("ts", mfrow=c(9,1))
plot(val1$t, val1$Tset, type="s")
Lnight <- list()
for(i in 1:nlags){
tmp <- to_equidist_and_lag(val1, i, nlags, X, sum(nhours))
plot(tmp$t, tmp$Tset, type="s")
Lnight[[i]] <- tmp
}
## ------------------------------------------------------------------------
## ------------------------------------------------------------------------
# Day values
X$Tset <- 22
val1 <- testseq(n=6, sum(nhours), starthour=0, Tmin=0, Tmax=1, seed=972)
# Round time to 30 min
val1$t <- round(val1$t)*2/2
# remove double time points
val1 <- val1[!duplicated(val1$t), ]
#
Lday <- list()
setpar("ts", mfrow=c(9,1))
plot(val1$t, val1$Tset, type="s")
for(i in 1:nlags){
tmp <- to_equidist_and_lag(val1, i, nlags, X, sum(nhours))
plot(tmp$t, tmp$Tset, type="s")
Lday[[i]] <- tmp
}
## ------------------------------------------------------------------------
## ------------------------------------------------------------------------
# Four Weekends
maxhours <- 4*2*48#round(sum(nhours)*2/5+2)
X <- testseq(n=6, maxhours, starthour=0, Tmin=17, Tmax=23, seed=873)
plot(X$t, X$Tset, type="s")
# left interval time points
X$t <- X$t - X$t[1]
X$t <- round(X$t)*2/2
#
val1 <- testseqPM(n=7, max(X$t), starthour=0, seed=972)
plot(val1$t, val1$Tset, type="s")
# Round time to 30 min
val1$t <- round(val1$t)*2/2
# remove double time points
val1 <- val1[!duplicated(val1$t), ]
#
Lweekend <- list()
setpar("ts", mfrow=c(9,1))
plot(val1$t, val1$Tset, type="s")
for(i in 1:nlags){
tmp <- to_equidist_and_lag(val1, i, nlags, X, maxhours)
plot(tmp$t, tmp$Tset, type="s")
Lweekend[[i]] <- tmp
}
## ------------------------------------------------------------------------
## ------------------------------------------------------------------------
# Combine them
i <- 1
Lfinal <- list()
par(mfrow=c(nlags,1))
for(i in 1:nlags){
Xnight <- Lnight[[i]]
Xday <- Lday[[i]]
Xweekend <- Lweekend[[i]]
#
L1 <- split(Xnight, Xnight$t %/% 12)
L2 <- split(Xday, Xday$t %/% 12)
L <- lapply(1:20, function(i){
rbind(L1[[i]], L2[[i]])
})
# Insert the two weekends
L3 <- split(Xweekend, Xweekend$t %/% 48)
Tset <- do.call("rbind", c(L[1:5], L3[1], L[6:10], L3[2], L[11:15], L3[3], L[16:20], L3[4]))[ ,2]
#
# Last Sunday at 17:00
tstartsunday <- trunc(tstart,units="days") - aslt(tstart)$wday*24*3600+17*3600
tseq <- seq(tstartsunday, by=1800, len=28*48)
X <- data.frame(t=tseq, Tset=Tset)
plot(X$t, X$Tset, type="s")
Lfinal[[i]] <- X
}
## ------------------------------------------------------------------------
X <- data.frame(t=X$t)
## Set RegI and RegP
##grep("^[HX_].*_spReg", idtags$VarName, value=TRUE)
## Larger I means less integration! thus its probably reciprocal: see 16. Nov. 2019 stepping.
X$HX_RC_spRegP_building_ltech <- 0.5 #rep(seq(0.3,by=0.1,len=6), each=2)
X$HX_RC_spRegI_building_ltech <- 80 #rep(seq(30,by=10,len=6), each=2)
X$HX_VE_spRegP_building_ltech <- 0.5
X$HX_VE_spRegI_building_ltech <- 80 #rep(seq(30,by=10,len=6), each=2)
## The set points
X[ ,nmTforwardRC] <- 60#rep(c(50,20), 6)
X[ ,nmTforwardVE] <- 45#c(val)#sample(seq(20,65,len=17))#c(65,60,55,50,45,40,35,30,25,30,35,40,45,50,55,60,65))#rep(c(50,20), 6)
## ------------------------------------------------------------------------
## Room set points
for(i in 1:length(nmRooms)){
X[ ,nmRooms[i]] <- Lfinal[[((i-1)%%nlags)+1]]$Tset
}
## ------------------------------------------------------------------------
plot_ts(X, c("RegP","RegI","spTforward"))
plot_ts(X, nmRooms)
Xres <- X
Xres <- Xres[in_range(tstart, Xres$t, tend), ]
plot_ts(Xres, nmRooms)
## Replace the column names with channelid
it <- which(names(Xres) == "t")
names(Xres)[-it] <- unlist(sapply(names(Xres)[-it], function(nm){ idtags$channelid[idtags$VarName == nm] }))
## Write directly in influxdb
# Delete all set points ahead in time
db <- "raw"
nms <- names(Xres)[-1]
for(i in 1:length(nms)){
influx_query(con=con, db=db, query=pst("DELETE FROM observations WHERE channelid = '",nms[i],"' AND direction = 'out' AND time > '", as.character(tstart), "'"), return_xts = FALSE)
}
#
write_values(Xres, idtags, con, colnamesTag = "channelid")
str(Xres)
plot_ts(Xres, patterns="*")
## Write the experiment to the incluxdb
library(onlineforecast)
set.seed(389)
# Start of the experiment
tstart <- ct("2021-03-15 16:00")
library(influxdbr)
## Source the files in the "functions" folder
files <- dir("functions",full.names=TRUE)
for(i in 1:length(files)) source(files[i])
## connection
##con <- influx_connection(host="localhost", port=8086)
con <- influx_connection(host="borgerskolen.centerdenmark.com", port=8086)
##influx_ping(con)
## ## Source the files in the "functions" folder
## files <- dir("functions",full.names=TRUE)
## for(i in 1:length(files)) source(files[i])
##
wd <- "../buildingctrl/main/"
source(paste0(wd,"lib/r/read_idtags.R"))
source(paste0(wd,"lib/r/write_values.R"))
## Find the output channels, to be able to check that the channelids are correct
idtags <- read_idtags(wd)
##
i <- grep("^Room_.*spTemperature__northq", idtags$VarName)
##
(nmRooms <- idtags[i, "VarName"])
(nmTforwardRC <- idtags[grep("^HX_RC_spTforward", idtags$VarName), "VarName"])
(nmTforwardVE <- idtags[grep("^HX_VE_spTforward", idtags$VarName), "VarName"])
##
## ------------------------------------------------------------------------
## Night values
## Generate room experiment, where all rooms are following same Tset, with \pm 1 PRBS making differences between the rooms
nmRooms
testseq <- function(n, nhours, starthour, Tmin, Tmax, lambda=1, plotit=TRUE, seed=42){
x <- prbs(n, seed, lambda)
t1 <- starthour + ((1:length(x))[abs(c(0,diff(x)))==1] / length(x)) * nhours
## Round to 5 min steps
t1 <- round(t1*12) / 12
val1 <- rep(c(Tmin,Tmax), length(t1))[1:length(t1)]
##
if(plotit){ plot(t1, val1, type="s") }
##
return(data.frame(t=t1, Tset=val1))
}
testseqPM <- function(n, nhours, starthour, seed, lambda=1, plotit=TRUE){
x <- prbs(n, 42, lambda)
ishift <- (1:length(x))[abs(c(0,diff(x)))==1]
X <- rbind(data.frame(ishift=ishift, vals=x[ishift]), data.frame(ishift=ishift, vals=-x[ishift]))
set.seed(seed)
X <- X[sample(1:nrow(X)), ]
t1 <- c(0, cumsum(X$ishift))
val1 <- X$vals
#
t1 <- starthour + t1 / max(t1) * nhours
t1 <- t1[-length(t1)]
##
if(plotit){
plot(t1, val1, type="s")
}
##
return(data.frame(t=t1, Tset=val1))
}
# The collective night changes
# Use n and nhours change the properties (lambda doesn't do anything)
nhours <- c(96,72,72)
val1 <- testseq(n=4, nhours[1], starthour=0, Tmin=19, Tmax=21, seed=87)
val2 <- testseq(n=4, nhours[2], starthour=sum(nhours[1]), Tmin=21, Tmax=23, seed=878)
val3 <- testseq(n=4, nhours[3], starthour=sum(nhours[1:2]), Tmin=21, Tmax=22, seed=278)
## Combined
X <- rbind(val1, val2, val3)
plot(X$t/24, X$Tset, type="s")
# left interval time points
X$t <- X$t - X$t[1]
# The individual changes
val1 <- testseqPM(n=6, sum(nhours), starthour=0, seed=279)
# Round time to hours
X$t <- round(X$t)
# Round time to 30 min
val1$t <- round(val1$t*2)/2
# remove double time points
val1 <- val1[!duplicated(val1$t), ]
to_equidist_and_lag <- function(val1, i, nlags, X, hoursmax){
tseq <- seq(0, hoursmax, by=0.5)
val1$t <- val1$t - val1$t[1]
val1[nrow(val1)+1, ] <- c(hoursmax, val1$Tset[nrow(val1)])
val1$t <- (val1$t + round(max(val1$t) * (i-1)/nlags)) %% max(val1$t)
# Insert end point
X[nrow(X)+1, ] <- c(hoursmax, X$Tset[nrow(X)])
# Make random delay of start of transient
iseq <- c(-1,-nrow(X))
X$t[iseq] <- X$t[iseq] + round(rexp(nrow(X)-2, rate=2)*2) / 2
# equidistant
X <- approx(X$t, X$Tset, tseq, method="constant")
val1[nrow(val1)+1, ] <- c(hoursmax, val1$Tset[nrow(val1)])
val1 <- approx(val1$t, val1$Tset, tseq, method="constant", rule=2)
return(data.frame(t=X$x, Tset=X$y+val1$y))
}
# How many different room lags
nlags <- 8
setpar("ts", mfrow=c(9,1))
plot(val1$t, val1$Tset, type="s")
Lnight <- list()
for(i in 1:nlags){
tmp <- to_equidist_and_lag(val1, i, nlags, X, sum(nhours))
plot(tmp$t, tmp$Tset, type="s")
Lnight[[i]] <- tmp
}
## ------------------------------------------------------------------------
## ------------------------------------------------------------------------
# Day values
X$Tset <- 22
val1 <- testseq(n=6, sum(nhours), starthour=0, Tmin=0, Tmax=1, seed=972)
# Round time to 30 min
val1$t <- round(val1$t)*2/2
# remove double time points
val1 <- val1[!duplicated(val1$t), ]
#
Lday <- list()
setpar("ts", mfrow=c(9,1))
plot(val1$t, val1$Tset, type="s")
for(i in 1:nlags){
tmp <- to_equidist_and_lag(val1, i, nlags, X, sum(nhours))
plot(tmp$t, tmp$Tset, type="s")
Lday[[i]] <- tmp
}
## ------------------------------------------------------------------------
## ------------------------------------------------------------------------
# Four Weekends
maxhours <- 4*2*48#round(sum(nhours)*2/5+2)
X <- testseq(n=6, maxhours, starthour=0, Tmin=17, Tmax=23, seed=873)
plot(X$t, X$Tset, type="s")
# left interval time points
X$t <- X$t - X$t[1]
X$t <- round(X$t)*2/2
#
val1 <- testseqPM(n=7, max(X$t), starthour=0, seed=972)
plot(val1$t, val1$Tset, type="s")
# Round time to 30 min
val1$t <- round(val1$t)*2/2
# remove double time points
val1 <- val1[!duplicated(val1$t), ]
#
Lweekend <- list()
setpar("ts", mfrow=c(9,1))
plot(val1$t, val1$Tset, type="s")
for(i in 1:nlags){
tmp <- to_equidist_and_lag(val1, i, nlags, X, maxhours)
plot(tmp$t, tmp$Tset, type="s")
Lweekend[[i]] <- tmp
}
## ------------------------------------------------------------------------
## ------------------------------------------------------------------------
# Combine them
i <- 1
Lfinal <- list()
par(mfrow=c(nlags,1))
for(i in 1:nlags){
Xnight <- Lnight[[i]]
Xday <- Lday[[i]]
Xweekend <- Lweekend[[i]]
#
L1 <- split(Xnight, Xnight$t %/% 12)
L2 <- split(Xday, Xday$t %/% 12)
L <- lapply(1:20, function(i){
rbind(L1[[i]], L2[[i]])
})
# Insert the two weekends
L3 <- split(Xweekend, Xweekend$t %/% 48)
Tset <- do.call("rbind", c(L[1:5], L3[1], L[6:10], L3[2], L[11:15], L3[3], L[16:20], L3[4]))[ ,2]
#
# Last Sunday at 17:00
tstartsunday <- trunc(tstart,units="days") - aslt(tstart)$wday*24*3600+17*3600
tseq <- seq(tstartsunday, by=1800, len=28*48)
X <- data.frame(t=tseq, Tset=Tset)
plot(X$t, X$Tset, type="s")
Lfinal[[i]] <- X
}
## ------------------------------------------------------------------------
X <- data.frame(t=X$t)
## Set RegI and RegP
##grep("^[HX_].*_spReg", idtags$VarName, value=TRUE)
## Larger I means less integration! thus its probably reciprocal: see 16. Nov. 2019 stepping.
X$HX_RC_spRegP_building_ltech <- 0.5 #rep(seq(0.3,by=0.1,len=6), each=2)
X$HX_RC_spRegI_building_ltech <- 80 #rep(seq(30,by=10,len=6), each=2)
X$HX_VE_spRegP_building_ltech <- 0.5
X$HX_VE_spRegI_building_ltech <- 80 #rep(seq(30,by=10,len=6), each=2)
## The set points
X[ ,nmTforwardRC] <- 60#rep(c(50,20), 6)
X[ ,nmTforwardVE] <- 45#c(val)#sample(seq(20,65,len=17))#c(65,60,55,50,45,40,35,30,25,30,35,40,45,50,55,60,65))#rep(c(50,20), 6)
## ------------------------------------------------------------------------
## Room set points
for(i in 1:length(nmRooms)){
X[ ,nmRooms[i]] <- Lfinal[[((i-1)%%nlags)+1]]$Tset
}
## ------------------------------------------------------------------------
plot_ts(X, c("RegP","RegI","spTforward"))
plot_ts(X, nmRooms)
Xres <- X
Xres <- Xres[in_range(tstart, Xres$t), ]
plot_ts(Xres, nmRooms)
## Replace the column names with channelid
it <- which(names(Xres) == "t")
names(Xres)[-it] <- unlist(sapply(names(Xres)[-it], function(nm){ idtags$channelid[idtags$VarName == nm] }))
## Write directly in influxdb
# Delete all set points ahead in time
db <- "raw"
nms <- names(Xres)[-1]
for(i in 1:length(nms)){
influx_query(con=con, db=db, query=pst("DELETE FROM observations WHERE channelid = '",nms[i],"' AND direction = 'out' AND time > '", as.character(tstart), "'"), return_xts = FALSE)
}
#
write_values(Xres, idtags, con, colnamesTag = "channelid")
str(Xres)
plot_ts(Xres, patterns="*")
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment