Model B discontinuity at 50% Humidity
Opened this issue · 1 comments
Dear all
here my example
ex<-data.table::CJ(Temperature=seq(-20,40,5),Humidity=seq(0,100,0.25) )
ex<-rbind(cbind(ex,Dewpoint=calcDewPoint(ex$Humidity, ex$Temperature,mode="A"),Model=c(rep("A",dim(ex)[1]))),
cbind(ex,Dewpoint=calcDewPoint(ex$Humidity, ex$Temperature,mode="B"),Model=c(rep("B",dim(ex)[1]))),
cbind(ex,Dewpoint=calcDewPoint(ex$Humidity, ex$Temperature,mode="C"),Model=c(rep("C",dim(ex)[1]))))
ex %>%
ggplot( aes(x=Humidity, y=Dewpoint, group=Temperature, color=Temperature)) +
geom_line() +
viridis::scale_color_viridis(discrete = FALSE, option="magma") +
theme(
legend.position="none",
plot.title = element_text(size=14)
) +
facet_wrap(~Model)
I've found a temporary solutions in the meanwhile @anadiedrichs will fix the package.
All credits belongs to @anadiedrichs , i've just give my 2 cents on top of her huge work.
There are two main problems:
- there was a mistake in 'else" equation where a temperature need to be considered as absolute temperature
- the two equations are not complementary but had a field of calculation that was overlapping 40% <= RH<= 100% and 0° < t < 30°C
So i've slightly modified the algorithm introducing those limitations and introducing the possibility of having the full range formula using equation 8
#Workaround calcDewPoint.B MOD
/ ```r
calcDewPoint.Bmod <- function(RH,temp, equation) {
dw <- 0
dw <- mapply(calculation, RH, temp, equation )
return(dw)
}
calculation <-function(RH,temp,equation)
{
dw <- 0
if(RH >= 40 & equation=="linear" ){dw <- (0.198 + 0.0017*temp) * RH + (0.84*temp) - 19.2}
else if(RH >= 40 & equation=="quadratic"){dw <- temp - ((100-RH)/5) * ((temp+273.15)/300)^2 -0.00135 * (RH - 84)^2 + 0.35}
else if(temp >= -40 & temp <= 50 & equation=="Log"){dw <- 243.04 * ( log(RH/100) + (17.625 * temp) / (243.04 + temp) ) / ( 17.625 - log(RH/100) - (17.625 * temp) / (243.04 + temp) )}
else dw<-NA
return(dw)
}
/ ```r
The new calcDewPoint.Bmod function documentation could be
Usage
calcDewPoint.Bmod(RH, temp, equation)
Arguments
RH [in percentage] relative humidity, an integer or double value between 0 and 100.
temp [°C] environmental temperature, an integer or double value between -20 and 60 °C
equation to choose within the following three options and associated limitations (out of those limits an NA will be produced:
linear Eq (20) page 230 field of calculation -> 40% <= RH<= 100% and 0° < t < 30°C
quadratic Eq (21) page 230 field of calculation -> 40% <= RH <= 100% and 0° < t <30°C
Log Eq 8 page 226 with Alduchov and Eskridge (1996) coefficents. These provide values for es with a relative error of < 0.4% over the range -40°C <= t <= 50°C
Eq 8 page 226 with Alduchov and Eskridge (1996) coefficents. These provide values for es with a relative error of < 0.4% over the range -40°C <= t <= 50°C
Eq (8) B1 * ( log(RH/100) + (A1 * temp) / (B1 + temp) ) / ( A1 - log(RH/100) - (A1 * temp) / (B1 + temp) )
243.04 * ( log(RH/100) + (17.625 * temp) / (243.04 + temp) ) / ( 17.625 - log(RH/100) - (17.625 * temp) / (243.04 + temp) )
A1= 17.625
B1 = 243.04 °C
C1= 610.94 Pa
Value
dew point value (double)
Here below the results
Code to generate the graph
ex<-data.table::CJ(Temperature=seq(-20,40,5),Humidity=seq(0,100,0.25) )
ex<-rbind(cbind(ex,Dewpoint=calcDewPoint(ex$Humidity, ex$Temperature,mode="A"),Model=c(rep("A",dim(ex)[1]))),
cbind(ex,Dewpoint=calcDewPoint(ex$Humidity, ex$Temperature,mode="B"),Model=c(rep("B_Original",dim(ex)[1]))),
cbind(ex,Dewpoint=calcDewPoint.Bmod(ex$Humidity, ex$Temperature, equation="linear"),Model=c(rep("BmodLinear",dim(ex)[1]))),
cbind(ex,Dewpoint=calcDewPoint.Bmod(ex$Humidity, ex$Temperature, equation="quadratic"),Model=c(rep("BmodQuadratic",dim(ex)[1]))),
cbind(ex,Dewpoint=calcDewPoint.Bmod(ex$Humidity, ex$Temperature, equation="Log"),Model=c(rep("BmodLog",dim(ex)[1]))),
cbind(ex,Dewpoint=calcDewPoint(ex$Humidity, ex$Temperature,mode="C"),Model=c(rep("C",dim(ex)[1]))))
ex %>%
ggplot( aes(x=Humidity, y=Dewpoint, group=Temperature, color=Temperature)) +
geom_line() +
viridis::scale_color_viridis(discrete = FALSE, option="magma") +
theme(
legend.position="none",
plot.title = element_text(size=14)
) +
facet_wrap(~Model)
in order to see deviations from humidity.to.dewpoint from weathermetrics package powered by Brooke Anderson https://github.com/geanders/weathermetrics/
that is based on equations from the source code for the US National Weather Service's online heat index calculator.
ex<-data.table::CJ(Temperature=seq(-20,60,1),Humidity=seq(0,100,5) )
ex<-rbind(cbind(ex,Dewpoint=calcDewPoint(ex$Humidity, ex$Temperature,mode="A")-weathermetrics::humidity.to.dewpoint(rh = ex$Humidity, t = ex$Temperature, temperature.metric = "celsius"),Model=c(rep("A",dim(ex)[1]))),
cbind(ex,Dewpoint=calcDewPoint(ex$Humidity, ex$Temperature,mode="B")-weathermetrics::humidity.to.dewpoint(rh = ex$Humidity, t = ex$Temperature, temperature.metric = "celsius"),Model=c(rep("B_Original",dim(ex)[1]))),
cbind(ex,Dewpoint=calcDewPoint.Bmod(ex$Humidity, ex$Temperature, equation="linear")-weathermetrics::humidity.to.dewpoint(rh = ex$Humidity, t = ex$Temperature, temperature.metric = "celsius"),Model=c(rep("BmodLinear",dim(ex)[1]))),
cbind(ex,Dewpoint=calcDewPoint.Bmod(ex$Humidity, ex$Temperature, equation="quadratic")-weathermetrics::humidity.to.dewpoint(rh = ex$Humidity, t = ex$Temperature, temperature.metric = "celsius"),Model=c(rep("BmodQuadratic",dim(ex)[1]))),
cbind(ex,Dewpoint=calcDewPoint.Bmod(ex$Humidity, ex$Temperature, equation="Log")-weathermetrics::humidity.to.dewpoint(rh = ex$Humidity, t = ex$Temperature, temperature.metric = "celsius"),Model=c(rep("BmodLog",dim(ex)[1]))) ,
cbind(ex,Dewpoint=calcDewPoint(ex$Humidity, ex$Temperature,mode="C")-weathermetrics::humidity.to.dewpoint(rh = ex$Humidity, t = ex$Temperature, temperature.metric = "celsius"),Model=c(rep("C",dim(ex)[1]))))
ex %>%
filter(Model!="B_Original" & Model!="C") %>%
ggplot( aes(x=Humidity, y=Dewpoint, group=Temperature, color=Temperature)) +
geom_line() +
viridis::scale_color_viridis(discrete = FALSE, option="magma") +
theme(
legend.position="none",
plot.title = element_text(size=14)
) +
facet_wrap(~Model)
ex %>%
filter(Model=="B_Original" ) %>%
ggplot( aes(x=Humidity, y=Dewpoint, group=Temperature, color=Temperature)) +
geom_line() +
viridis::scale_color_viridis(discrete = FALSE, option="magma") +
theme(
legend.position="none",
plot.title = element_text(size=14)
) +
facet_wrap(~Model)
ex %>%
filter(Model=="C" ) %>%
ggplot( aes(x=Humidity, y=Dewpoint, group=Temperature, color=Temperature)) +
geom_line() +
viridis::scale_color_viridis(discrete = FALSE, option="magma") +
theme(
legend.position="none",
plot.title = element_text(size=14)
) +
facet_wrap(~Model)