Possible bug in thermo.relative_humidity_from_specific_humidity

[mpvginde]

What happened?

When starting from a given pressure, specific humidity and temperature, using earthkit.meteo.thermo functions to calculate specific humidity --> dewpoint temperature --> relative humidity and specific humidity --> relaitve humidity results in different results for relative humidity.

What are the steps to reproduce the bug?

from earthkit.meteo import thermo

p = 101325
q = 0.00035
t = 247.5
d = thermo.dewpoint_from_specific_humidity(q=q,p=p)

rh_from_dp = thermo.relative_humidity_from_dewpoint(t=t,td=d)
rh_from_q = thermo.relative_humidity_from_specific_humidity(t=t,q=q,p=p)
print(f"Relative humidity from dewpoint: {rh_from_dp}%")
print(f"Relative humidity from specific humidty: {rh_from_q}%")

dp_from_rh_from_dp = thermo.dewpoint_from_relative_humidity(t=t,r=rh_from_dp)
dp_from_rh_from_q = thermo.dewpoint_from_relative_humidity(t=t,r=rh_from_q)
print(f"Original dewpoint: {d} K")
print(f"Dewpoint from RH from dewpoint: {dp_from_rh_from_dp} K")
print(f"Dewpoint from RH from specific humidity: {dp_from_rh_from_q} K")

### Version

0.1.1.dev

### Platform (OS and architecture)

ATOS

### Relevant log output

_No response_

### Accompanying data

_No response_

### Organisation

RMIB

[sandorkertesz]

@mpvginde, thank your for reporting this issue.

My reply below is not a solution rather an explanation.

This happens when the temperature is < 0C. If you try e.g. t = 277.5 your example will work. This is down to the saturation vapour pressure formula used in the computations, see:

https://earthkit-meteo.readthedocs.io/en/latest/_api/meteo/thermo/array/thermo/index.html#meteo.thermo.array.thermo.saturation_vapour_pressure

By default it uses the phase="mixed" option (the other values are "water" and "ice"), which for t > 0C is the same as phase="water". Some of the methods allow specifying the phase keyword, but not the relative/specific humidity ones, which always use the default saturation vapour pressure.

The dewpoint is computed by inverting this formula to find the temperature from the saturation vapour pressure. However, this is always done with the phase="water" option. So this results in the discrepancy for t < 0C values.

We will review these computations and come up with a solution.

[mpvginde]

Thanks @sandorkertesz, For the explanation.

So if I understand correct, in case of t <= 0C it is best to go from specific hum. to dewpoint to relative hum. As this will take into account the mixed phase.