From f51823686ed9d4bec0bbb1498b8c9a1d55a76ec5 Mon Sep 17 00:00:00 2001
From: sbiri <sbiri@noc.ac.uk>
Date: Fri, 27 Aug 2021 09:04:41 +0100
Subject: [PATCH] Do not calculate lhf if a measure of humidity is not input
---
AirSeaFluxCode.py | 28 +++++++++++++++++++---------
get_init.py | 2 +-
hum_subs.py | 5 ++---
3 files changed, 22 insertions(+), 13 deletions(-)
diff --git a/AirSeaFluxCode.py b/AirSeaFluxCode.py
index ace0647..671982c 100644
--- a/AirSeaFluxCode.py
+++ b/AirSeaFluxCode.py
@@ -128,9 +128,9 @@ def AirSeaFluxCode(spd, T, SST, lat=None, hum=None, P=None, hin=18, hout=10,
39. relative humidity (rh)
40. thickness of the viscous layer (delta)
41. lv latent heat of vaporization (Jkg−1)
- 40. flag ("n": normal, "o": out of nominal range,
+ 42. flag ("n": normal, "o": out of nominal range,
"u": u10n<0, "q":q10n<0
- "m": missing,
+ "m": missing,
"l": Rib<-0.5 or Rib>0.2 or z/L>1000,
"r" : rh>100%,
"i": convergence fail at n)
@@ -180,9 +180,7 @@ def AirSeaFluxCode(spd, T, SST, lat=None, hum=None, P=None, hin=18, hout=10,
np.round(np.nanmedian(qair), 7))
if (np.all(np.isnan(qsea)) or np.all(np.isnan(qair))):
print("qsea and qair cannot be nan")
- if (hum == None):
- rh = np.ones(sst.shape)*80
- elif (hum[0] == 'rh'):
+ if ((hum[0] == 'rh') or (hum[0] == 'no')):
rh = hum[1]
elif (hum[0] == 'Td'):
Td = hum[1] # dew point temperature (K)
@@ -193,8 +191,11 @@ def AirSeaFluxCode(spd, T, SST, lat=None, hum=None, P=None, hin=18, hout=10,
rh = 100*esd/es
flag = np.empty(spd.shape, dtype="object")
flag[:] = "n"
- flag = np.where(np.isnan(spd+T+SST+hum[1]+P+Rs+Rl), "m", flag)
- flag = np.where(rh > 100, "r", flag)
+ if (hum[0] == 'no'):
+ flag = np.where(np.isnan(spd+T+SST+P+Rs+Rl), "m", flag)
+ else:
+ flag = np.where(np.isnan(spd+T+SST+hum[1]+P+Rs+Rl), "m", flag)
+ flag = np.where(rh > 100, "r", flag)
dt = Ta - sst
dq = qair - qsea
@@ -530,6 +531,14 @@ def AirSeaFluxCode(spd, T, SST, lat=None, hum=None, P=None, hin=18, hout=10,
if (((cskin == 0) and (wl == 0)) and
(np.all(Rl == 370) and np.all(Rs == 150))):
Rl, Rs, Rnl = Rl*np.nan, Rs*np.nan, Rnl*np.nan
+ # Do not calculate lhf if a measure of humidity is not input
+ if (hum[0] == 'no'):
+ latent = np.ones(sst.shape)*np.nan
+ qsr = np.ones(sst.shape)*np.nan
+ q10n = np.ones(sst.shape)*np.nan
+ qref = np.ones(sst.shape)*np.nan
+ qair = np.ones(sst.shape)*np.nan
+ rh = np.ones(sst.shape)*np.nan
res = np.zeros((41, len(spd)))
res[0][:] = tau
@@ -577,8 +586,9 @@ def AirSeaFluxCode(spd, T, SST, lat=None, hum=None, P=None, hin=18, hout=10,
if (out == 0):
res[:, ind] = np.nan
# set missing values where data have non acceptable values
- res = np.asarray([np.where(q10n < 0, np.nan,
- res[i][:]) for i in range(41)])
+ if (hum[0] != 'no'):
+ res = np.asarray([np.where(q10n < 0, np.nan,
+ res[i][:]) for i in range(41)])
res = np.asarray([np.where(u10n < 0, np.nan,
res[i][:]) for i in range(41)])
# output with pandas
diff --git a/get_init.py b/get_init.py
index 144714c..7d2170b 100644
--- a/get_init.py
+++ b/get_init.py
@@ -108,7 +108,7 @@ def get_init(spd, T, SST, lat, hum, P, Rl, Rs, cskin, skin, wl, gust, L, tol,
lat = np.ones(spd.shape)*np.copy(lat)
if (hum == None):
RH = np.ones(SST.shape)*80
- hum = ['rh', RH]
+ hum = ['no', RH]
else:
hum = hum
if ((np.all(P == None)) or np.all(np.isnan(P))):
diff --git a/hum_subs.py b/hum_subs.py
index 9eff6a4..76777bd 100644
--- a/hum_subs.py
+++ b/hum_subs.py
@@ -372,15 +372,14 @@ def get_hum(hum, T, sst, P, qmeth):
specific humidity over sea surface
"""
- if (hum[0] not in ['rh', 'q', 'Td']):
+ if (hum[0] not in ['rh', 'q', 'Td', 'no']):
sys.exit("unknown humidity input")
qair, qsea = np.nan, np.nan
- elif (hum[0] == 'rh'):
+ elif ((hum[0] == 'rh') or (hum[0] == 'no')):
RH = hum[1]
if (np.all(RH < 1)):
sys.exit("input relative humidity units should be \%")
qair, qsea = np.nan, np.nan
- # RH = np.where(RH > 100, np.nan, RH) # ensure RH <=100
qsea = qsat_sea(sst, P, qmeth)/1000 # surface water q (kg/kg)
qair = qsat_air(T, P, RH, qmeth)/1000 # q of air (kg/kg)
elif (hum[0] == 'q'):
--
GitLab