Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Model B discontinuity at 50% Humidity #68

Open
ScymnusRIP opened this issue Mar 17, 2023 · 1 comment
Open

Model B discontinuity at 50% Humidity #68

ScymnusRIP opened this issue Mar 17, 2023 · 1 comment

Comments

@ScymnusRIP
Copy link

ScymnusRIP commented Mar 17, 2023

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)

What's happen to model B?
image

@ScymnusRIP ScymnusRIP changed the title Model B discontinuity at 50°C Model B discontinuity at 50% Humidity Mar 17, 2023
@ScymnusRIP
Copy link
Author

ScymnusRIP commented Mar 18, 2023

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:

  1. there was a mistake in 'else" equation where a temperature need to be considered as absolute temperature
  2. 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)


image

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)

image

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)

image

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)

image

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

1 participant