From 36e23fff83547e364daea22868d4895b5e0923b8 Mon Sep 17 00:00:00 2001 From: sbiri <sbiri@noc.ac.uk> Date: Tue, 21 Jul 2020 12:34:44 +0100 Subject: [PATCH] will be in hum_subs --- VaporPressure.py | 345 ----------------------------------------------- 1 file changed, 345 deletions(-) delete mode 100644 VaporPressure.py diff --git a/VaporPressure.py b/VaporPressure.py deleted file mode 100644 index 366e588..0000000 --- a/VaporPressure.py +++ /dev/null @@ -1,345 +0,0 @@ -import numpy as np -import logging -""" - Calculate the saturation vapor pressure - - Example: Calculate the vapor pressure over liquid water using the WMO - formula - temp = FINDGEN(100) - 100. - P = VaporPressure(temp, 'liquid', 'WMO') - - For temperatures above 0 deg C the vapor pressure over liquid water - is calculated. - - The optional parameter 'liquid' changes the calculation to vapor pressure - over liquid water over the entire temperature range. - The formula to be used can be selected with the appropriate keyword - - The current default fomulas are Hyland and Wexler for liquid and - Goff Gratch for ice. - - Parameters - ---------- - temp : float - Temperature in [deg C] - phase : str - 'liquid' : Calculate vapor pressure over liqiud water or - 'ice' : Calculate vapor pressure over ice - meth : str - formula to be used - MartiMauersberger : vaporpressure formula from Marti Mauersberger - MagnusTetens : vaporpressure formula from Magnus Tetens - GoffGratch : vaporpressure formula from Goff Gratch - Buck_original : vaporpressure formula from Buck (original) - Buck_manual : vaporpressure formula from the Buck manual - CIMO : vaporpressure formula recommended by CIMO - Vaisala : vaporpressure formula from Vaisala - WMO_Goff : vaporpressure formula from WMO, which should have been Goff - WMO2000 : vaporpressure formula from WMO (2000) containing a typo - Wexler : vaporpressure formula from Wexler (1976) - Sonntag : vaporpressure formula from Sonntag (1994) - Bolton : vaporpressure formula from Bolton (1980) - HylandWexler : vaporpressure formula from Hyland and Wexler (1983) - IAPWS : vaporpressure formula from IAPWS (2002) - Preining : vaporpressure formula from Preining (2002) - MurphyKoop : vaporpressure formula from Murphy and Koop (2005) - - Returns - ------- - P : float - Saturation vapor pressure in [hPa] - - - Author - ------ - Holger Voemel - NCAR / EOL - PO BOX 3000 - Boulder, CO 80303 - USA - Ported to Python and modified by S. Biri - - COPYRIGHT - --------- - Copyright (c) 2015, Holger Voemel -""" -def VaporPressure(temp, P, phase, meth): - logging.basicConfig(filename='VaporPressure.log', - format='%(asctime)s %(message)s',level=logging.info) - - Psat = np.zeros(temp.size)*np.nan - if (np.nanmin(temp) > 200): # if Ta in Kelvin convert to Celsius - temp = temp-273.16 - T = np.copy(temp)+273.16 # Most formulas use T in [K] - # Formulas using [C] use the variable temp - # Default uses Hyland and Wexler over liquid. While this may not be - # the best formula, it is consistent with what Vaisala uses in their - # system - # Calculate saturation pressure over liquid water - if (phase == 'liquid'): - if (meth == 'HylandWexler' or meth == ''): - logging.info("Source Hyland, R. W. and A. Wexler, Formulations \ - for the Thermodynamic Properties of the saturated \ - Phases of H2O from 173.15K to 473.15K, \ - ASHRAE Trans, 89(2A), 500-519, 1983.") - Psat = np.exp(-0.58002206e4/T+0.13914993e1-0.48640239e-1*T + - 0.41764768e-4*np.power(T, 2) - - 0.14452093e-7*np.power(T, 3) + - 0.65459673e1*np.log(T))/100 - if (meth == 'Hardy'): - logging.info("Source Hardy, B., 1998, ITS-90 Formulations \ - for Vapor Pressure, Frostpoint Temperature, \ - Dewpoint Temperature, and Enhancement Factors \ - in the Range -100 to +100°C, The Proceedings of \ - the Third International Symposium on Humidity & \ - Moisture, London, England") - Psat = np.exp(-2.8365744e3/np.power(T, 2)-6.028076559e3/T + - 1.954263612e1-2.737830188e-2*T + - 1.6261698e-5*np.power(T, 2) + - 7.0229056e-10*np.power(T, 3) - - 1.8680009e-13*np.power(T, 4) + - 2.7150305*np.log(T))/100 - if (meth == 'Preining'): - logging.info("Source : Vehkamaeki, H., M. Kulmala, I. Napari, \ - K. E. J. Lehtinen, C.Timmreck, M. Noppel, and \ - A. Laaksonen (2002), J. Geophys. Res., 107, \ - doi:10.1029/2002JD002184.") - Psat = np.exp(-7235.424651/T+77.34491296+5.7113e-3*T - - 8.2*np.log(T))/100 - if (meth == 'Wexler'): - logging.info("Wexler, A., Vapor Pressure Formulation for Water \ - in Range 0 to 100 C. A Revision, Journal of Research \ - of the National Bureau of Standards - A. Physics and \ - Chemistry, September - December 1976, Vol. 80A, \ - Nos.5 and 6, 775-785") - Psat = np.exp(-0.29912729e4*np.power(T, -2) - - 0.60170128e4*np.power(T, -1) + - 0.1887643854e2-0.28354721e-1*T + - 0.17838301e-4*np.power(T, 2) - - 0.84150417e-9*np.power(T, 3) + - 0.44412543e-12*np.power(T, 4) + # This line was corrected from '-' to '+' following the original citation. (HV 20140819). The change makes only negligible difference - 2.858487*np.log(T))/100 - if ((meth == 'GoffGratch') or (meth == 'MartiMauersberger')): - # logging.info("Marti and Mauersberger don't have a vapor pressure \ - # curve over liquid. Using Goff Gratch instead; \ - # Goff Gratch formulation Source : Smithsonian \ - # Meteorological Tables, 5th edition, p. 350, 1984 \ - # From original source: Goff & Gratch (1946), p. 107) \ - # Goff Gratch formulation; Source : Smithsonian \ - # Meteorological Tables, 5th edition, p. 350, 1984 \ - # From original source: Goff and Gratch (1946), p. 107") - Ts = 373.16 # steam point temperature in K - ews = 1013.246 # saturation pressure at steam point temperature - Psat = np.power(10, -7.90298*(Ts/T-1)+5.02808*np.log10(Ts/T) - - 1.3816e-7*(np.power(10, 11.344*(1-T/Ts))-1) + - 8.1328e-3*(np.power(10, -3.49149*(Ts/T-1))-1) + - np.log10(ews)) - if (meth == 'CIMO'): - logging.info("Source: Annex 4B, Guide to Meteorological \ - Instruments and Methods of Observation, \ - WMO Publication No 8, 7th edition, Geneva, 2008. \ - (CIMO Guide)") - Psat = (6.112*np.exp(17.62*temp/(243.12+temp)) * - (1.0016+3.15e-6*P-0.074/P)) - if (meth == 'MagnusTetens'): - logging.info("Source: Murray, F. W., On the computation of \ - saturation vapor pressure, J. Appl. Meteorol., \ - 6, 203-204, 1967.") - Psat = np.power(10, 7.5*(temp)/(temp+237.5)+0.7858) - # Murray quotes this as the original formula and - Psat = 6.1078*np.exp(17.269388*temp/(temp+237.3)) - # this as the mathematical aquivalent in the form of base e. - if (meth == 'Buck'): - logging.info("Bucks vapor pressure formulation based on \ - Tetens formula. Source: Buck, A. L., \ - New equations for computing vapor pressure and \ - enhancement factor, J. Appl. Meteorol., 20, \ - 1527-1532, 1981.") - Psat = (6.1121*np.exp(17.502*temp/(240.97+temp)) * - (1.0007+(3.46e-6*P))) - if (meth == 'Buck2'): - logging.info("Bucks vapor pressure formulation based on \ - Tetens formula. Source: Buck Research, \ - Model CR-1A Hygrometer Operating Manual, Sep 2001") - Psat = (6.1121*np.exp((18.678-(temp)/234.5)*(temp)/(257.14+temp)) * - (1+1e-4*(7.2+P*(0.0320)+5.9e-6*np.power(T, 2)))) - if (meth == 'WMO'): - logging.info("Intended WMO formulation, originally published by \ - Goff (1957) incorrectly referenced by WMO technical \ - regulations, WMO-NO 49, Vol I, General \ - Meteorological Standards and Recommended Practices, \ - App. A, Corrigendum Aug 2000. and incorrectly \ - referenced by WMO technical regulations, \ - WMO-NO 49, Vol I, General Meteorological Standards \ - and Recommended Practices, App. A, 1988.") - Ts = 273.16 # triple point temperature in K - Psat = np.power(10, 10.79574*(1-Ts/T)-5.028*np.log10(T/Ts) + - 1.50475e-4*(1-np.power(10, -8.2969*(T/Ts-1))) + - 0.42873e-3*(np.power(10, -4.76955*(1-Ts/T))-1) + - 0.78614) - if (meth == 'WMO2000'): - logging.info("WMO formulation, which is very similar to Goff \ - Gratch. Source : WMO technical regulations, \ - WMO-NO 49, Vol I, General Meteorological Standards \ - and Recommended Practices, App. A, \ - Corrigendum Aug 2000.") - Ts = 273.16 # triple point temperature in K - Psat = np.power(10, 10.79574*(1-Ts/T)-5.028*np.log10(T/Ts) + - 1.50475e-4*(1-np.power(10, -8.2969*(T/Ts-1))) + - 0.42873e-3*(np.power(10, -4.76955*(1.-Ts/T))-1) + - 0.78614) - if (meth == 'Sonntag'): - logging.info("Source: Sonntag, D., Advancements in the field of \ - hygrometry, Meteorol. Z., N. F., 3, 51-66, 1994.") - Psat = np.exp(-6096.9385*np.power(T, -1)+16.635794 - - 2.711193e-2*T+1.673952e-5*np.power(T, 2) + - 2.433502*np.log(T))#*(1.0016+P*3.15e-6-0.074/P) - if (meth == 'Bolton'): - logging.info("Source: Bolton, D., The computation of equivalent \ - potential temperature, Monthly Weather Report, 108, \ - 1046-1053, 1980. equation (10)") - Psat = 6.112*np.exp(17.67*temp/(temp+243.5)) - if (meth == 'IAPWS'): - logging.info("Source: Wagner W. and A. Pruss (2002), The IAPWS \ - formulation 1995 for the thermodynamic properties \ - of ordinary water substance for general and \ - scientific use, J. Phys. Chem. Ref. Data, 31(2), \ - 387-535. This is the 'official' formulation from \ - the International Association for the Properties of \ - Water and Steam The valid range of this formulation \ - is 273.16 <= T <= 647.096 K and is based on the \ - ITS90 temperature scale.") - Tc = 647.096 # K : Temperature at the critical point - Pc = 22.064e4 # hPa : Vapor pressure at the critical point - nu = (1-T/Tc) - a1, a2, a3 = -7.85951783, 1.84408259, -11.7866497 - a4, a5, a6 = 22.6807411, -15.9618719, 1.80122502 - Psat = (Pc*np.exp(Tc/T*(a1*nu+a2*np.power(nu, 1.5) + - a3*np.power(nu, 3)+a4*np.power(nu, 3.5) + - a5*np.power(nu, 4)+ a6*np.power(nu, 7.5)))) - if (meth == 'MurphyKoop'): - logging.info("Source : Murphy and Koop, Review of the vapour \ - pressure of ice and supercooled water for \ - atmospheric applications, Q. J. R. Meteorol. \ - Soc (2005), 131, pp. 1539-1565.") - Psat = np.exp(54.842763-6763.22/T-4.210*np.log(T)+0.000367*T + - np.tanh(0.0415*(T-218.8))*(53.878-1331.22/T - - 9.44523*np.log(T)+0.014025*T))/100 - # Calculate saturation pressure over ice ---------------------------------- - elif (phase == 'ice'): - logging.info("Default uses Goff Gratch over ice. There is little \ - ambiguity in the ice saturation curve. Goff Gratch \ - is widely used.") - if (meth == 'MartiMauersberger'): - logging.info("Source : Marti, J. and K Mauersberger, A survey and \ - new measurements of ice vapor pressure at \ - temperatures between 170 and 250 K, GRL 20, \ - 363-366, 1993.") - Psat = np.power(10, -2663.5/T+12.537)/100 - if (meth == 'HylandWexler'): - logging.info("Source Hyland, R. W. and A. Wexler, Formulations \ - for the Thermodynamic Properties of the saturated \ - Phases of H2O from 173.15K to 473.15K, ASHRAE Trans,\ - 89(2A), 500-519, 1983.") - Psat = np.exp(-0.56745359e4/T+0.63925247e1-0.96778430e-2*T + - 0.62215701e-6*np.power(T, 2) + - 0.20747825e-8*np.power(T, 3) - - 0.9484024e-12*np.power(T, 4) + - 0.41635019e1*np.log(T))/100 - if (meth == 'Wexler'): - logging.info("Wexler, A., Vapor pressure formulation for ice, \ - Journal of Research of the National Bureau of \ - Standards-A. 81A, 5-20, 1977.") - Psat = np.exp(-0.58653696e4*np.power(T, -1)+0.22241033e2 + - 0.13749042e-1*T-0.34031775e-4*np.power(T, 2) + - 0.26967687e-7*np.power(T, 3) + - 0.6918651*np.log(T))/100 - if (meth == 'Hardy'): - logging.info("Source Hardy, B., 1998, ITS-90 Formulations for \ - Vapor Pressure, Frostpoint Temperature, Dewpoint \ - Temperature, and Enhancement Factors in the Range \ - -100 to +100°C, The Proceedings of the Third \ - International Symposium on Humidity & Moisture, \ - London, England. These coefficients are updated to \ - ITS90 based on the work by Bob Hardy at Thunder \ - Scientific: http://www.thunderscientific.com/\ - tech_info/reflibrary/its90formulas.pdf \ - The difference to the older ITS68 coefficients used \ - by Wexler is academic.") - Psat = np.exp(-0.58666426e4*np.power(T, -1)+0.2232870244e2 + - 0.139387003e-1*T-0.34262402e-4*np.power(T, 2) + - 0.27040955e-7*np.power(T, 3) + - 0.67063522e-1*np.log(T))/100 - if (meth == 'GoffGratch' or meth == '' or meth == 'IAPWS'): - logging.info("IAPWS does not provide a vapor pressure formulation \ - over ice use Goff Gratch instead.\ - Source : Smithsonian Meteorological Tables, \ - 5th edition, p. 350, 1984") - ei0 = 6.1071 # mbar - T0 = 273.16 # triple point in K - Psat = np.power(10, -9.09718*(T0/T-1)-3.56654*np.log10(T0/T) + - 0.876793*(1-T/T0)+np.log10(ei0)) - if (meth == 'MagnusTetens'): - logging.info("Source: Murray, F. W., On the computation of \ - saturation vapor pressure, J. Appl. Meteorol., 6, \ - 203-204, 1967.") - Psat = np.power(10, 9.5*temp/(265.5+temp)+0.7858) - # Murray quotes this as the original formula and - Psat = 6.1078*np.exp(21.8745584*(T-273.16)/(T-7.66)) - # this as the mathematical aquivalent in the form of base e. - if (meth == 'Buck'): - logging.info("Bucks vapor pressure formulation based on Tetens \ - formula. Source: Buck, A. L., New equations for \ - computing vapor pressure and enhancement factor, \ - J. Appl. Meteorol., 20, 1527-1532, 1981.") - Psat = (6.1115*np.exp(22.452*temp/(272.55+temp)) * - (1.0003+(4.18e-6*P))) - if (meth == 'Buck2'): - logging.info("Bucks vapor pressure formulation based on Tetens \ - formula. Source: Buck Research, Model CR-1A \ - Hygrometer Operating Manual, Sep 2001") - Psat = (6.1115*np.exp((23.036-temp/333.7)*temp/(279.82+temp)) * - (1+1e-4*(2.2+P*(0.0383+6.4e-6*np.power(T, 2))))) - if (meth == 'CIMO'): - logging.info("Source: Annex 4B, Guide to Meteorological \ - Instruments and Methods of Observation, \ - WMO Publication No 8, 7th edition, Geneva, 2008. \ - (CIMO Guide)") - Psat = (6.112*np.exp(22.46*temp/(272.62+temp)) * - (1.0016+3.15e-6*P-0.074/P)) - if (meth == 'WMO' or meth == 'WMO2000'): - logging.info("There is no typo issue in the WMO formulations for \ - ice. WMO formulation, which is very similar to Goff \ - Gratch. Source : WMO technical regulations, \ - WMO-NO 49, Vol I, General Meteorological Standards \ - and Recommended Practices, Aug 2000, App. A.") - T0 = 273.16 # triple point temperature in K - Psat = np.power(10, -9.09685*(T0/T-1)-3.56654*np.log10(T0/T) + - 0.87682*(1-T/T0)+0.78614) - if (meth == 'Sonntag'): - logging.info("Source: Sonntag, D., Advancements in the field of \ - hygrometry, Meteorol. Z., N. F., 3, 51-66, 1994.") - Psat = np.exp(-6024.5282*np.power(T, -1)+24.721994+1.0613868e-2*T - - 1.3198825e-5*np.power(T, 2)-0.49382577*np.log(T)) - if (meth == 'MurphyKoop'): - logging.info("Source : Murphy and Koop, Review of the vapour \ - pressure of ice and supercooled water for \ - atmospheric applications, Q. J. R. Meteorol. Soc \ - (2005), 131, pp. 1539-1565.") - Psat = np.exp(9.550426-5723.265/T+3.53068*np.log(T) - - 0.00728332*T)/100 - s = np.where(temp > 0) - if (s.size[0] >= 1): - logging.info("Independent of the formula used for ice, use \ - Hyland Wexler (water) for temperatures above freezing\ - (see above). Source Hyland, R. W. and A. Wexler, \ - Formulations for the Thermodynamic Properties of the \ - saturated Phases of H2O from 173.15K to 473.15K, \ - ASHRAE Trans, 89(2A), 500-519, 1983.") - Psat_w = np.exp(-0.58002206e4/T+0.13914993e1-0.48640239e-1*T + - 0.41764768e-4*np.power(T, 2) - - 0.14452093e-7*np.power(T, 3) + - 0.65459673e1*np.log(T))/100 - Psat[s] = Psat_w[s] - - return Psat -- GitLab