Converter from RH to q using standard atmosheric pressure and altitude

metno · Nov 5, 2024 · 09e9afc · 09e9afc
1 parent 94ae3c5
commit 09e9afc
Show file tree

Hide file tree

Showing 3 changed files with 26 additions and 5 deletions.
diff --git a/pysurfex/cfg/config.yml b/pysurfex/cfg/config.yml
@@ -251,7 +251,18 @@ QA:
           mslp:
             name: air_pressure_at_sea_level
             units: "Pa"
-
+        rh2q_z:
+          rh:
+            name: relative_humidity_2m
+            level: 2
+            units: "%"
+          t:
+            name: air_temperature_2m
+            level: 2
+            units: "K"
+          altitude:
+            name: altitude
+            units: "m"
     grib1:
       converter:
         none:

diff --git a/pysurfex/cmd_parsing.py b/pysurfex/cmd_parsing.py
@@ -176,7 +176,7 @@ def parse_args_create_forcing(argv):
         type=str,
         help="Converter function to specific humidity",
         default="none",
-        choices=["none", "rh2q", "rh2q_mslp"],
+        choices=["none", "rh2q", "rh2q_mslp", "rh2q_z"],
     )
 
     group_ps = parser.add_argument_group("PS", description="Surface air pressure [Pa]")

diff --git a/pysurfex/read.py b/pysurfex/read.py
@@ -186,6 +186,12 @@ def __init__(self, name, initial_time, defs, conf, fileformat):
                 fileformat, defs, conf[self.name]["altitude"]
             )
             self.pres = self.create_variable(fileformat, defs, conf[self.name]["mslp"])
+        elif name == "rh2q_z":
+            self.r_h = self.create_variable(fileformat, defs, conf[self.name]["rh"])
+            self.temp = self.create_variable(fileformat, defs, conf[self.name]["t"])
+            self.altitude = self.create_variable(
+                fileformat, defs, conf[self.name]["altitude"]
+            )
         elif name == "windspeed" or name == "winddir":
             self.x_wind = self.create_variable(fileformat, defs, conf[self.name]["x"])
             self.y_wind = self.create_variable(fileformat, defs, conf[self.name]["y"])
@@ -336,7 +342,7 @@ def read_time_step(self, geo, validtime, cache):
             elif self.name == "winddir":
                 field = np.mod(90 - np.rad2deg(np.arctan2(field_y, field_x)), 360)
 
-        elif self.name == "rh2q" or self.name == "rh2q_mslp":
+        elif self.name == "rh2q" or self.name == "rh2q_mslp" or self.name == "rh2q_z":
             """
             ZES = 6.112 * exp((17.67 * (ZT - 273.15)) / ((ZT - 273.15) + 243.5))
             ZE = ZRH * ZES
@@ -345,8 +351,12 @@ def read_time_step(self, geo, validtime, cache):
             """
             field_r_h = self.r_h.read_variable(geo, validtime, cache)  # %
             field_temp = self.temp.read_variable(geo, validtime, cache)  # In K
-            field_pres = self.pres.read_variable(geo, validtime, cache)  # In Pa
-            if self.name == "rh2q_mslp":
+            if self.name == "rh2q_z":
+                field_pres = np.array(field_r_h.shape)
+                field_pres.fill(101325.0)
+            else:
+                field_pres = self.pres.read_variable(geo, validtime, cache)  # In Pa
+            if self.name == "rh2q_mslp" or self.name == "rh2q_z":
                 field_altitude = self.altitude.read_variable(
                     geo, validtime, cache
                 )  # In m