library(rPython)
## Warning: package 'rPython' was built under R version 3.1.3
## Loading required package: RJSONIO
python.exec( "
import math
def comfPMV(ta, tr, vel, rh, met, clo, wme):
#returns [pmv, ppd]
#ta, air temperature (C)
#tr, mean radiant temperature (C)
#vel, relative air velocity (m/s)
#rh, relative humidity (%) Used only this way to input humidity level
#met, metabolic rate (met)
#clo, clothing (clo)
#wme, external work, normally around 0 (met)
pa = rh * 10 * math.exp(16.6536 - 4030.183 / (ta + 235))
icl = 0.155 * clo #thermal insulation of the clothing in M2K/W
m = met * 58.15 #metabolic rate in W/M2
w = wme * 58.15 #external work in W/M2
mw = m - w #internal heat production in the human body
if (icl <= 0.078):
fcl = 1 + (1.29 * icl)
else:
fcl = 1.05 + (0.645 * icl)
#heat transf. coeff. by forced convection
hcf = 12.1 * math.sqrt(vel)
taa = ta + 273
tra = tr + 273
tcla = taa + (35.5 - ta) / (3.5 * icl + 0.1)
p1 = icl * fcl
p2 = p1 * 3.96
p3 = p1 * 100
p4 = p1 * taa
p5 = (308.7 - 0.028 * mw) + (p2 * math.pow(tra / 100, 4))
xn = tcla / 100
xf = tcla / 50
eps = 0.00015
n = 0
while abs(xn - xf) > eps:
xf = (xf + xn) / 2
hcn = 2.38 * math.pow(abs(100.0 * xf - taa), 0.25)
if (hcf > hcn):
hc = hcf
else:
hc = hcn
xn = (p5 + p4 * hc - p2 * math.pow(xf, 4)) / (100 + p3 * hc)
n += 1
if (n > 150):
print 'Max iterations exceeded'
return 1
tcl = 100 * xn - 273
#heat loss diff. through skin
hl1 = 3.05 * 0.001 * (5733 - (6.99 * mw) - pa)
#heat loss by sweating
if mw > 58.15:
hl2 = 0.42 * (mw - 58.15)
else:
hl2 = 0
#latent respiration heat loss
hl3 = 1.7 * 0.00001 * m * (5867 - pa)
#dry respiration heat loss
hl4 = 0.0014 * m * (34 - ta)
#heat loss by radiation
hl5 = 3.96 * fcl * (math.pow(xn, 4) - math.pow(tra / 100, 4))
#heat loss by convection
hl6 = fcl * hc * (tcl - ta)
ts = 0.303 * math.exp(-0.036 * m) + 0.028
pmv = ts * (mw - hl1 - hl2 - hl3 - hl4 - hl5 - hl6)
ppd = 100.0 - 95.0 * math.exp(-0.03353 * pow(pmv, 4.0)
- 0.2179 * pow(pmv, 2.0))
return pmv
" )
library(data.table)
## Warning: package 'data.table' was built under R version 3.1.3
sensors <- fread("~/Downloads/data-for-carlos/sensor_info2.csv")
samples <- fread("~/Downloads/data-for-carlos/all_sensor_data.csv",showProgress=FALSE)
library(lubridate)
## Warning: package 'lubridate' was built under R version 3.1.3
##
## Attaching package: 'lubridate'
##
## The following objects are masked from 'package:data.table':
##
## hour, mday, month, quarter, wday, week, yday, year
samples <- samples[,.(Date=ymd_hms(datetime),sensor_id,reading)]
dt <- merge(sensors,samples,by="sensor_id")
west <- dt[building=="KWWest"]
east <- dt[building=="KWEast"]
illima <- dt[building=="Illima"]
west <- dcast(west,Date ~ sensor_id,value.var=c("reading"))
west <- west[,c("Date","4","3"),with=FALSE]
west <- west[complete.cases(west)]
colnames(west) <- c("Date","temp","relh")
east <- dcast(east,Date ~ sensor_id,value.var=c("reading"))
illima <- dcast(illima,Date ~ sensor_id,value.var=c("reading"))
#python.call("comfPMV",36.6,36.6,0.1,75,1.1,0.5,0)
#ta <- west[,"4",with=FALSE]
f2c <- function(tempF){
tempC <- 5.0 / 9.0 * (tempF - 32)
return(tempC)
}
ta <- f2c(west[,temp])
tr <- ta
#vel <- 0.1
vel <- array(0.1 ,nrow(west))
#rh <- west[,"3",with=FALSE]
rh <- west[,relh]
#met <- 1.1
met <- array(1.1,nrow(west))
#clo <- 0.5
clo <- array(0.5,nrow(west))
#wme <- 0
wme <- array(0,nrow(west))
a <- array(NA,nrow(west))
for (i in 1:nrow(west)){
a[i] <- python.call("comfPMV",ta[i],tr[i],vel[i],rh[i],met[i],clo[i],wme[i])
}
west$pmv <- a
west$iscomfortable <- ifelse(west$pmv > -0.5 & west$pmv < 0.5,"comfortable","uncomfortable")
library(dygraphs)
## Warning: package 'dygraphs' was built under R version 3.1.3
dygraph(west[,.(Date,pmv)], main = "Comfort Level") %>% dyRangeSelector()