diff --git a/.gitignore b/.gitignore index 0c6c9d0..a348e50 100644 --- a/.gitignore +++ b/.gitignore @@ -1,5 +1 @@ -/.pydevproject -/.project -/dist/ -/build/ -/.settings/ +/__pycache__/ diff --git a/README.md b/README.md index d3442fd..2bb6aad 100644 --- a/README.md +++ b/README.md @@ -1,83 +1,93 @@ # NAME -penmon - Implementation of weather station class in Python that supports +penmon - Implementation of weather station class in Python that supports Penman-Monteith ETo Equation. -# INSTALL +## INSTALL - $ pip install penmon +```shell +pip install git+ +``` -# USAGE +## USAGE - import penmon as pm - - ### create a station class with known location and elevation - station = pm.Station(latitude=41.42, altitude=109) - station.anemometer_height = 10 +getting a day instance for August 16th 2022 - ### getting a day instance for August 16th - day = station.day_entry(238, - temp_min = 19.5, - temp_max = 25.6, - wind_speed = 2.5, - humidity_mean = 65, - radiation_s =2 5.6 - ) - print("ETo for this day is", day.eto() ) +```py3 +import penmon as pm -# DESCRIPTION +### create a station class with known location and elevation +station = pm.Station(latitude=41.42, altitude=109) +station.anemometer_height = 10 -Full implementation of *Penman-Monteith* *ETo* equation based on +day = station.day_entry( + "2022-08-16", + temp_min = 19.5, + temp_max = 25.6, + wind_speed = 2.5, + humidity_mean = 65, + radiation_s =2 5.6 +) +print("ETo for this day is", day.eto() ) +``` + +### NOTE: radiation_s is given in MJ m-2 day-1. Divide by 11.6 for W m-2 + +## DESCRIPTION + +Full implementation of *Penman-Monteith* *ETo* equation based on [UAN-FAO Irrigation and Drainage Paper 56][1]. -*Penman-Monteith* equation is used to calculate reference crop evapotranspiration *ETo* -for a given location using available climate data. This method provides many ways of estimating -missing climate data. +*Penman-Monteith* equation is used to calculate reference crop evapotranspiration +*ETo* for a given location using available climate data. This method provides +many ways of estimating missing climate data. Following are the least data required to estimate ETo (But the more data you provide the better the accuracy gets): - * latitude - * elevation - * date - * daily minimum temperature - * daily maximum temperature +- latitude +- elevation +- date +- daily minimum temperature +- daily maximum temperature It can do this by making certain assumptions about the climate. These assumptions can be fine-tuned. *Climate*-class is responsible for setting these assumptions. We'll talk more about it later (see Climate class below) -# OVERVIEW +## OVERVIEW -To calculate *ETo*, including intermediate atmospheric data you first need to +To calculate *ETo*, including intermediate atmospheric data you first need to define an instance of a **Station** with a known *latitude* and *altitude*. Then you request the station to create an instance of a **DayEntry**, which represents -a single day with a known date. We then set whatever data we know about that particular -day, and ask the **day** to calculate information that we do not know, including *ETo*. +a single day with a known date. We then set whatever data we know about that particular +day, and ask the **day** to calculate information that we do not know, including +*ETo*. -# BEWARE +## BEWARE -This is pre-release version of the library and intended for review only. API of -the class may change in future releases. Do not assume backwards compatability +This is pre-release version of the library and intended for review only. API of +the class may change in future releases. Do not assume backwards compatability in future releases. Consults **CHANGES** file before upgrading! -# Station CLASS +## Station CLASS -## CREATE A STATION instance +### CREATE A STATION instance - import penmon as pm + import penmon as pm - station = pm.Station(latitude=41.42, altitude=109) + station = pm.Station(latitude=41.42, altitude=109) -*latitude* must be a float between -90 and 90. *altitudu* must be a positive -integer. These arguments values are of utmost importance. Please make sure these -data are as accurate as you can get them be. *latitude* is used to calculate -**sunset hour angle** (Eq. 25) and **Extraterrestrial radiation** (Eq. 21), +*latitude* must be a float between -90 and 90. *altitudu* must be a positive +integer. These arguments values are of utmost importance. Please make sure these +data are as accurate as you can get them be. *latitude* is used to calculate +**sunset hour angle** (Eq. 25) and **Extraterrestrial radiation** (Eq. 21), which in turn, along with date, is used to calculate solar radiation! -*altitude* is used to calculate **atmospheric pressure**, which in turn is used to -calculate **psychrometric constant**, which in turn is used to calculate **vapour pressure**, -which is used to calculate **net longwave radiation**. As you can see, these +*altitude* is used to calculate **atmospheric pressure**, which in turn is used to +calculate **psychrometric constant**, which in turn is used to calculate +**vapour pressure**, +which is used to calculate **net longwave radiation**. As you can see, these very innocent looking numbers are pretty much backbone of the whole equation. Show them respect! @@ -86,123 +96,135 @@ Show them respect! Above *station* assumes its anemometer height is 2m above ground surface. If it's not, you can set the height as: - station.anemometer_height = 10 - + station.anemometer_height = 10 + Now all the wind_speed information is assumed to be wind speed at 10m height. This is important information, since ETo is calculated with speed close the crop surface, which is 2m. Library uses logarithmic algorithm -to convert the data accordingly. Again, important setting! Shoud you wish to +to convert the data accordingly. Again, important setting! Shoud you wish to access calculated wind speed at 2m use *wind_speed_2m()* method: - day.wind_speed=2.0; - u2m = day.wind_speed_2m(); - +```py3 +day.wind_speed=2.0; +u2m = day.wind_speed_2m(); +``` + In the above example *u2m* returns *2.0* if the anemometer was set to 2 meters. If it was set to 10m, it returns 1.5. If it was set to 50 meters, it reads 1.2 m/s. - + ## STATION CLIMATE Station also makes certain assumptions about its climate. You can set this by creating a new climate instance (see **Climate** class) and set is as: - humid_climate = pm.Climate().humid().coastal().strong_winds() - arid_climate = pm.Climate().arid().interior().moderare_winds() - station.climate = humid_climate +```py3 +humid_climate = pm.Climate().humid().coastal().strong_winds() +arid_climate = pm.Climate().arid().interior().moderare_winds() +station.climate = humid_climate +``` By default it assumes we are in *arid, interior location, with moderate winds*. -If it's in arid climate, it makes certain assumption about the dew point temperature. +If it's in arid climate, it makes certain assumption about the dew point temperature. This temperature will be used to calculate relative humidity if humidity data -is missing. It deduces dew_point temperature by looking at the temp_min of the record. +is missing. It deduces dew_point temperature by looking at the temp_min of the record. In particulare, for arid location it substracts 2.0 degrees from temp_min. In humid location it treats temp_min as temp_dew. In the following example we set dew_point temperature 4.0 below temp_min - climate=pm.Climate() - # above is the same as saying: - climate=pm.Climate().arid().interior().moderate_winds() - - climate.dew_point_difference=4.0 - - station.climate=climate; - - # from now on if humidity data is missing it substtracts 4.0 degrees - from temp_min to take a guess at temp_dew - +```py3 +climate=pm.Climate() +# above is the same as saying: +climate=pm.Climate().arid().interior().moderate_winds() + +climate.dew_point_difference=4.0 + +station.climate=climate; + +# from now on if humidity data is missing it substtracts 4.0 degrees +from temp_min to take a guess at temp_dew +``` + ## REFERENCE CROP It assumes it will be calculating ETo for a refernce crop. According -to the original paper reference crop is assumed to be grass of 0.12m height, -aerodynamic resistance of 70 and albedo of 0.23. It you wish to calculate ETo +to the original paper reference crop is assumed to be grass of 0.12m height, +aerodynamic resistance of 70 and albedo of 0.23. It you wish to calculate ETo for a different reference crop you may do so as: - different_crop = pm.Crop(resistance_a=208, albedo=0.23, height=0.12) - station.ref_crop = different_crop - -## ITERATE THROUGH INDIVIDUAL ENTRIES - -Based on the above example this station has no climate data available at this -moment. But if it were you would've been able to iterate through these records + different_crop = pm.Crop(resistance_a=208, albedo=0.23, height=0.12) + station.ref_crop = different_crop + +### ITERATE THROUGH INDIVIDUAL ENTRIES + +Based on the above example this station has no climate data available at this +moment. But if it were you would've been able to iterate through these records like following: - for a_day in station.days: - # do stuff with a_day + for a_day in station.days: -# DayEntry class +- do stuff with a_day + +### DayEntry class Once we have station data available we work with a day at a time. We first need to get a single day, identified by a day number: - day = station.day_entry(238) - + day = station.day_entry(238) + *day* is an instance of **DayEntry** class. *238* above represents *August 26th* - it is 238th day of the year. Day number can only be an integer in 1-366 range. It also supports a date string: - day = station.day_entry("2020-08-26") - day.day_number # returns 238 - -If you have to pass a date string that has a different template than "%Y-%m-%d", + day = station.day_entry("2020-08-26") + day.day_number # returns 238 + +If you have to pass a date string that has a different template than "%Y-%m-%d", you can pass your template string to the method as follows. Above example assumes following default *date_template* value: - day = station.day_entry("2020-08-26", date_template="%Y-%m-%d") - +```py3 +day = station.day_entry("2020-08-26", date_template="%Y-%m-%d") +``` + To learn more about date template placeholders refer to *strptime()* manual either -from *datetime* module, or by refering to your system's *strptime* manual +from *datetime* module, or by refering to your system's *strptime* manual (available in all linux/unix machines). Based on this day number alone library is able to calculate inverse -of the relative Sun-Earth distance (Eq. 23), solar declination (Eq. 24), -which is used to calculate possible daylight hours for that particular day for -that particular latitude. It's amazing how much information we can deduce based +of the relative Sun-Earth distance (Eq. 23), solar declination (Eq. 24), +which is used to calculate possible daylight hours for that particular day for +that particular latitude. It's amazing how much information we can deduce based on this single number! -Once you have *day* instance now all the fun begins! +Once you have *day* instance now all the fun begins! Everyday has following attributes: - - day_number (whatever we passed to the constructor) - - station (set automatically by station) - - # following are required to be set: - - temp_min - - temp_max - - # following are better to be set: - - temp_dew - - wind_speed - - humidity_min - - humidity_max - - radiation_s - Solar Radiation - - # following are optional, if above attributes are known - - radiation_a - - temp_mean - - temp_dry - dry bulb temperature - - temp_wet - wet bulb temperature - +- day_number (whatever we passed to the constructor) +- station (set automatically by station) + +### The following variables are required to be set + +- temp_min +- temp_max + +### The following are better to be set + +- temp_dew +- wind_speed +- humidity_min +- humidity_max +- radiation_s - Solar Radiation + +### The following are optional, if above attributes are known + +- radiation_a +- temp_mean +- temp_dry - dry bulb temperature +- temp_wet - wet bulb temperature + ## INTERMEDIATE CALCULATIONS - + Before setting any of the above attributes following information is available for us. These information do not use any recorded data, but only uses mostly astronomical calculations. @@ -217,27 +239,27 @@ is just an alias to it. Value returned is a pressure in *kPa*. If you wish to convert it to mercury scale multiply it by *7.50* to get *mm Hg*, or *0.295* to get *in Hg*. -### day.latent_heat_of_vaporization() +### day.latent_heat_of_vaporization -Returns 2.45 +- 2.45 -### day.specific_heat() +### day.specific_heat -Returns 0.001013 +- 0.001013 ### day.psychrometric_constant() For the elevation in this example returns 0.0665. -### day.saturation_vapour_pressure(T) +### day.saturation_vapour_pressure(temp) -Given *T* - temperature returns saturation vapour pressure. For example for T=25.0 -returns 3.168 +Given *temp* - temperature returns saturation vapour pressure. For example for T=25.0 +returns 3.168 -### day.slope_of_saturation_vapour_pressure(T) +### day.slope_of_saturation_vapour_pressure(temp) -Given *T* - temperature returns slope of saturation vapour pressure curve. For example -for T=25.0 returns 0.188684 +Given *temp* - temperature returns slope of saturation vapour pressure curve. +E.g. for temp=25.0 returns 0.188684 ### day.relative_sun_distance() @@ -245,10 +267,10 @@ Returns inverse of relative Earth-Sun distance. For this example returns 0.981 ### day.solar_declination() -For this example returns 0.172. This is angle in radians. To convert this to +For this example returns 0.172. This is angle in radians. To convert this to degrees multiply by 180 and devide to *math.pi()*. For example, 0.172 rad is the same as 9.8 degrees. So Sun's declination was 9.8 degrees north relative to equatorial -plane. If the value were negative it would've meant the Sun is declined to the +plane. If the value were negative it would've meant the Sun is declined to the south of the equatorial plane. ### day.sunset_hour_angle() @@ -256,17 +278,17 @@ south of the equatorial plane. For this example, returns 1.725. See *day.solar_declination* to convert this to degrees. -### day.R_a() +### day.r_a() Returns extraterrestrial radiation in MJ/m2/day. For this given example returns 40.3. -### day.R_a_in_mm() +### day.ra_to_mm() The same as above but in mm. 16.4 for this example. ### day.daylight_hours() -Possible daylight hours for this day. +Possible daylight hours for this day. ### day.solar_radiation() @@ -277,14 +299,14 @@ in MJ/m2/day Same as above, but returns in mm. -### day.R_so() +### day.clear_sky_radiation() Calculates Clear-Sky solar radiation. -### day.R_ns() +### day.net_solar_radiation() -Given hours of direct sunlight (in *day.sunshine_hours* attribute) calculates -net solar radiation (or shortwave radiation) Takes into account albedo of the crop. +Given hours of direct sunlight (in *day.sunshine_hours* attribute) calculates +net solar radiation (or shortwave radiation) Takes into account albedo of the crop. ### day.soil_heat_flux() @@ -295,44 +317,42 @@ Returns 0.00 for daily calculations. At this point we still cannot calculate *ETo*. some vital information are still missing. Let's record some data for this day: - day.temp_min = 19.5 - day.tmep_max = 25.6 + day.temp_min = 19.5 + day.tmep_max = 25.6 Now we can calculate *ETo* for the first time: - day.eto() # returns 3.2mm - + day.eto() # returns 3.2mm + To calculate ETo with the current recorded data the library did lots of assumptions -and calculations based on empirical data. Let's help it further by recording +and calculations based on empirical data. Let's help it further by recording other important information: - day.temp_dew = 15.0 - day.eto() # returns 3.58mm + day.temp_dew = 15.0 + day.eto() # returns 3.58mm + + day.wind_speed=2.50 + day.eto()# returns 3.82mm - day.wind_speed=2.50 - day.eto()# returns 3.82mm - Recording solar radiation gets us the most accurate ETo: - day.radiation_s = 25.0 - day.eto() # returns 5.04m - -# TODO AND ISSUES + day.radiation_s = 25.0 + day.eto() # returns 5.04m +## TODO AND ISSUES See: [Issues at github.com][3] - -# SEE ALSO - +## SEE ALSO *libpenmon* - port of the current module into C++. See [github.com/sherzodr/libpenmon][2] - -# AUTHOR +## AUTHOR Sherzod Ruzmetov - + Revised by Benjamin Schmidt + [1]: http://www.fao.org/3/X0490E/x0490e00.htm [2]: https://github.com/sherzodr/libpenmon - [3]: https://github.com/sherzodr/penmon/issues \ No newline at end of file + [3]: https://github.com/sherzodr/penmon/issues + diff --git a/eto.py b/eto.py index 7e2beee..21dbd28 100644 --- a/eto.py +++ b/eto.py @@ -1,10 +1,10 @@ """ Penman-Monteith Equation implementation in Python. -Full implementation of Penman-Monteith ETo equation based on UAN-FAO +Full implementation of Penman-Monteith ETo equation based on UAN-FAO [Irrigation and Drainage Paper 56](http://www.fao.org/3/X0490E/x0490e00.htm) -Penman-Monteith equation is used to calculate reference crop evapotranspiration (ETo) +Penman-Monteith equation is used to calculate reference crop evapotranspiration (ETo) for a given location using available climate data. This method provides many ways of estimating missing climate data using minimal data. @@ -97,7 +97,7 @@ def __init__(self, latitude, longitude, altitude, anemometer_height=2, timezone_ self.anemometer_height = anemometer_height self.climate = Climate() self.ref_crop = Crop() - + def day_entry(self, day_number, date_template="%Y-%m-%d", temp_min=None, temp_max=None, @@ -108,13 +108,13 @@ def day_entry(self, day_number, date_template="%Y-%m-%d", sunshine_hours=None ): """ - Given a day number (integer type from 1-366) returns a **StationDay*** instance for + Given a day number (integer type from 1-366) returns a **StationDay*** instance for that day. Logs the day in *days* attribute of the **Station()** class. If it receives a string it expects it to be in "yyyy-mm-dd" format, in which case it parses the string into **datetime** and calculates day number - If your date format is different than assumed, you can adjust *date_template* + If your date format is different than assumed, you can adjust *date_template* as the second parameter. For example, following all three lines are identical day = station.day_entry(229) @@ -266,12 +266,12 @@ class TimeEntry: """ Represents a single day retrieved from the Station. - This class is usually not instantiated directly. It's instantniated by the + This class is usually not instantiated directly. It's instantniated by the **Station()**'s day_entry() method, passing all reuqired state data. - Since bulk of Penman-Moneith is concerned with a daily ETo **StationDay** is - heart of the module. Penman-Monteith equatoin is implemented within the - methods of **StationDay**. + Since bulk of Penman-Moneith is concerned with a daily ETo **StationDay** is + heart of the module. Penman-Monteith equatoin is implemented within the + methods of **StationDay**. All meteorological data are stored within this class instance. """ @@ -336,7 +336,7 @@ def __init__(self, time, station): def wind_speed_2m(self): """ - Returns wind speed at 2m height. + Returns wind speed at 2m height. If this information is already logged, returns as is. If anemometer of the Station is located higher and wind speed information is available it @@ -395,7 +395,7 @@ def latent_heat_of_vaporization(self): def specific_heat(self): """ - constant: 1.013*10**(-3) + constant: 1.013*10**(-3) """ return 1.013 * 10 ** (-3) @@ -562,7 +562,7 @@ def daylight_hours(self): def solar_radiation_in_mm(self): """ - Alias to *solar_radiation(n)* but converts the output to mm equivalent, + Alias to *solar_radiation(n)* but converts the output to mm equivalent, rounded to 1 decimal. """ rs = self.solar_radiation() @@ -582,7 +582,7 @@ def R_so(self): def R_ns(self): """ - Net solar or net shortwave radiation. Uses Crop's albedo in calculations. (Eq. 38). + Net solar or net shortwave radiation. Uses Crop's albedo in calculations. (Eq. 38). Return radiation in MJ/m2/day """ @@ -645,10 +645,10 @@ def RH(self, T): self.saturation_vapour_pressure(T)), 3) def RH_mean(self): - + if self.humidity_mean != None: return self.humidity_mean - + if self.temp_min and self.temp_max: return int(round(( self.RH(self.temp_min) + self.RH(self.temp_max) ) / 2, 0)) @@ -666,7 +666,7 @@ def air_conductance_coefficient(self): def eto_hargreaves(self): """ - ETo estimating using Hargreaves euqation. If wind and humidty information is + ETo estimating using Hargreaves euqation. If wind and humidty information is available, or can be estimated thsi equation is not recommended. ( Eq. 52 ) """ Tmean = (self.temp_max + self.temp_min) / 2 @@ -884,8 +884,8 @@ def air_conductance_coefficient(self): return 37 # return the official approximation class Climate: - """ - Represents a default climate according to *UN-FAO Paper 56*. If + """ + Represents a default climate according to *UN-FAO Paper 56*. If module has to make any assumptions regarding the climate it consults this class for any clues. If you wish to not use any assumptions and rely soleley on logged station data (if such is available) you may set @@ -1023,9 +1023,9 @@ def describe(self): humidity = "humid, or sub-humid region" return f""" -Climate can be described as located in {location} area, {humidity}. -Average wind speeds are estimated at {self.average_wind_speed}. -Dew point temperature is usually runs {self.dew_point_difference}C lower than +Climate can be described as located in {location} area, {humidity}. +Average wind speeds are estimated at {self.average_wind_speed}. +Dew point temperature is usually runs {self.dew_point_difference}C lower than day's minimal observed temperature""" @@ -1046,7 +1046,7 @@ def __init__(self, resistance_a=208, albedo=0.23, height=0.12): class StationDay(DayEntry): """ - Left here for backwards-compatability with earlier versions of + Left here for backwards-compatability with earlier versions of the library """ diff --git a/penmon/__init__.py b/penmon/__init__.py new file mode 100644 index 0000000..01e3a9d --- /dev/null +++ b/penmon/__init__.py @@ -0,0 +1,2 @@ +from .eto import (CHECK_RADIATION_RANGE, CHECK_SUNSHINE_HOURS_RANGE, Climate, + Crop, DayEntry, Station) diff --git a/src/penmon/eto.py b/penmon/eto.py similarity index 61% rename from src/penmon/eto.py rename to penmon/eto.py index bee4ead..0b98f1e 100644 --- a/src/penmon/eto.py +++ b/penmon/eto.py @@ -1,877 +1,853 @@ -""" -Penman-Monteith Equation implementation in Python. - -Full implementation of Penman-Monteith ETo equation based on UAN-FAO -[Irrigation and Drainage Paper 56](http://www.fao.org/3/X0490E/x0490e00.htm) - -Penman-Monteith equation is used to calculate reference crop evapotranspiration (ETo) -for a given location using available climate data. This method provides many ways of estimating -missing climate data using minimal data. - -Homepage of the project: https://github.com/sherzodr/penmon - -""" - -import math -import datetime as dt - -CHECK_RADIATION_RANGE = True -CHECK_SUNSHINE_HOURS_RANGE = True - -def is_number(s): - try: - float(s) - return True - except ValueError: - return False - - -class Station: - """ Class that implements a weather station at a known latitude and elevation.""" - - def __init__(self, latitude, altitude, anemometer_height=2): - """ - Required parameters: - - :param latitude: latitude of the location in decimal format. For southern - hemisphere negative number must be used - :type latitude: float - - :param altitude: altitude (elevation) of the location in meters - :type altitude: int - - :param anemometer_height=2: height of the anemometer (wind-speed) - measuring device - :type anemometer_height: int - - Following are additional attributes that you can get/set on this station - after instantiation: - - * latitude_rad - latitude in radian, alculated based on latitude - * days - dictionary of days recorded (or calculated) by this station - * climate - set to default **Climate()** instance - * ref_crop - instance of **Crop** class, which sets default chracteristics - of the reference crop according to the paper. - - Should you wish to change assumes Climate and Crop characteristics - you can do so after the object is innitialized, like so: - - station=Station(41.42, 109) - station.ref_crop = Crop(albedo=0.25, height=0.35) - """ - - if not type(latitude) is float: - raise TypeError("latitude must be a float") - - if latitude < -90.0 or latitude > 90.0: - raise Exception("latitude must be between -90.0 and 90.0") - - if not type(altitude) is int: - raise TypeError("altitude must be an integer") - - if altitude < 0: - raise Exception("'altitude' must be above 0") - - self.latitude = latitude - self.altitude = altitude - self.latitude_rad = round((math.pi / 180 * self.latitude), 3) - self.days = {} - - # setting default parameters for the station - self.anemometer_height = anemometer_height - self.climate = Climate() - self.ref_crop = Crop() - - def day_entry(self, day_number, date_template="%Y-%m-%d", - temp_min=None, - temp_max=None, - temp_mean=None, - wind_speed=None, - humidity_mean=None, - humidity_min=None, - humidity_max=None, - radiation_s=None, - sunshine_hours=None - ): - """ - Given a day number (integer type from 1-366) returns a **StationDay*** instance for - that day. Logs the day in *days* attribute of the **Station()** class. - - If it receives a string it expects it to be in "yyyy-mm-dd" format, in which case - it parses the string into **datetime** and calculates day number - - If your date format is different than assumed, you can adjust *date_template* - as the second parameter. For example, following all three lines are identical - - day = station.day_entry(229) - day = station.day_entry("2020-08-16") - day = station.day_entry('08/16/2020', '%m/%d/%Y') - - You can pass the following named-parameters to the method: - - - temp_min - - temp_max - - wind_speed - - radiation_s - - sunshine_hours - - If *radiation_s* and *sunshine_hours* is out of range for this location - for this date (based on solar declination, sun-distance and daylight hours) - raises ValueError exception. - """ - - if type(day_number) is str: - try: - dt1 = dt.datetime.strptime(day_number, date_template) - except ValueError: - raise ValueError( - "Date must be in YYYY-MM-DD format (ex: 2020-09-28)") - - dt0 = dt.datetime(dt1.year, 1, 1) - dates_delta = dt1 - dt0 - day_number = dates_delta.days + 1 - - if not type(day_number) is int: - raise TypeError("'day_number' must be an integer") - - if not (day_number >= 1 and day_number <= 366): - raise Exception("'day_number' must be between in the range 1-366") - - day = DayEntry(day_number, self) - - self.days[day_number] = day - day.temp_min = temp_min - day.temp_max = temp_max - day.temp_mean = temp_mean - day.humidity_min=humidity_min - day.humidity_max=humidity_max - day.humidity_mean = humidity_mean - day.wind_speed = wind_speed - - if radiation_s: - if CHECK_RADIATION_RANGE: - if radiation_s <= day.R_so(): - day.radiation_s = radiation_s - else: - raise ValueError("Raditaion out of range") - else: - day.radiation_s = radiation_s - - if sunshine_hours: - if CHECK_SUNSHINE_HOURS_RANGE: - if sunshine_hours <= day.daylight_hours(): - day.sunshine_hours = sunshine_hours - else: - raise ValueError("Sunshine hours out of range") - else: - day.sunshine_hours = sunshine_hours - - return day - - get_day = day_entry - - def atmospheric_pressure(self): - """ - Calculates atmospheric pressure *in kPa* based on station's altitude. (Eq. 7) - """ - return round(101.3 * ((293 - 0.0065 * self.altitude) / 293) ** 5.26, 2) - - def describe(self): - """ - Describes the station and all its assumptions in human-friendly text - """ - return self - - -class DayEntry: - """ - Represents a single day retrieved from the Station. - - This class is usually not instantiated directly. It's instantniated by the - **Station()**'s day_entry() method, passing all reuqired state data. - - Since bulk of Penman-Moneith is concerned with a daily ETo **StationDay** is - heart of the module. Penman-Monteith equatoin is implemented within the - methods of **StationDay**. - - All meteorological data are stored within this class instance. - """ - - def __init__(self, day_number, station): - """ - *day_number* and *station* are two required arguments passed to - instantiate the class. Day number must be between 1-366. It represents a - single day in a year, *1* meaning January 1st, and 365 (or 366) being - Decement 31st. *station* must be instance of **Station** class. It reads - such important information from the station as its location, altitude, - climate conditions. Without this information it's impossible to - calculate solar radiation and humidity data. - - Following attributes of the class are available. They can be both set - and get. - - - day_number - - station - references **Station** class. - - temp_min - - temp_max - - temp_mean - - temp_dew - - temp_dry - - temp_wet - - humidity_mean - - humidity_min - - humidity_max - - vapour_pressure - - wind_speed - - radiation_s - - stephan_boltzman_constant - - climate - convenient reference to station.climate - - sunshine_hours - - - """ - - self.day_number = day_number - self.station = station - - self.temp_min = None - self.temp_max = None - self.temp_mean = None - self.humidity_mean = None - self.humidity_min = None - self.humidity_max = None - self.wind_speed = None - self.radiation_s = None - self.temp_dew = None - self.temp_dry = None - self.temp_wet = None - self.climate = station.climate - self.stephan_boltzmann_constant = 4.903 * (10 ** -9) - self.vapour_pressure = None - self.sunshine_hours = None - - def wind_speed_2m(self): - """ - Returns wind speed at 2m height. - - If this information is already logged, returns as is. If anemometer of - the Station is located higher and wind speed information is available it - converts this information to wind speed as 2ms based on logarithimc - conversion (Eq. 47) - - If wind speed was not logged for this date, and if climate is known - tries to rely on Climate data to estimate average wind speed - """ - - # if wind speed at 2m height is given, return it - if self.wind_speed and (self.station.anemometer_height == 2): - return self.wind_speed - - # if wind speed at height different than 2m is given, calculate wind - # speed at 2m - if self.wind_speed and self.station.anemometer_height != 2: - return round(self.wind_speed * (4.87 / - math.log(67.8 * self.station.anemometer_height - 5.42)), 1) - - # if we reach this far no wind information is available to work with. we - # consult if station has any climatic data, in which case we try to - # deduce information off of that: - if self.station.climate: - return self.station.climate.average_wind_speed - - return None - - def dew_point(self): - """ - If *temp_dew* attribute is logged returns as is. If this data was not - logged, but *temp_min* data is available tries to estimate *temp_dew* - based on Station's Climate. If either is not possible returns *None*. - """ - if self.temp_dew: - return self.temp_dew - - if self.temp_min and self.climate: - return self.temp_min - self.climate.dew_point_difference - - def atmospheric_pressure(self): - """ - Calculates atmospheric pressure *in kPa* based on station's altitude. (Eq. 7) - """ - return self.station.atmospheric_pressure() - - def latent_heat_of_vaporization(self): - """ - constant *2.45* - """ - return 2.45 - - def specific_heat(self): - """ - constant: 1.013*10**(-3) - """ - return 1.013 * 10 ** (-3) - - def psychrometric_constant(self): - """ - Calculates psychrometric constant based on Station's altitude (and - atmospheric pressure). (Eq. 8) - """ - return round(0.665 * 10 ** (-3) * self.atmospheric_pressure(), 6) - - def Tmean(self): - """ - If *temp_mean* is logged returns is as is. If *temp_min* and *temp_max* - are available computes *Tmean* based on these data. If none are - available returns *None*. (Eq. 9) - """ - if self.temp_mean: - return self.temp_mean - - if self.temp_max and self.temp_min: - return ((self.temp_max + self.temp_min) / 2) - - return None - - def saturation_vapour_pressure(self, T): - """ - Calculates saturation vapour pressure for a given temperature. (Eq. 11) - """ - return round((0.6108 * (2.7183 ** (17.27 * T / (T + 237.3)))), 3) - - def mean_saturation_vapour_pressure(self): - """ - Given *temp_max* and *temp_min* calculates mean saturation vapour pressure. (Eq. 12) - """ - if self.temp_max and self.temp_min: - vp_max = self.saturation_vapour_pressure(self.temp_max) - vp_min = self.saturation_vapour_pressure(self.temp_min) - return (vp_max + vp_min) / 2 - - if self.temp_mean: - return self.saturation_vapour_pressure(self.temp_mean) - - def slope_of_saturation_vapour_pressure(self, T): - """ - Calculates slope of the saturation vapour pressure curve for a given - temperature. It's the required information to calculate ETo. (Eq. 13) - """ - return round((4098 * (0.6108 * 2.7183 ** (17.27 * T / (T + 237.3)))) - / ((T + 237.3) ** 2), 6) - - def actual_vapour_pressure(self): - """ - Attepmts to calculate vapour pressure based on several available weather - data. - - If *temp_dry* and *temp_wet* data are logged (psychrometric data) uses - (Eq. 15) to calculate actual vapour pressure. If only *temp_dew* - information is logged uses (Eq. 14) to calculate actual vapour pressure. - If *humidity_max* and *humidity_min* are logged uses (Eq. 17) to - calculate vapour pressure. If only *humidity_max* is known uses (Eq. 18) - to estimate actual vapour pressure. If only *humidity_mean* is known - uses (Eq. 19) to estimate actual vapour pressure. - """ - if self.vapour_pressure: - return self.vapour_pressure - - if self.temp_dry and self.temp_wet: - vp_wet = self.saturation_vapour_pressure(self.temp_wet) - psychrometric_constant = self.psychrometric_constant() - return round(vp_wet - psychrometric_constant * - (self.temp_dry - self.temp_wet), 3) - - if self.humidity_max and self.humidity_min and self.temp_max and self.temp_min: - vp_min = self.saturation_vapour_pressure(self.temp_min) - vp_max = self.saturation_vapour_pressure(self.temp_max) - return round((vp_min * (self.humidity_max / 100) + - vp_max * (self.humidity_min / 100)) / 2, 3) - - if self.humidity_max and self.temp_min: - vp_min = self.saturation_vapour_pressure(self.temp_min) - return round(vp_min * (self.humidity_max / 100), 3) - - if self.humidity_mean and self.temp_max and self.temp_min: - vp_min = self.saturation_vapour_pressure(self.temp_min) - vp_max = self.saturation_vapour_pressure(self.temp_max) - return round((self.humidity_mean / 100) * ((vp_max + vp_min) / 2), 3) - - if self.dew_point(): - return round(self.saturation_vapour_pressure(self.dew_point()), 3) - - def vapour_pressure_deficit(self): - if self.temp_min and self.temp_max: - vp_min = self.saturation_vapour_pressure(self.temp_min) - vp_max = self.saturation_vapour_pressure(self.temp_max) - actual_vp = self.actual_vapour_pressure() - return round(((vp_min + vp_max) / 2) - actual_vp, 3) - - def relative_sun_distance(self): - """ - Eq. 23 - """ - return round(1 + 0.033 * math.cos((2 * math.pi / 365) * self.day_number), 3) - - def solar_declination(self): - """ - Eq. 24 - """ - return round(0.409 * math.sin((2 * math.pi / 365) * self.day_number - 1.39), 3) - - def X(self): - """ - Eq. 27 - """ - x = (1 - math.tan(self.station.latitutde_radians) * - math.tan(self.solar_declination())) - if x <= 0: - x = 0.00001 - return x - - def sunset_hour_angle(self): - """ - Eq. 25 - """ - - return round(math.acos(-1 * math.tan(self.station.latitude_rad) * - math.tan(self.solar_declination())), 3) - - def R_a(self): - """ - Extraterrestrial radiation for daily periods.( Eq. 21 ). - """ - - return round( - 24 * 60 / math.pi * 0.0820 * self.relative_sun_distance() * - ( - ( - self.sunset_hour_angle() * math.sin(self.station.latitude_rad) * - math.sin(self.solar_declination()) - ) + - ( - math.cos(self.station.latitude_rad) * - math.cos(self.solar_declination()) * - math.sin(self.sunset_hour_angle()) - ) - ), - 1) - - def R_a_in_mm(self): - """ - Same as R_a(), but returns the result in mm-equivalents - """ - return round(self.R_a() * 0.408, 1) - - def daylight_hours(self): - """ - Eq. 34 - """ - return round((24 / math.pi) * self.sunset_hour_angle(), 1) - - # Rs - def solar_radiation(self): - """ - If *radiation_s* is logged, returns the value as is. If *sunshine_hours* - attribute of the day class is set returns solar radiation amount in mJ/m2/day. - To convert this number to W/m2 multiply multiply it by 11.57 or divide by - 0.0864. Uses Angstrom equation (Eq. 35). - - If climate data is available, and climate is *island* location and - station elevation is between 0-100m it uses simplified (Eq. 51). This - equation does not use temperature data, but just solar radiation and a - coefficient. - - If station elevation is higher than 100m and/or location is not island - it uses (Eq. 50) that calculates solar radiation by using temperature - data along with a *krs* constant. - - If climate is not known it assumes **n=N**, meaning daily sunshine hours - is the same as daylight hours for the given day and location. - """ - - if self.radiation_s: - # We need to make sure that solar radiation if set, is not - # larger than clear-sky solar radiation - if CHECK_RADIATION_RANGE and (self.radiation_s > self.R_so()): - raise ValueError( - "Solar radiation out ot range. Rso=" + str(self.R_so())) - return self.radiation_s - - n = self.sunshine_hours - if n == None: - # if we are in island location we refer to equation 51 in UAN-FAO - # paper 56 - if (self.station.climate and self.station.climate.island_location - and self.station.altitude < 100): - Ra = self.R_a() - return round((0.7 * Ra) - 4, 1) - - if self.station.climate and self.temp_min and self.temp_max: - # We assume caller has only temperature informaiton, and no - # information on overcast conditions. So we resort to Hargreaves - # and Samani's radiation formula: - climate = self.station.climate - Ra = self.R_a() - krs = 0.16 - if climate.coastal_location: - krs = 0.19 - elif climate.interior_location: - krs = 0.16 - - return round(krs * math.sqrt(self.temp_max - self.temp_min) * Ra, 1) - - else: - n = self.daylight_hours() - - if n and not is_number(n): - raise TypeError("'n' must be a number") - - if n < 0: - raise ValueError("Observed daylight hours cannot be less than 0") - - # n cannot be more than N, which is available daylight hours - if (n > self.daylight_hours()) and CHECK_SUNSHINE_HOURS_RANGE: - raise ValueError("Daylight hours out of range") - - a_s = 0.25 - b_s = 0.50 - N = self.daylight_hours() # this is the maximum possible sunshine duration - Ra = self.R_a() - - return round((a_s + b_s * n / N) * Ra, 1) - - def solar_radiation_in_mm(self): - """ - Alias to *solar_radiation(n)* but converts the output to mm equivalent, - rounded to 1 decimal. - """ - rs = self.solar_radiation() - return round(rs * 0.408, 1) - - # clear-skype solar radiation - def R_so(self): - """ - Calculates clear sky radiation when n=N. Uses (Eq. 36) for elevations - below 100m. Above 100m uses (Eq. 37) - """ - - if self.station.altitude < 100: - return round((0.25 + 0.50) * self.R_a(), 1) - - return round((0.75 + (2 * 10**(-5)) * self.station.altitude) * self.R_a(), 1) - - def R_ns(self): - """ - Net solar or net shortwave radiation. Uses Crop's albedo in calculations. (Eq. 38). - Return radiation in MJ/m2/day - - """ - ref_crop = self.station.ref_crop - - return round((1 - ref_crop.albedo) * self.solar_radiation(), 1) - - def R_nl(self): - """ - Net longwave radiation. (Eq. 39) - """ - - if not (self.temp_max and self.temp_min): - raise Exception( - "Net longwave radiation cannot be calculated without min/max temperature") - - TmaxK = self.temp_max + 273.16 - TminK = self.temp_min + 273.16 - ea = self.actual_vapour_pressure() - rs = self.solar_radiation() - rso = self.R_so() - - sb_constant = self.stephan_boltzmann_constant - return round(sb_constant * ((TmaxK ** 4 + TminK ** 4) / 2) * - (0.34 - 0.14 * math.sqrt(ea)) * - (1.35 * (rs / rso) - 0.35), 1) - - def net_radiation(self): - """ - Net Radiation. (Eq. 40) - """ - ns = self.R_ns() - - try: - nl = self.R_nl() - except Exception as e: - raise(str(e)) - - if (not ns is None) and (not nl is None): - return round(ns - nl, 1) - - def net_radition_in_mm(self): - """ - Same as *net_radiation()*, except returns results in mm-equivalents - """ - net_radition = self.net_radiation() - if net_radition: - return round(net_radition * 0.408, 1) - - def RH(self, T): - - if not is_number(T): - raise TypeError("Number is expected") - - """ - Calculates relative humidity of the air for certain temperature using vapour pressure - """ - - return round(100 * (self.actual_vapour_pressure() / - self.saturation_vapour_pressure(T)), 3) - - def RH_mean(self): - """ - If possible returns mean relative humidity for the day. If humidity_mean was logged in the station - returns it as is. - If min/max humidity values were logged at the station, computes - the mean of the two values - If none was logged, but min/max temperature values were logged it attempts to calculate - relative humidity through saturation vapour pressure: - """ - if self.humidity_mean != None: - return self.humidity_mean - - if (self.humidity_min != None) and (self.humidity_max != None): - return round( ( self.humidity_max + self.humidity_min ) / 2, 0 ) - - if ( self.temp_min != None ) and ( self.temp_max != None ): - return int(round(( self.RH(self.temp_min) + self.RH(self.temp_max) ) / 2, 0)) - - return None - - - def soil_heat_flux(self): - """ - Soil heat flux. Returns 0.00 (daily coefficient) - """ - return 0.00 - - def eto_hargreaves(self): - """ - ETo estimating using Hargreaves euqation. If wind and humidty information is - available, or can be estimated thsi equation is not recommended. ( Eq. 52 ) - """ - Tmean = (self.temp_max + self.temp_min) / 2 - return round(0.0023 * (Tmean + 17.8) * - (self.temp_max - self.temp_min) ** 0.5 * self.R_a(), 2) - - def eto(self): - """ - Eq. 6 - """ - - # if we cannot get wind speed data we revert to Hargreaves formula. - # Which is not ideal! This can happen only if user removed default 'climate' - # reference - if not self.wind_speed_2m(): - return self.eto_hargreaves() - - if self.Tmean() == None: - raise Exception( - "Cannot calculate eto(): temp_mean (mean temperature) is missing") - - try: - net_radiation = self.net_radiation() - except Exception as e: - raise(str(e)) - - Tmean = self.Tmean() - slope_of_vp = self.slope_of_saturation_vapour_pressure(Tmean) - G = self.soil_heat_flux() - u2m = self.wind_speed_2m() - eto_nominator = ( - ( 0.408 * slope_of_vp * (net_radiation - G) ) + - self.psychrometric_constant() * (900 / (Tmean + 273)) * u2m * - self.vapour_pressure_deficit() - ) - - eto_denominator = slope_of_vp + self.psychrometric_constant() * (1 + 0.34 * u2m) - return round(eto_nominator / eto_denominator, 2) - - -class Climate: - """ - Represents a default climate according to *UN-FAO Paper 56*. If - module has to make any assumptions regarding the climate it consults - this class for any clues. If you wish to not use any assumptions and - rely soleley on logged station data (if such is available) you may set - Station's *climate* attribute to *None*. - - station = Station(latitude=-20.5, altitude=200) - station.climate = None - - If you want to set a new climate: - - humid_climate = Climate().humid().coastal().moderate_winds() - - station = Station(latitude=-20.5, altitude=200) - station.climate = humid_climate - """ - - def __init__(self): - """ - Accepts no arguments. Default initialization is as follows: - - - interior_location - - arid_climate - - dew_point_difference = 2 - - average_wind_speed = 2.0 m/s - - k_rs = 0.16 - - - To affect these values use respected methods documented below. - """ - self.interior_location = True - self.coastal_location = False - self.island_location = False - - self.arid_climate = True - self.humid_climate = False - - # Assining default values for the climatic condition to be able to - # calculate missing data accurately - self.dew_point_difference = 2 - self.average_wind_speed = 2.0 - self.k_rs = 0.16 - - def light_winds(self): - """ - Sets *average_wind_speed* to 0.5m/s - """ - self.average_wind_speed = 0.5 - return self - - def moderate_winds(self): - """ - Sets *average_wind_speed* to 2.0m/s - """ - self.average_wind_speed = 2 - return self - - def strong_winds(self): - """ - Sets average_wind_speed to 4m/s - """ - self.average_wind_speed = 4 - return self - - def arid(self): - """ - Sets *arid_climate*, sets *dew_point_difference* to 2 - """ - self.arid_climate = True - self.humid_climate = False - self.dew_point_difference = 2 - return self - - def humid(self): - """ - Sets *humid_climate*, *dew_point_difference* to 0 - """ - self.arid_climate = False - self.humid_climate = True - self.dew_point_difference = 0 - return self - - def interior(self): - """ - Sets *interior_location*, *k_rs* coefficient to *0.16* - """ - self.interior_location = True - self.coastal_location = False - self.island_location = False - self.k_rs = 0.16 - return self - - def coastal(self): - """ - Sets *coastal_location*, *k_rs* to *0.19* - """ - self.interior_location = False - self.coastal_location = True - self.island_location = False - self.k_rs = 0.19 - return self - - def island(self): - """ Sets *island_location*. Sets *k_rs* to 0.19 """ - self.interior_location = False - self.coastal_location = False - self.island_location = True - self.k_rs = 0.19 - return self - - """ - def set_location(self, location): - if location == "coastal": - return self.coastal() - - if location == "interior": - return self.interior() - - if location == "island": - return self.island() - - return self - """ - - def describe(self): - location = "" - if self.interior_location: - location = "interior" - else: - location = "coastal" - - humidity = "" - if self.arid_climate: - humidity = "arid, or semi-arid region" - else: - humidity = "humid, or sub-humid region" - - return f""" -Climate can be described as located in {location} area, {humidity}. -Average wind speeds are estimated at {self.average_wind_speed}. -Dew point temperature is usually runs {self.dew_point_difference}C lower than -day's minimal observed temperature""" - - -class Crop: - """ Represents reference crop as assumed by Penman-Monteith equation.""" - - def __init__(self, resistance_a=208, albedo=0.23, height=0.12): - """ - Reference crop is instantiated with the following assumptions: - - resistance_a = 208 (aerodynamic resistance) - - albedo = 0.23 - - height - 0.12 - """ - self.resistance_a = resistance_a - self.albedo = albedo - self.height = height - - -class StationDay(DayEntry): - """ - Left here for backwards-compatability with earlier versions of - the library - """ - - pass - - -class MonthEntry: - def __init__(self): - pass - - -class WeekEntry: - def __init__(self): - pass - - -class HourEntry: - def __init__(self): - pass +""" +Penman-Monteith Equation implementation in Python. + +Full implementation of Penman-Monteith ETo equation based on UAN-FAO +[Irrigation and Drainage Paper 56](http://www.fao.org/3/X0490E/x0490e00.htm) + +Penman-Monteith equation is used to calculate reference crop evapotranspiration (ETo) +for a given location using available climate data. This method provides many ways of estimating +missing climate data using minimal data. + +Homepage of the project: https://github.com/sherzodr/penmon + +""" + +import datetime as dt +from numpy import nan +import math + +CHECK_RADIATION_RANGE = True +CHECK_SUNSHINE_HOURS_RANGE = True + + +class Station: + """Class that implements a weather station at a known latitude and elevation.""" + + def __init__(self, latitude, altitude, anemometer_height=2): + """ + Instantiate Station object from latitude, altitude, and anemometer height. + + Required parameters: + + :param latitude: latitude of the location in decimal format. For southern + hemisphere negative number must be used + :type latitude: float + + :param altitude: altitude (elevation) of the location in meters + :type altitude: int + + :param anemometer_height=2: height of the anemometer (wind-speed) + measuring device + :type anemometer_height: int + + Following are additional attributes that you can get/set on this station + after instantiation: + + * latitude_rad - latitude in radian, alculated based on latitude + * days - dictionary of days recorded (or calculated) by this station + * climate - set to default **Climate()** instance + * ref_crop - instance of **Crop** class, which sets default chracteristics + of the reference crop according to the paper. + + Should you wish to change assumes Climate and Crop characteristics + you can do so after the object is innitialized, like so: + + station=Station(41.42, 109) + station.ref_crop = Crop(albedo=0.25, height=0.35) + """ + if type(latitude) is not float: + raise TypeError("'latitude' must be a float") + + if latitude < -90.0 or latitude > 90.0: + raise ValueError("'latitude' must be between -90.0 and 90.0") + + if type(altitude) is not int: + raise TypeError("'altitude' must be an integer") + + if altitude < 0: + raise ValueError("'altitude' must be above 0") + + self.latitude = latitude + self.altitude = altitude + self.latitude_rad = round((math.pi / 180 * self.latitude), 3) + self.days = {} + + # setting default parameters for the station + self.anemometer_height = anemometer_height + self.climate = Climate() + self.ref_crop = Crop() + + def day_entry( + self, + day_number, + date_template="%Y-%m-%d", + temp_min=None, + temp_max=None, + temp_mean=None, + wind_speed=None, + humidity_mean=None, + humidity_min=None, + humidity_max=None, + radiation_s=None, + sunshine_hours=None, + ): + """ + Calculate the evapotranspiration for a given day or day of year. + + Given a day number (integer type from 1-366) returns a **StationDay*** instance for + that day. Logs the day in *days* attribute of the **Station()** class. + + If it receives a string it expects it to be in "yyyy-mm-dd" format, in which case + it parses the string into **datetime** and calculates day number + + If your date format is different than assumed, you can adjust *date_template* + as the second parameter. For example, following all three lines are identical + + day = station.day_entry(229) + day = station.day_entry("2020-08-16") + day = station.day_entry('08/16/2020', '%m/%d/%Y') + + You can pass the following named-parameters to the method: + + - temp_min + - temp_max + - wind_speed + - radiation_s + - sunshine_hours + + If *radiation_s* and *sunshine_hours* is out of range for this location + for this date (based on solar declination, sun-distance and daylight hours) + raises ValueError exception. + """ + if type(day_number) is str: + try: + dt1 = dt.datetime.strptime(day_number, date_template) + except ValueError as e: + raise ValueError("Date must be in YYYY-MM-DD format (ex: 2020-09-28)") from e + + dt0 = dt.datetime(dt1.year, 1, 1) + dates_delta = dt1 - dt0 + day_number = dates_delta.days + 1 + + if type(day_number) is not int: + raise TypeError("'day_number' must be an integer") + + if day_number < 1 or day_number > 366: + raise ValueError("'day_number' must be between in the range 1-366") + + day = DayEntry(day_number, self) + + self.days[day_number] = day + day.temp_min = temp_min + day.temp_max = temp_max + day.temp_mean = temp_mean + day.humidity_min = humidity_min + day.humidity_max = humidity_max + day.humidity_mean = humidity_mean + day.wind_speed = wind_speed + + if radiation_s: + if ( + CHECK_RADIATION_RANGE + and radiation_s <= day.clear_sky_rad() + or not CHECK_RADIATION_RANGE + ): + day.radiation_s = radiation_s + elif radiation_s is nan: + pass + else: + raise ValueError("Radiation out of range") + if sunshine_hours: + if ( + CHECK_SUNSHINE_HOURS_RANGE + and sunshine_hours <= day.daylight_hours() + or not CHECK_SUNSHINE_HOURS_RANGE + ): + day.sunshine_hours = sunshine_hours + else: + raise ValueError("Sunshine hours out of range") + return day + + get_day = day_entry + + def atmospheric_pressure(self): + """Calculate atmospheric pressure *in kPa* based on station's altitude (Eq. 7).""" + return round(101.3 * ((293 - 0.0065 * self.altitude) / 293) ** 5.26, 2) + + def describe(self): + """Describe the station and all its assumptions in human-friendly text.""" + return self + + +class DayEntry: + """ + Represents a single day retrieved from the Station. + + This class is usually not instantiated directly. It's instantniated by the + **Station()**'s day_entry() method, passing all reuqired state data. + + Since bulk of Penman-Moneith is concerned with a daily ETo **StationDay** is + heart of the module. Penman-Monteith equatoin is implemented within the + methods of **StationDay**. + + All meteorological data are stored within this class instance. + """ + + def __init__(self, day_number, station): + """ + Instantiate a DayEntry-object. + + *day_number* and *station* are two required arguments passed to + instantiate the class. Day number must be between 1-366. It represents a + single day in a year, *1* meaning January 1st, and 365 (or 366) being + Decement 31st. *station* must be instance of **Station** class. It reads + such important information from the station as its location, altitude, + climate conditions. Without this information it's impossible to + calculate solar radiation and humidity data. + + Following attributes of the class are available. They can be both set + and get. + + - day_number + - station - references **Station** class. + - temp_min + - temp_max + - temp_mean + - temp_dew + - temp_dry + - temp_wet + - humidity_mean + - humidity_min + - humidity_max + - vapour_pressure + - wind_speed + - radiation_s + - stephan_boltzman_constant + - climate - convenient reference to station.climate + - sunshine_hours + - + """ + self.day_number = day_number + self.station = station + self.temp_min = None + self.temp_max = None + self.temp_mean = None + self.humidity_mean = None + self.humidity_min = None + self.humidity_max = None + self.wind_speed = None + self.radiation_s = None + self.temp_dew = None + self.temp_dry = None + self.temp_wet = None + self.climate = station.climate + self.stephan_boltzmann_constant = 4.903 * (10**-9) + self.vapour_pressure = None + self.sunshine_hours = None + self.specific_heat = 1.013 * 10 ** (-3) + self.latent_heat_of_vaporization = 2.45 + + def wind_speed_2m(self): + """ + Calculate the wind speed at 2m height. + + If this information is already logged, returns as is. If anemometer of + the Station is located higher and wind speed information is available it + converts this information to wind speed as 2ms based on logarithimc + conversion (Eq. 47) + + If wind speed was not logged for this date, and if climate is known + tries to rely on Climate data to estimate average wind speed + + if wind speed at 2m height is given, return it + if wind speed at height different than 2m is given, calculate wind + speed at 2m + """ + if self.wind_speed: + if self.station.anemometer_height == 2: + return self.wind_speed + + return round( + self.wind_speed * (4.87 / math.log(67.8 * self.station.anemometer_height - 5.42)), 1 + ) + + # if we reach this far no wind information is available to work with. we + # consult if station has any climatic data, in which case we try to + # deduce information off of that: + if self.station.climate: + return self.station.climate.average_wind_speed + + return None + + def dew_point(self): + """ + Calculate dew point temperature. + + If *temp_dew* attribute is logged returns as is. + If this data was not logged, but *temp_min* data is available tries to estimate *temp_dew* + based on Station's Climate. If either is not possible returns *None*. + """ + if self.temp_dew: + return self.temp_dew + + if self.temp_min and self.climate: + return self.temp_min - self.climate.dew_point_difference + return False + + def atmospheric_pressure(self): + """Calculate atmospheric pressure *in kPa* based on station's altitude (Eq. 7).""" + return self.station.atmospheric_pressure() + + def psychrometric_constant(self): + """ + Calculate psychrometric constant. + + The computation is based on Station's altitude (and atmospheric pressure, Eq. 8). + """ + return round(0.665 * 10 ** (-3) * self.atmospheric_pressure(), 6) + + def interpolate_temp_mean(self): + """ + Interpolate the mean temperature from max and min. + + If *temp_mean* is logged returns is as is. If *temp_min* and *temp_max* + are available computes *temp_mean* based on these data. If none are + available returns *None*. (Eq. 9) + """ + if self.temp_mean: + return self.temp_mean + + if self.temp_max and self.temp_min: + return (self.temp_max + self.temp_min) / 2 + + return None + + def saturation_vapour_pressure(self, temp): + """Calculate saturation vapour pressure for a given temperature (Eq. 11).""" + return round((0.6108 * (2.7183 ** (17.27 * temp / (temp + 237.3)))), 3) + + def mean_saturation_vapour_pressure(self): + """ + Calculate mean saturation vapour pressure. + + Given *temp_max* and *temp_min* calculates mean saturation vapour pressure. (Eq. 12) + """ + if self.temp_max and self.temp_min: + vp_max = self.saturation_vapour_pressure(self.temp_max) + vp_min = self.saturation_vapour_pressure(self.temp_min) + return (vp_max + vp_min) / 2 + + if self.temp_mean: + return self.saturation_vapour_pressure(self.temp_mean) + + def slope_of_saturation_vapour_pressure(self, temp): + """ + Calculate slope of the saturation vapour pressure curve for a given temperature. + + It's the required information to calculate ETo. (Eq. 13) + """ + return round( + (4098 * (0.6108 * 2.7183 ** (17.27 * temp / (temp + 237.3)))) / ((temp + 237.3) ** 2), 6 + ) + + def actual_vapour_pressure(self): + """ + Attempt to calculate vapour pressure based on several available weather data. + + If *temp_dry* and *temp_wet* data are logged (psychrometric data) uses + (Eq. 15) to calculate actual vapour pressure. If only *temp_dew* + information is logged uses (Eq. 14) to calculate actual vapour pressure. + If *humidity_max* and *humidity_min* are logged uses (Eq. 17) to + calculate vapour pressure. If only *humidity_max* is known uses (Eq. 18) + to estimate actual vapour pressure. If only *humidity_mean* is known + uses (Eq. 19) to estimate actual vapour pressure. + """ + if self.vapour_pressure: + return self.vapour_pressure + + if self.temp_dry and self.temp_wet: + vp_wet = self.saturation_vapour_pressure(self.temp_wet) + psychrometric_constant = self.psychrometric_constant() + return round(vp_wet - psychrometric_constant * (self.temp_dry - self.temp_wet), 3) + + if self.humidity_max and self.humidity_min and self.temp_max and self.temp_min: + vp_min = self.saturation_vapour_pressure(self.temp_min) + vp_max = self.saturation_vapour_pressure(self.temp_max) + return round( + (vp_min * (self.humidity_max / 100) + vp_max * (self.humidity_min / 100)) / 2, 3 + ) + + if self.humidity_max and self.temp_min: + vp_min = self.saturation_vapour_pressure(self.temp_min) + return round(vp_min * (self.humidity_max / 100), 3) + + if self.humidity_mean and self.temp_max and self.temp_min: + vp_min = self.saturation_vapour_pressure(self.temp_min) + vp_max = self.saturation_vapour_pressure(self.temp_max) + return round((self.humidity_mean / 100) * ((vp_max + vp_min) / 2), 3) + + if self.dew_point(): + return round(self.saturation_vapour_pressure(self.dew_point()), 3) + + def vapour_pressure_deficit(self): + """Calculate vapour pressure deficit.""" + if self.temp_min and self.temp_max: + vp_min = self.saturation_vapour_pressure(self.temp_min) + vp_max = self.saturation_vapour_pressure(self.temp_max) + actual_vp = self.actual_vapour_pressure() + return round(((vp_min + vp_max) / 2) - actual_vp, 3) + + def relative_sun_distance(self): + """Eq. 23.""" + return round(1 + 0.033 * math.cos((2 * math.pi / 365) * self.day_number), 3) + + def solar_declination(self): + """Eq. 24.""" + return round(0.409 * math.sin((2 * math.pi / 365) * self.day_number - 1.39), 3) + + def x_distance_from_sun(self): + """Eq. 27.""" + x = 1 - math.tan(self.station.latitude_rad) * math.tan(self.solar_declination()) + if x <= 0: + x = 0.00001 + return x + + def sunset_hour_angle(self): + """Eq. 25.""" + return round( + math.acos( + -1 * math.tan(self.station.latitude_rad) * math.tan(self.solar_declination()) + ), + 3, + ) + + def r_a(self): + """Extraterrestrial radiation for daily periods ( Eq. 21 ).""" + return round( + 24 + * 60 + / math.pi + * 0.0820 + * self.relative_sun_distance() + * ( + ( + self.sunset_hour_angle() + * math.sin(self.station.latitude_rad) + * math.sin(self.solar_declination()) + ) + + ( + math.cos(self.station.latitude_rad) + * math.cos(self.solar_declination()) + * math.sin(self.sunset_hour_angle()) + ) + ), + 1, + ) + + def ra_to_mm(self): + """Transform as r_a to mm-equivalents.""" + return round(self.r_a() * 0.408, 1) + + def daylight_hours(self): + """Calculate daylight hours from Eq. 34.""" + return round((24 / math.pi) * self.sunset_hour_angle(), 1) + + # Rs + def solar_radiation(self): + """ + Calculate solar radiation. + + If *radiation_s* is logged, returns the value as is. If *sunshine_hours* + attribute of the day class is set returns solar radiation amount in mJ/m2/day. + To convert this number to W/m2 multiply multiply it by 11.57 or divide by + 0.0864. Uses Angstrom equation (Eq. 35). + + If climate data is available, and climate is *island* location and + station elevation is between 0-100m it uses simplified (Eq. 51). This + equation does not use temperature data, but just solar radiation and a + coefficient. + + If station elevation is higher than 100m and/or location is not island + it uses (Eq. 50) that calculates solar radiation by using temperature + data along with a *krs* constant. + + If climate is not known it assumes **n=N**, meaning daily sunshine hours + is the same as daylight hours for the given day and location. + """ + if self.radiation_s: + # We need to make sure that solar radiation if set, is not + # larger than clear-sky solar radiation + if CHECK_RADIATION_RANGE and (self.radiation_s > self.clear_sky_rad()): + raise ValueError(f"Solar radiation out ot range. Rso={str(self.clear_sky_rad())}") + return self.radiation_s + + n = self.sunshine_hours + if n is None: + # if we are in island location we refer to equation 51 in UAN-FAO + # paper 56 + if ( + self.station.climate + and self.station.climate.island_location + and self.station.altitude < 100 + ): + ra = self.r_a() + return round((0.7 * ra) - 4, 1) + + if self.station.climate and self.temp_min and self.temp_max: + # We assume caller has only temperature informaiton, and no + # information on overcast conditions. So we resort to Hargreaves + # and Samani's radiation formula: + climate = self.station.climate + ra = self.r_a() + krs = 0.16 + if climate.coastal_location: + krs = 0.19 + elif climate.interior_location: + krs = 0.16 + + return round(krs * math.sqrt(self.temp_max - self.temp_min) * ra, 1) + + else: + n = self.daylight_hours() + + if n and not isinstance(n, (int, float)): + raise TypeError("'n' must be a number") + + if n < 0: + raise ValueError("Observed daylight hours cannot be less than 0") + + # n cannot be more than N, which is available daylight hours + if (n > self.daylight_hours()) and CHECK_SUNSHINE_HOURS_RANGE: + raise ValueError("Daylight hours out of range") + + a_s = 0.25 + b_s = 0.50 + max_daylight_hours = self.daylight_hours() # this is the maximum possible sunshine duration + ra = self.r_a() + + return round((a_s + b_s * n / max_daylight_hours) * ra, 1) + + def solar_radiation_in_mm(self): + """ + Convert solar radiation to mm. + + Alias to *solar_radiation(n)* but converts the output to mm equivalent, + rounded to 1 decimal. + """ + rs = self.solar_radiation() + return round(rs * 0.408, 1) + + # clear-skype solar radiation + def clear_sky_rad(self): + """ + Calculate clear sky radiation when n=N. + + Uses (Eq. 36) for elevations below 100m. Above 100m uses (Eq. 37). + """ + if self.station.altitude < 100: + return round((0.25 + 0.50) * self.r_a(), 1) + + return round((0.75 + (2 * 10 ** (-5)) * self.station.altitude) * self.r_a(), 1) + + def net_solar_rad(self): + """ + Calculate net solar or net shortwave radiation. + + Uses Crop's albedo in calculations. (Eq. 38). + Return radiation in MJ/m2/day + """ + ref_crop = self.station.ref_crop + return round((1 - ref_crop.albedo) * self.solar_radiation(), 1) + + def net_longwave_rad(self): + """Calculate net longwave radiation (Eq. 39).""" + if not (self.temp_max and self.temp_min): + raise AttributeError( + "Net longwave radiation cannot be calculated without min/max temperature" + ) + + temp_max_k = self.temp_max + 273.16 + temp_min_k = self.temp_min + 273.16 + ea = self.actual_vapour_pressure() + rs = self.solar_radiation() + rso = self.clear_sky_rad() + + sb_constant = self.stephan_boltzmann_constant + return round( + sb_constant + * ((temp_max_k**4 + temp_min_k**4) / 2) + * (0.34 - 0.14 * math.sqrt(ea)) + * (1.35 * (rs / rso) - 0.35), + 1, + ) + + def net_radiation(self): + """Calculate net radiation (Eq. 40).""" + ns = self.net_solar_rad() + + try: + nl = self.net_longwave_rad() + except BaseException as e: + raise (str(e)) from e + + if (ns is not None) and (nl is not None): + return round(ns - nl, 1) + + def net_radiation_to_mm(self): + """Convert *net_radiation()* to mm.""" + if net_radiation := self.net_radiation(): + return round(net_radiation * 0.408, 1) + + def relative_humidity(self, T): + """Calculate relative humidity of air given temperature using vapour pressure.""" + if not isinstance(T, (int, float)): + raise TypeError("Number is expected") + return round(100 * (self.actual_vapour_pressure() / self.saturation_vapour_pressure(T)), 3) + + def relative_humidity_mean(self): + """ + Calculate the mean of relative humidity. + + If possible returns mean relative humidity for the day. + If humidity_mean was logged in the station returns it as is. + If min/max humidity values were logged at the station, computes the mean of the two values + If none was logged, but min/max temperature values were logged it attempts to calculate + relative humidity through saturation vapour pressure: + """ + if self.humidity_mean is not None: + return self.humidity_mean + + if (self.humidity_min is not None) and (self.humidity_max is not None): + return round((self.humidity_max + self.humidity_min) / 2, 0) + + if (self.temp_min is not None) and (self.temp_max is not None): + return int( + round( + (self.relative_humidity(self.temp_min) + self.relative_humidity(self.temp_max)) + / 2, + 0, + ) + ) + + return None + + def soil_heat_flux(self): + """Calculate soil heat flux. + + Returns 0.00 (daily coefficient) + """ + return 0.00 + + def eto_hargreaves(self): + """ + Eto estimating using Hargreaves equation. + + If wind and humidty information is available, or can be estimated by this equation + is not recommended. ( Eq. 52 ) + """ + temp_mean = self.interpolate_temp_mean() + return round( + 0.0023 + * (temp_mean + 17.8) + * (self.temp_max - self.temp_min) ** 0.5 + * self.r_a(), + 2, + ) + + def eto(self): + """Calculate evapotranspiration from Eq. 6.""" + # if we cannot get wind speed data we revert to Hargreaves formula. + # Which is not ideal! This can happen only if user removed default 'climate' + # reference + if not self.wind_speed_2m(): + return self.eto_hargreaves() + + if self.interpolate_temp_mean() is None: + raise AttributeError("Cannot calculate eto(): temp_mean (mean temperature) is missing") + + try: + net_radiation = self.net_radiation() + except Exception as e: + raise (str(e)) from e + + temp_mean = self.interpolate_temp_mean() + slope_of_vp = self.slope_of_saturation_vapour_pressure(temp_mean) + ground_heat_flux = self.soil_heat_flux() + u2m = self.wind_speed_2m() + eto_nominator = ( + 0.408 * slope_of_vp * (net_radiation - ground_heat_flux) + ) + self.psychrometric_constant() * ( + 900 / (temp_mean + 273) + ) * u2m * self.vapour_pressure_deficit() + + eto_denominator = slope_of_vp + self.psychrometric_constant() * (1 + 0.34 * u2m) + return round(eto_nominator / eto_denominator, 2) + + +class Climate: + """ + Represents a default climate according to *UN-FAO Paper 56*. + + If module has to make any assumptions regarding the climate it consults + this class for any clues. If you wish to not use any assumptions and + rely soleley on logged station data (if such is available) you may set + Station's *climate* attribute to *None*. + + station = Station(latitude=-20.5, altitude=200) + station.climate = None + + If you want to set a new climate: + + humid_climate = Climate().humid().coastal().moderate_winds() + + station = Station(latitude=-20.5, altitude=200) + station.climate = humid_climate + """ + + def __init__(self): + """ + Instantiate Climate-class object. + + Accepts no arguments. Default initialization is as follows: + + - interior_location + - arid_climate + - dew_point_difference = 2 + - average_wind_speed = 2.0 m/s + - k_rs = 0.16 + + + To affect these values use respected methods documented below. + """ + self.interior_location = True + self.coastal_location = False + self.island_location = False + + self.arid_climate = True + self.humid_climate = False + + # Assining default values for the climatic condition to be able to + # calculate missing data accurately + self.dew_point_difference = 2 + self.average_wind_speed = 2.0 + self.k_rs = 0.16 + + def light_winds(self): + """Set *average_wind_speed* to 0.5m/s.""" + self.average_wind_speed = 0.5 + return self + + def moderate_winds(self): + """Set *average_wind_speed* to 2.0m/s.""" + self.average_wind_speed = 2 + return self + + def strong_winds(self): + """Set average_wind_speed to 4m/s.""" + self.average_wind_speed = 4 + return self + + def arid(self): + """Set *arid_climate* and *dew_point_difference* to 2.""" + self.arid_climate = True + self.humid_climate = False + self.dew_point_difference = 2 + return self + + def humid(self): + """Set *humid_climate* and *dew_point_difference* to 0.""" + self.arid_climate = False + self.humid_climate = True + self.dew_point_difference = 0 + return self + + def interior(self): + """Set *interior_location*, *k_rs* coefficient to *0.16*.""" + self.interior_location = True + self.coastal_location = False + self.island_location = False + self.k_rs = 0.16 + return self + + def coastal(self): + """Set *coastal_location*, *k_rs* to *0.19*.""" + self.interior_location = False + self.coastal_location = True + self.island_location = False + self.k_rs = 0.19 + return self + + def island(self): + """Set *island_location* and *k_rs* to 0.19.""" + self.interior_location = False + self.coastal_location = False + self.island_location = True + self.k_rs = 0.19 + return self + + """ + def set_location(self, location): + if location == "coastal": + return self.coastal() + + if location == "interior": + return self.interior() + + if location == "island": + return self.island() + + return self + """ + + def __str__(self): + """Get string representation of object.""" + location = "interior" if self.interior_location else "coastal" + if self.arid_climate: + humidity = "arid, or semi-arid region" + else: + humidity = "humid, or sub-humid region" + + return f""" + Climate can be described as located in {location} area, {humidity}. + Average wind speeds are estimated at {self.average_wind_speed}. + Dew point temperature is usually runs {self.dew_point_difference}C lower than + day's minimal observed temperature + """ + + +class Crop: + """Represents reference crop as assumed by Penman-Monteith equation.""" + + def __init__(self, resistance_a=208, albedo=0.23, height=0.12): + """ + Instantiate crop class. + + Contains the following assumptions (default values): + - resistance_a = 208 (aerodynamic resistance) + - albedo = 0.23 + - height - 0.12 + """ + self.resistance_a = resistance_a + self.albedo = albedo + self.height = height + + +class StationDay(DayEntry): + """Left here for backwards-compatability with earlier versions of the library.""" + + pass + + +class MonthEntry: + def __init__(self): + pass + + +class WeekEntry: + def __init__(self): + pass + + +class HourEntry: + def __init__(self): + pass diff --git a/pyproject.toml b/pyproject.toml index e74c9ec..f0525a1 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -1,6 +1,129 @@ +[tool.poetry] +name = "penmon" +version = "0.5.0" +description = "Implements weather station class in Python that calculates ETo (reference crop's evapotranspiration) based on UN-FAO Irrigation and Drainage Paper 56" +authors = ["Sherzod RUZMETOV "] +license = "MIT" +readme = "README.md" + +[tool.poetry.dependencies] +python = ">=3.7" + +[tool.poetry.dev-dependencies] +black = "^22.8.0" +flake8 = "^5.0.4" +pylint = "^2.15.2" + [build-system] -requires = [ - "setuptools>=42", - "wheel" -] -build-backend = "setuptools.build_meta" +requires = ["poetry-core"] +build-backend = "poetry.core.masonry.api" + +[tool.mypy] +python_executable = "python" + +[tool.pylint.main] +fail-under = 10 +ignore = ["CVS"] +ignore-patterns = ["^\\.#"] +jobs = 1 + +limit-inference-results = 100 +persistent = true + +py-version = "3.7" + +suggestion-mode = true + + +[tool.pylint.basic] +argument-naming-style = "snake_case" + +attr-naming-style = "snake_case" + +bad-names = ["foo", "bar", "baz", "toto", "tutu", "tata"] +class-attribute-naming-style = "any" +class-const-naming-style = "UPPER_CASE" +class-naming-style = "PascalCase" +class-rgx = "[A-Z][a-z]+" +const-naming-style = "UPPER_CASE" +docstring-min-length = -1 +function-naming-style = "snake_case" +good-names = ["ax", "ds", "df", "i", "j", "k", "ex", "Run", "_"] +inlinevar-naming-style = "any" +method-naming-style = "snake_case" +module-naming-style = "snake_case" +no-docstring-rgx = "^_" +property-classes = ["abc.abstractproperty"] +variable-naming-style = "snake_case" + +[tool.pylint.classes] +defining-attr-methods = ["__init__", "__new__", "setUp", "__post_init__"] +exclude-protected = ["_asdict", "_fields", "_replace", "_source", "_make"] +valid-classmethod-first-arg = ["cls"] +valid-metaclass-classmethod-first-arg = ["cls"] + +[tool.pylint.design] +# Maximum number of arguments for function / method. +max-args = 5 +# Maximum number of attributes for a class (see R0902). +max-attributes = 7 +# Maximum number of boolean expressions in an if statement (see R0916). +max-bool-expr = 5 +# Maximum number of branch for function / method body. +max-branches = 12 +# Maximum number of locals for function / method body. +max-locals = 15 +# Maximum number of parents for a class (see R0901). +max-parents = 7 +# Maximum number of public methods for a class (see R0904). +max-public-methods = 20 +# Maximum number of return / yield for function / method body. +max-returns = 6 +# Maximum number of statements in function / method body. +max-statements = 50 +# Minimum number of public methods for a class (see R0903). +min-public-methods = 2 + +[tool.pylint.exceptions] +# Exceptions that will emit a warning when caught. +overgeneral-exceptions = ["BaseException", "Exception"] + +[tool.pylint.format] +# Regexp for a line that is allowed to be longer than the limit. +ignore-long-lines = "^\\s*(# )??$" +indent-after-paren = 4 +indent-string = " " +max-line-length = 100 +max-module-lines = 1000 + +[tool.pylint.logging] +logging-modules = ["logging"] + +[tool.pylint.miscellaneous] +# List of note tags to take in consideration, separated by a comma. +notes = ["FIXME", "TODO", "NOTE"] + +[tool.pylint.refactoring] +max-nested-blocks = 5 +never-returning-functions = ["sys.exit", "argparse.parse_error"] + +[tool.pylint.reports] +evaluation = "max(0, 0 if fatal else 10.0 - ((float(5 * error + warning + refactor + convention) / statement) * 10))" +score = true + +[tool.pylint.similarities] +ignore-comments = true +ignore-docstrings = true +ignore-imports = true +ignore-signatures = true +min-similarity-lines = 4 + +[tool.flake8] +max-line-length = 100 +max-complexity = 10 + +[tool.black] +line-length = 100 + +[tool.ruff] +line-length = 100 diff --git a/setup.cfg b/setup.cfg index ce5894b..c1821fb 100644 --- a/setup.cfg +++ b/setup.cfg @@ -1,6 +1,6 @@ [metadata] name = penmon -version = 1.5 +version = 0.5.0 author = Sherzod RUZMETOV author_email = sherzodr@gmail.com description = Implementation of weather station that calculates daily ETo for a reference crop using Penman-Monteith equation @@ -17,7 +17,13 @@ classifiers = Topic :: Scientific/Engineering :: Hydrology [options] -python_requires = >=3.6 +platforms = any +python_requires = >=3.7 +packages = find: +include_package_data = True [options.packages.find] -where = src +where = penmon + +[bdist_wheel] +universal = 1 diff --git a/src/.gitignore b/src/.gitignore deleted file mode 100644 index 87c0554..0000000 --- a/src/.gitignore +++ /dev/null @@ -1 +0,0 @@ -/penmon.egg-info/ diff --git a/src/penmon/.gitignore b/src/penmon/.gitignore deleted file mode 100644 index a348e50..0000000 --- a/src/penmon/.gitignore +++ /dev/null @@ -1 +0,0 @@ -/__pycache__/ diff --git a/src/penmon/__init__.py b/src/penmon/__init__.py deleted file mode 100644 index dfe1304..0000000 --- a/src/penmon/__init__.py +++ /dev/null @@ -1,2 +0,0 @@ -from penmon.eto import Station, DayEntry, Climate, Crop -from penmon.eto import CHECK_RADIATION_RANGE, CHECK_SUNSHINE_HOURS_RANGE \ No newline at end of file diff --git a/tests/__init__.py b/tests/__init__.py index d77cfe5..8876372 100644 --- a/tests/__init__.py +++ b/tests/__init__.py @@ -2,4 +2,4 @@ import os import sys -sys.path.insert(0, os.path.abspath( os.path.join(os.path.dirname(__file__), '../src/') )) \ No newline at end of file +sys.path.insert(0, os.path.abspath(os.path.join(os.path.dirname(__file__), "../src/"))) diff --git a/tests/bug_rs_3_20.py b/tests/bug_rs_3_20.py index 675efc3..3d2f45b 100644 --- a/tests/bug_rs_3_20.py +++ b/tests/bug_rs_3_20.py @@ -1,36 +1,35 @@ -''' - - -Problem dataset #1 - -date 1/11/1979 -latitude 37.3113 -altitude 263 -temp_max 3.1 -temp_min 1.247 -wind_speed 2.13654 -humidity_mean 0.853413 -solar_radiation 3.21356 - -''' -import unittest -import penmon as pm - -class Test(unittest.TestCase): - - - def test_bug(self): - station = pm.Station(37.31, 263) - day = station.day_entry("1979-01-11") - day.temp_min=1.247 - day.temp_max=3.10 - day.humidity_mean=85.34 - day.wind_speed=2.13654 - day.radiation_s=3.21356 - self.assertTrue(day.net_radiation(), day.net_radiation()) - self.assertTrue(day.eto(), day.eto()) - - -if __name__ == "__main__": - #import sys;sys.argv = ['', 'Test.testName'] - unittest.main() +""" + + +Problem dataset #1 + +date 1/11/1979 +latitude 37.3113 +altitude 263 +temp_max 3.1 +temp_min 1.247 +wind_speed 2.13654 +humidity_mean 0.853413 +solar_radiation 3.21356 + +""" +from penmon import eto as pm +import unittest + + +class Test(unittest.TestCase): + def test_bug(self): + station = pm.Station(37.31, 263) + day = station.day_entry("1979-01-11") + day.temp_min = 1.247 + day.temp_max = 3.10 + day.humidity_mean = 85.34 + day.wind_speed = 2.13654 + day.radiation_s = 3.21356 + self.assertTrue(day.net_radiation(), day.net_radiation()) + self.assertTrue(day.eto(), day.eto()) + + +if __name__ == "__main__": + # import sys;sys.argv = ['', 'Test.testName'] + unittest.main() diff --git a/tests/error_handling.py b/tests/error_handling.py index d06a5ec..02c67bb 100644 --- a/tests/error_handling.py +++ b/tests/error_handling.py @@ -1,74 +1,75 @@ -''' -Created on Sep 30, 2020 - -@author: sherzodr -''' -import unittest -import penmon as pm - -class Test(unittest.TestCase): - - def test_smoke(self): - station = pm.Station(latitude=41.42, altitude=109) - self.assertIsInstance(station, pm.Station, "Smoke test passed") - - def test_type_error(self): - try: - pm.Station(latitude=10, altitude=109) - except TypeError: - self.assertTrue(True, "Exception was expected and raised") - else: - self.assertTrue(False, "Exception was expected but was NOT raised") - - def test_latitude_range_error(self): - try: - pm.Station(latitude=100.0, altitude=100) - except: - self.assertTrue(True, "Exception was expected and raised") - else: - self.assertTrue(False, "Exception was expected but was NOT raised") - - def test_altitude_range_error(self): - try: - pm.Station(latitude=41.42, altitude=-1) - except: - self.assertTrue(True, "Exception was expected and raised") - else: - self.assertTrue(False, "Exception was expected but NOT raised") - - def test_day_number_type(self): - station = pm.Station(latitude=41.42, altitude=109) - try: - station.day_entry(365.0) - except TypeError: - self.assertTrue(True, "Exception was expected and raised") - else: - self.assertTrue(False, "Exception was expected but was NOT raised") - - def test_day_number_range(self): - station = pm.Station(latitude=41.42, altitude=109) - - try: - station.day_entry(367) - except: - self.assertTrue(True, "Exception was expected and raised") - else: - self.assertTrue(False, "Exception was expected but was NOT raised") - - def test_immature_eto(self): - station = pm.Station(41.42, 109) - day = station.day_entry(238, temp_mean=25.00) - self.assertTrue(day.temp_mean != None, "temp_mean was set") - self.assertEqual(day.temp_mean, day.Tmean()) - - # following code should raise an exception: - try: - day.eto() - except Exception as e: - self.assertTrue(True, str(e)) - else: - self.assertTrue(False, "Exception was expected") - -if __name__ == "__main__": - # import sys;sys.argv = ['', 'Test.testName'] - unittest.main() +""" +Created on Sep 30, 2020 + +@author: sherzodr +""" +from penmon import eto as pm +import unittest + + +class Test(unittest.TestCase): + def test_smoke(self): + station = pm.Station(latitude=41.42, altitude=109) + self.assertIsInstance(station, pm.Station, "Smoke test passed") + + def test_type_error(self): + try: + pm.Station(latitude=10, altitude=109) + except TypeError: + self.assertTrue(True, "Exception was expected and raised") + else: + self.assertTrue(False, "Exception was expected but was NOT raised") + + def test_latitude_range_error(self): + try: + pm.Station(latitude=100.0, altitude=100) + except: + self.assertTrue(True, "Exception was expected and raised") + else: + self.assertTrue(False, "Exception was expected but was NOT raised") + + def test_altitude_range_error(self): + try: + pm.Station(latitude=41.42, altitude=-1) + except: + self.assertTrue(True, "Exception was expected and raised") + else: + self.assertTrue(False, "Exception was expected but NOT raised") + + def test_day_number_type(self): + station = pm.Station(latitude=41.42, altitude=109) + try: + station.day_entry(365.0) + except TypeError: + self.assertTrue(True, "Exception was expected and raised") + else: + self.assertTrue(False, "Exception was expected but was NOT raised") + + def test_day_number_range(self): + station = pm.Station(latitude=41.42, altitude=109) + + try: + station.day_entry(367) + except: + self.assertTrue(True, "Exception was expected and raised") + else: + self.assertTrue(False, "Exception was expected but was NOT raised") + + def test_immature_eto(self): + station = pm.Station(41.42, 109) + day = station.day_entry(238, temp_mean=25.00) + self.assertTrue(day.temp_mean != None, "temp_mean was set") + self.assertEqual(day.temp_mean, day.interpolate_temp_mean()) + + # following code should raise an exception: + try: + day.eto() + except Exception as e: + self.assertTrue(True, str(e)) + else: + self.assertTrue(False, "Exception was expected") + + +if __name__ == "__main__": + # import sys;sys.argv = ['', 'Test.testName'] + unittest.main() diff --git a/tests/et0_calc.py b/tests/et0_calc.py index 99f74fa..7db3aba 100644 --- a/tests/et0_calc.py +++ b/tests/et0_calc.py @@ -1,15 +1,17 @@ -import unittest, penmon as pm - -class Test(unittest.TestCase): - - def test_daylight_hours(self): - station = pm.Station(41.42, 109) - day = station.day_entry(135) - day.temp_min = 19.5 - day.temp_max = 28 - - self.assertEqual(day.daylight_hours(), 14.3, "daylighth_hours") - -if __name__ == "__main__": - unittest.main() - \ No newline at end of file +import unittest + +from penmon import eto as pm + + +class Test(unittest.TestCase): + def test_daylight_hours(self): + station = pm.Station(41.42, 109) + day = station.day_entry(135) + day.temp_min = 19.5 + day.temp_max = 28 + + self.assertEqual(day.daylight_hours(), 14.3, "daylighth_hours") + + +if __name__ == "__main__": + unittest.main() diff --git a/tests/fao56.py b/tests/fao56.py index 1f512fc..0877325 100644 --- a/tests/fao56.py +++ b/tests/fao56.py @@ -1,4 +1,4 @@ -''' +""" Created on Sep 29, 2020 @author: sherzodr @@ -9,9 +9,10 @@ Examples in the original paper does not use floating points as consistently. I had to standardadize floating point arithmetic for the purposes of consistency. -''' +""" +from penmon.eto import Climate, DayEntry, Station import unittest -from penmon import Station, Climate, DayEntry + class Test(unittest.TestCase): @@ -31,20 +32,19 @@ def test_station_day(self): """Testing station's 'day()' method""" day = self.station.day_entry(130) self.assertEqual(day.day_number, 130, "Station.day.day_number") - #print(type(day)) - self.assertIsInstance(day, DayEntry, "staion.day is instance of " ) + # print(type(day)) + self.assertIsInstance(day, DayEntry, "staion.day is instance of ") def test_atmospheric_pressure(self): station = Station(41.42, 1800) station.climate = Climate() day = station.day_entry(130) self.assertEqual(round(day.atmospheric_pressure(), 1), 81.8, "Atmosphertic pressure") - + def test_station_atmospheric_pressure(self): station = Station(41.42, 1800) self.assertEqual(round(station.atmospheric_pressure(), 1), 81.8, "Atmosphertic pressure") - def test_psychrometric_constant(self): station = Station(41.42, 1800) station.climate = Climate() @@ -53,33 +53,45 @@ def test_psychrometric_constant(self): def test_specific_heat(self): day = self.station.day_entry(130) - self.assertEqual(day.specific_heat(), 1.013 * 10 ** (-3)) + self.assertEqual(day.specific_heat, 1.013 * 10 ** (-3)) def test_latent_heat_of_vaporization(self): day = self.station.day_entry(130) - self.assertEqual(day.latent_heat_of_vaporization(), 2.45) + self.assertEqual(day.latent_heat_of_vaporization, 2.45) def test_mean_saturation_vp(self): day_130 = self.day_130 day_130.temp_max = 24.5 day_130.temp_min = 15 - self.assertEqual(day_130.mean_saturation_vapour_pressure(), 2.39, "mean_saturation_vapour_pressure()") - + self.assertEqual( + day_130.mean_saturation_vapour_pressure(), + 2.39, + "mean_saturation_vapour_pressure()", + ) + def test_RH(self): - day = self.station.day_entry(130) + day = self.station.day_entry(130) day.temp_dew = 19.5 - self.assertEqual(round(day.RH(35), 2), 40.32, "We can calculate relative humidity for a given T if dew point is known") + self.assertEqual( + round(day.relative_humidity(35), 2), + 40.32, + "We can calculate relative humidity for a given T if dew point is known", + ) def test_mean_saturation_vp2(self): day = self.station.day_entry(130) day.temp_mean = 19.75 - self.assertEqual(day.mean_saturation_vapour_pressure(), 2.302, "mean_saturation_vapour_pressure using just Tmean") - + self.assertEqual( + day.mean_saturation_vapour_pressure(), + 2.302, + "mean_saturation_vapour_pressure using just Tmean", + ) + def test_slope_of_saturation_vp(self): day = self.station.day_entry(130) slope = day.slope_of_saturation_vapour_pressure(24.5) self.assertEqual(slope, 0.183837, "slope of vapour pressure curve") - + def test_actual_vapour_pressure_psychrometric(self): station = Station(41.42, 1200) station.climate = Climate() @@ -89,7 +101,11 @@ def test_actual_vapour_pressure_psychrometric(self): self.assertEqual(day.atmospheric_pressure(), 87.9, 87.9) self.assertEqual(day.saturation_vapour_pressure(19.5), 2.267) self.assertEqual(day.temp_dew, None, "No dew point temperature known") - self.assertEqual(day.actual_vapour_pressure(), 1.91, "Actual vapour pressure from psychrometric data") + self.assertEqual( + day.actual_vapour_pressure(), + 1.91, + "Actual vapour pressure from psychrometric data", + ) def test_actual_vapour_pressure_dew(self): station = Station(41.42, 1200) @@ -99,7 +115,11 @@ def test_actual_vapour_pressure_dew(self): self.assertEqual(day.temp_dew, None, "Dew temp is not known") day.temp_dew = 17.0 self.assertEqual(day.temp_dew, 17.0, "Dew temp is set") - self.assertEqual(day.actual_vapour_pressure(), 1.938, "Actual vapour pressure with no psychrometric data") + self.assertEqual( + day.actual_vapour_pressure(), + 1.938, + "Actual vapour pressure with no psychrometric data", + ) def test_actual_vapour_pressure_humidity(self): @@ -117,10 +137,10 @@ def test_actual_vapour_pressure_humidity(self): sat_vap_pres_max = day.saturation_vapour_pressure(day.temp_max) self.assertEqual(sat_vap_pres_max, 3.168) - + actual_vp = day.actual_vapour_pressure() self.assertEqual(actual_vp, 1.702) - + def test_actual_vapour_pressure_humidity_max_temp_min(self): station = Station(41.42, 1200) @@ -129,7 +149,11 @@ def test_actual_vapour_pressure_humidity_max_temp_min(self): day.humidity_max = 82 actual_vp = day.actual_vapour_pressure() - self.assertEqual(actual_vp, 1.692, "Actual vapour pressure with only humidity_max and temp_min") + self.assertEqual( + actual_vp, + 1.692, + "Actual vapour pressure with only humidity_max and temp_min", + ) def test_actual_vapour_pressure_humidity_mean(self): station = Station(41.42, 1200) @@ -140,11 +164,15 @@ def test_actual_vapour_pressure_humidity_mean(self): day.humidity_mean = 68 actual_vp = day.actual_vapour_pressure() - self.assertEqual(actual_vp, 1.779, "Actual vapour pressure with only humidity_mean, temp_min and temp_max") - + self.assertEqual( + actual_vp, + 1.779, + "Actual vapour pressure with only humidity_mean, temp_min and temp_max", + ) + def test_vapour_pressure_deficit(self): station = Station(41.42, 1200) - + day = station.day_entry(130) day.temp_min = 18 day.temp_max = 25 @@ -163,7 +191,7 @@ def test_latitude_rad(self): def test_relative_sun_distance(self): dr = self.station.day_entry(246).relative_sun_distance() self.assertEqual(dr, 0.985) - + def test_solar_declination(self): declination = self.station.day_entry(246).solar_declination() self.assertEqual(declination, 0.120) @@ -174,11 +202,11 @@ def test_sunset_hour_angle(self): self.assertEqual(sunset_angle, 1.527) def test_r_a(self): - r_a = Station(-20.0, 1200).day_entry(246).R_a() + r_a = Station(-20.0, 1200).day_entry(246).r_a() self.assertEqual(r_a, 32.2) def test_r_a_in_mm(self): - r_a = Station(-20.0, 1200).day_entry(246).R_a_in_mm() + r_a = Station(-20.0, 1200).day_entry(246).ra_to_mm() self.assertEqual(r_a, 13.10) def test_daylight_hours(self): @@ -186,24 +214,24 @@ def test_daylight_hours(self): self.assertEqual(hours, 11.7) def test_solar_radiation(self): - day=Station(-22.90,1200).day_entry(135) + day = Station(-22.90, 1200).day_entry(135) day.sunshine_hours = 7.10 radiation = day.solar_radiation() self.assertEqual(radiation, 14.4) - + def test_clear_sky_solar_radiation(self): day = Station(-22.90, 0).day_entry(135) - day.sunshine_hours=7.1 + day.sunshine_hours = 7.1 solar_radiation = day.solar_radiation() self.assertEqual(solar_radiation, 14.4) - - clear_sky_radiation = day.R_so() + + clear_sky_radiation = day.clear_sky_rad() self.assertEqual(clear_sky_radiation, 18.8) - + def test_net_shortwave_radiation(self): day = Station(-22.90, 1200).day_entry(135) - day.sunshine_hours=7.1 - r_ns = day.R_ns() + day.sunshine_hours = 7.1 + r_ns = day.net_solar_rad() self.assertEqual(r_ns, 11.1) def test_net_longwave_radiation(self): @@ -215,10 +243,10 @@ def test_net_longwave_radiation(self): vp = day.actual_vapour_pressure() self.assertEqual(vp, 2.1) - - day.sunshine_hours = 7.1; - r_nl = day.R_nl() + day.sunshine_hours = 7.1 + + r_nl = day.net_longwave_rad() self.assertEqual(r_nl, 3.3) def test_net_radiation(self): @@ -226,18 +254,18 @@ def test_net_radiation(self): day.temp_max = 25.1 day.temp_min = 19.1 day.vapour_pressure = 2.1 - day.sunshine_hours=7.1 + day.sunshine_hours = 7.1 net_radiation = day.net_radiation() self.assertEqual(net_radiation, 7.8) - + def test_wind_speed2m(self): day = Station(-22.90, 1200).day_entry(135) - day.wind_speed = 5 - + day.wind_speed = 5 + self.assertEqual(day.wind_speed_2m(), 5) - day.station.anemometer_height= 10 + day.station.anemometer_height = 10 day.wind_speed = 3.2 self.assertEqual(day.wind_speed_2m(), 2.4) @@ -245,11 +273,11 @@ def test_solar_radiation_from_temp(self): day = Station(45.72, 200).day_entry(196) day.temp_max = 26.6 day.temp_min = 14.8 - - ra = day.R_a() + + ra = day.r_a() self.assertEqual(ra, 40.6) - solar_radiation = day.solar_radiation() + solar_radiation = day.solar_radiation() self.assertEqual(solar_radiation, 22.3) self.assertEqual(day.solar_radiation_in_mm(), 9.1) @@ -261,19 +289,19 @@ def test_net_radiation_without_radiation_data(self): day.temp_min = 25.6 day.temp_max = 34.8 - + self.assertEqual(Station(13.73, 2).latitude_rad, 0.24) - ra = day.R_a() + ra = day.r_a() self.assertEqual(ra, 38.0) net_radiation = day.net_radiation() self.assertEqual(net_radiation, 14.0) - + def test_solar_radiation_in_island(self): day = Station(41.42, 10).day_entry(105) day.station.climate.island() - + rs = day.solar_radiation() self.assertEqual(rs, 20.0) @@ -289,9 +317,9 @@ def test_solar_radiation_in_island(self): def test_eto_hargreaves(self): day = Station(41.42, 109).day_entry(295) - day.temp_min = 19.5; + day.temp_min = 19.5 day.temp_max = 26.5 - + eto = day.eto_hargreaves() self.assertEqual(eto, 4.97) @@ -300,7 +328,7 @@ def test_eto(self): day.temp_min = 19.5 day.temp_max = 36.5 day.wind_speed = 2 - #day.humidity_mean = 60 + # day.humidity_mean = 60 self.assertEqual(day.slope_of_saturation_vapour_pressure(23), 0.169921) self.assertEqual(day.net_radiation(), 16.1) @@ -313,6 +341,7 @@ def test_eto(self): self.assertEqual(eto, 6.98) + if __name__ == "__main__": - #import sys;sys.argv = ['', 'Test.testName'] - unittest.main() \ No newline at end of file + # import sys;sys.argv = ['', 'Test.testName'] + unittest.main() diff --git a/tests/get_day_new_syntax.py b/tests/get_day_new_syntax.py index 221644f..83e1036 100644 --- a/tests/get_day_new_syntax.py +++ b/tests/get_day_new_syntax.py @@ -1,31 +1,34 @@ -''' -Created on Oct 5, 2020 - -@author: sherzodr -''' -import unittest -import penmon as pm - -class Test(unittest.TestCase): - - station=pm.Station(41.42, 109) - - def test_name(self): - #day_228 = self.station.day_entry(228, temp_min=19.5, temp_max=25.6, humidity_mean=60, wind_speed=2.35) - - day_228 = self.station.day_entry(228) - day_228.temp_min=19.5 - day_228.temp_max=25.6 - day_228.humidity_mean=60 - day_228.wind_speed=2.35 - - self.assertEqual(day_228.eto(), 3.94) - - - def test_new_day_entry(self): - day_228_2=self.station.day_entry(228, temp_min=19.5, temp_max=25.6, humidity_mean=60, wind_speed=2.35) - self.assertEqual(day_228_2.eto(), 3.94) - -if __name__ == "__main__": - #import sys;sys.argv = ['', 'Test.testName'] - unittest.main() \ No newline at end of file +""" +Created on Oct 5, 2020 + +@author: sherzodr +""" +from penmon import eto as pm +import unittest + + +class Test(unittest.TestCase): + + station = pm.Station(41.42, 109) + + def test_name(self): + # day_228 = self.station.day_entry(228, temp_min=19.5, temp_max=25.6, humidity_mean=60, wind_speed=2.35) + + day_228 = self.station.day_entry(228) + day_228.temp_min = 19.5 + day_228.temp_max = 25.6 + day_228.humidity_mean = 60 + day_228.wind_speed = 2.35 + + self.assertEqual(day_228.eto(), 3.94) + + def test_new_day_entry(self): + day_228_2 = self.station.day_entry( + 228, temp_min=19.5, temp_max=25.6, humidity_mean=60, wind_speed=2.35 + ) + self.assertEqual(day_228_2.eto(), 3.94) + + +if __name__ == "__main__": + # import sys;sys.argv = ['', 'Test.testName'] + unittest.main() diff --git a/tests/rh_mean.py b/tests/rh_mean.py index e7971ff..11f1396 100644 --- a/tests/rh_mean.py +++ b/tests/rh_mean.py @@ -1,30 +1,27 @@ -''' +""" Created on Jan 26, 2021 @author: sherzodr -''' +""" +from penmon import eto as pm import unittest -import penmon class Test(unittest.TestCase): - - def testName(self): - station = penmon.Station(41.42, 109) + station = pm.Station(41.42, 109) qarmish = station.day_entry("2020-07-04", temp_min=15, temp_max=30) - rh_mean = qarmish.RH_mean() - self.assertEqual(rh_mean, 62, "RH_mean()") + rh_mean = qarmish.relative_humidity_mean() + self.assertEqual(rh_mean, 62, "relative_humidity_mean() wrong calc.") - humidity_min = qarmish.RH(15) - humidity_max = qarmish.RH(30) - humidity_mean = int(round(( humidity_max + humidity_min ) / 2, 0)) + humidity_min = qarmish.relative_humidity(15) + humidity_max = qarmish.relative_humidity(30) + humidity_mean = int(round((humidity_max + humidity_min) / 2, 0)) self.assertEqual(humidity_mean, 62, "Humidity is 62") - if __name__ == "__main__": - #import sys;sys.argv = ['', 'Test.testName'] - unittest.main() \ No newline at end of file + # import sys;sys.argv = ['', 'Test.testName'] + unittest.main() diff --git a/tests/rh_mean_vs_minmax.py b/tests/rh_mean_vs_minmax.py index 7685740..6369e66 100644 --- a/tests/rh_mean_vs_minmax.py +++ b/tests/rh_mean_vs_minmax.py @@ -1,41 +1,54 @@ -''' +""" Created on Mar 31, 2021 @author: sherzodr -''' -import unittest -import penmon +""" import math +from penmon import eto as pm +import unittest + +# penmon.eto.CHECK_RADIATION_RANGE=False -#penmon.eto.CHECK_RADIATION_RANGE=False class Test(unittest.TestCase): - - station = penmon.Station(37.85, 100, anemometer_height=10) - day180 = station.get_day(180, temp_min=15, temp_max=35, humidity_min=30, humidity_max=70, wind_speed=2, radiation_s=30) + + station = pm.Station(37.85, 100, anemometer_height=10) + day180 = station.get_day( + 180, + temp_min=15, + temp_max=35, + humidity_min=30, + humidity_max=70, + wind_speed=2, + radiation_s=30, + ) def testDecl(self): day180 = self.day180 - self.assertEqual(day180.solar_declination(), 0.405, "solar declination is 0.405 rad") - self.assertEqual(round(day180.solar_declination() * 180 / math.pi, 2), 23.20, "solar declination is 23.23 degrees") - - + self.assertEqual(day180.solar_declination(), 0.405, "solar declination is 0.405 rad") + self.assertEqual( + round(day180.solar_declination() * 180 / math.pi, 2), + 23.20, + "solar declination is 23.23 degrees", + ) + def testHourAngle(self): self.assertEqual(self.day180.sunset_hour_angle(), 1.911, "Sunset hour angle is 1.91") - - def testETo(self): - self.assertEqual(self.day180.eto(), 6.79, "eto equals 6.79mm") - - - def testSunDistance(self): - self.assertEqual(self.day180.relative_sun_distance(), 0.967, "inverse of relative sun distance") + def testETo(self): + self.assertEqual(self.day180.eto(), 6.79, "eto equals 6.79mm") - def testRs(self): - self.assertEqual(self.day180.R_a(), 41.7, "Ra is 41.7") + def testSunDistance(self): + self.assertEqual( + self.day180.relative_sun_distance(), + 0.967, + "inverse of relative sun distance", + ) + def testRs(self): + self.assertEqual(self.day180.r_a(), 41.7, "Ra is 41.7") if __name__ == "__main__": - #import sys;sys.argv = ['', 'Test.testName'] - unittest.main() \ No newline at end of file + # import sys;sys.argv = ['', 'Test.testName'] + unittest.main() diff --git a/tests/test_date_parse.py b/tests/test_date_parse.py index dba16c6..ab9a103 100644 --- a/tests/test_date_parse.py +++ b/tests/test_date_parse.py @@ -1,34 +1,35 @@ -''' -Created on Sep 30, 2020 - -@author: sherzodr -''' -import unittest -import penmon as pm - - -class Test(unittest.TestCase): - - station = pm.Station(41.42, 109) - - def test_date_parse(self): - day = self.station.day_entry(229) - self.assertEqual(day.day_number, 229) - - def test_first_date(self): - day=self.station.day_entry("2020-01-01") - self.assertEqual(day.day_number, 1) - - def test_isodate_parse(self): - - day=self.station.day_entry("2020-08-16") - self.assertEqual(day.day_number, 229) - - def test_parse_us_date(self): - - day=self.station.day_entry("08/16/2020", date_template="%m/%d/%Y") - self.assertEqual(day.day_number, 229) - -if __name__ == "__main__": - #import sys;sys.argv = ['', 'Test.testName'] - unittest.main() \ No newline at end of file +""" +Created on Sep 30, 2020 + +@author: sherzodr +""" +from penmon import eto as pm +import unittest + + +class Test(unittest.TestCase): + + station = pm.Station(41.42, 109) + + def test_date_parse(self): + day = self.station.day_entry(229) + self.assertEqual(day.day_number, 229) + + def test_first_date(self): + day = self.station.day_entry("2020-01-01") + self.assertEqual(day.day_number, 1) + + def test_isodate_parse(self): + + day = self.station.day_entry("2020-08-16") + self.assertEqual(day.day_number, 229) + + def test_parse_us_date(self): + + day = self.station.day_entry("08/16/2020", date_template="%m/%d/%Y") + self.assertEqual(day.day_number, 229) + + +if __name__ == "__main__": + # import sys;sys.argv = ['', 'Test.testName'] + unittest.main() diff --git a/tests/test_solar_range.py b/tests/test_solar_range.py index ef04cfa..80259b1 100644 --- a/tests/test_solar_range.py +++ b/tests/test_solar_range.py @@ -1,78 +1,86 @@ -''' -Created on Oct 6, 2020 - -@author: sherzodr -''' -import unittest -import penmon as pm - - -class Test(unittest.TestCase): - - def test_solar_radiation_range(self): - station = pm.Station(41.42, 109) - - self.assertTrue(pm.CHECK_RADIATION_RANGE == True) - - try: - day = station.day_entry("2020-08-16", temp_min=19.8, temp_max=29.9, radiation_s=27.5) - self.assertTrue(day.eto()) - except Exception as e: - self.assertTrue(False, "Out of range ValueError caught: " + str(e)) - else: - self.assertTrue(True, "Out of range ValueError NOT caught") - - def test_daylight_hours(self): - station = pm.Station(41.42, 109) - - try: - day = station.day_entry("2020-08-16", temp_min=19.8, temp_max=29.9, radiation_s=35) - self.assertTrue(day.eto()) - except ValueError: - self.assertTrue(True, "out of range ValueError caught") - else: - self.assertTrue(False, "out of range ValueError NOT caught") - day = station.day_entry("2020-08-16", temp_min=19.8, temp_max=29.9, radiation_s=27.5) - - def test_no_solar_range_test(self): - - self.assertTrue(pm.eto.CHECK_RADIATION_RANGE, "CHECK_RADIATION_RANGE defaults to True") - - pm.eto.CHECK_RADIATION_RANGE = False - - self.assertTrue(pm.eto.CHECK_RADIATION_RANGE == False, "CHECK_RADIATION_RANGE defaults to True") - - try: - day = pm.Station(41.42, 109).day_entry("2020-08-16", temp_min=19.8, temp_max=29.9, radiation_s=35) - self.assertTrue(day.eto() > 0) - except: - self.assertTrue(False, "out of range exception raised") - else: - self.assertTrue(True, "out of range exception was not raised, as planned") - - pm.CHECK_RADIATION_RANGE = True - - def test_check_sunshine_hours_range(self): - try: - day = pm.Station(41.42, 109).day_entry("2019-12-21", sunshine_hours=15, temp_min=19.8, temp_max=29.9) - self.assertTrue(day.eto()) - except Exception as e: - self.assertTrue(True, str(e)) - else: - self.assertTrue(False, "Out of range exception was supposed to be raised") - - def test_check_sunshine_hours_range_error_ignore(self): - pm.eto.CHECK_SUNSHINE_HOURS_RANGE = False - - try: - day = pm.Station(41.42, 109).day_entry("2019-12-21", temp_min=-5, temp_max=5.50, sunshine_hours=15) - self.assertTrue(day.eto()) - except Exception as e: - self.assertTrue(False, str(e)) - else: - self.assertTrue(True, "Out of range exception was NOT supposed to be raised") - pm.eto.CHECK_SUNSHINE_HOURS_RANGE = True - -if __name__ == "__main__": - # import sys;sys.argv = ['', 'Test.test_solar_radiation_range'] - unittest.main() +""" +Created on Oct 6, 2020 + +@author: sherzodr +""" +import unittest + +from penmon import eto as pm + + +class Test(unittest.TestCase): + def test_solar_radiation_range(self): + station = pm.Station(41.42, 109) + + self.assertTrue(pm.CHECK_RADIATION_RANGE is True) + + eto = False + try: + day = station.day_entry("2020-08-16", temp_min=19.8, temp_max=29.9, radiation_s=27.5) + self.assertTrue(day.eto()) + except Exception as e: + self.assertIsInstance(e, ValueError, f"Out of range ValueError caught: {str(e)}") + else: + self.assertTrue(False, "Out of range ValueError NOT caught") + + def test_daylight_hours(self): + station = pm.Station(41.42, 109) + + try: + day = station.day_entry("2020-08-16", temp_min=19.8, temp_max=29.9, radiation_s=35) + self.assertTrue(day.eto()) + except ValueError: + self.assertTrue(True, "out of range ValueError caught") + else: + self.assertTrue(False, "out of range ValueError NOT caught") + day = station.day_entry("2020-08-16", temp_min=19.8, temp_max=29.9, radiation_s=27.5) + + def test_no_solar_range_test(self): + + self.assertTrue(pm.CHECK_RADIATION_RANGE, "CHECK_RADIATION_RANGE defaults to True") + + pm.CHECK_RADIATION_RANGE = False + + self.assertTrue(not pm.CHECK_RADIATION_RANGE, "CHECK_RADIATION_RANGE defaults to True") + + try: + day = pm.Station(41.42, 109).day_entry( + "2020-08-16", temp_min=19.8, temp_max=29.9, radiation_s=35 + ) + self.assertTrue(day.eto() > 0) + except Exception: + self.assertTrue(False, "out of range exception raised") + else: + self.assertTrue(True, "out of range exception was not raised, as planned") + + pm.CHECK_RADIATION_RANGE = True + + def test_check_sunshine_hours_range(self): + try: + day = pm.Station(41.42, 109).day_entry( + "2019-12-21", sunshine_hours=15, temp_min=19.8, temp_max=29.9 + ) + self.assertTrue(day.eto()) + except Exception as e: + self.assertTrue(True, str(e)) + else: + self.assertTrue(False, "Out of range exception was supposed to be raised") + + def test_check_sunshine_hours_range_error_ignore(self): + pm.CHECK_SUNSHINE_HOURS_RANGE = False + + try: + day = pm.Station(41.42, 109).day_entry( + "2019-12-21", temp_min=-5, temp_max=5.50, sunshine_hours=15 + ) + self.assertTrue(day.eto()) + except Exception as e: + self.assertTrue(False, str(e)) + else: + self.assertTrue(True, "Out of range exception was NOT supposed to be raised") + pm.CHECK_SUNSHINE_HOURS_RANGE = True + + +if __name__ == "__main__": + # import sys;sys.argv = ['', 'Test.test_solar_radiation_range'] + unittest.main() diff --git a/tests/validate_results.py b/tests/validate_results.py index 70aa336..19fd27a 100644 --- a/tests/validate_results.py +++ b/tests/validate_results.py @@ -1,42 +1,47 @@ -import unittest -from penmon import Station - -class Test(unittest.TestCase): - - def test_et0(self): - station = Station(41.42, 109) - - data = { - 228: [19.8, 29.9, 4.49], - 229:[15.1, 29.9, 5.13], - 230:[13.4, 29.6, 5.24], - 231:[12.7, 31.1, 5.62], - 232:[16.2, 32.8, 5.66], - 233:[16.9, 33.6, 5.72], - 234:[16.0, 35.5, 6.27] - } - - for day_number, wdata in data.items(): - day = station.day_entry(day_number) - day.temp_min = wdata[0] - day.temp_max = wdata[1] - day.wind_speed = 2 - self.assertEqual(day.eto(), wdata[2], f"day {day_number}: Tmin={wdata[0]}, Tmax={wdata[1]}, ETo={wdata[2]}") - - def test_with_no_assumptions(self): - - station = Station(41.42, 109) - station.climate = None - - day = station.day_entry(228) # august 16 - day.temp_min = 19.8 - day.temp_max = 29.9 - day.radiation_s = 264 * 0.0864 - day.humidity_min = 37 - day.humidity_max = 88 - day.wind_speed = 2 - - self.assertEqual(day.eto(), 5.3) - -if __name__ == "__main__": - unittest.main() +from penmon.eto import Station +import unittest + + +class Test(unittest.TestCase): + def test_et0(self): + station = Station(41.42, 109) + + data = { + 228: [19.8, 29.9, 4.49], + 229: [15.1, 29.9, 5.13], + 230: [13.4, 29.6, 5.24], + 231: [12.7, 31.1, 5.62], + 232: [16.2, 32.8, 5.66], + 233: [16.9, 33.6, 5.72], + 234: [16.0, 35.5, 6.27], + } + + for day_number, wdata in data.items(): + day = station.day_entry(day_number) + day.temp_min = wdata[0] + day.temp_max = wdata[1] + day.wind_speed = 2 + self.assertEqual( + day.eto(), + wdata[2], + f"day {day_number}: Tmin={wdata[0]}, Tmax={wdata[1]}, ETo={wdata[2]}", + ) + + def test_with_no_assumptions(self): + + station = Station(41.42, 109) + station.climate = None + + day = station.day_entry(228) # august 16 + day.temp_min = 19.8 + day.temp_max = 29.9 + day.radiation_s = 264 * 0.0864 + day.humidity_min = 37 + day.humidity_max = 88 + day.wind_speed = 2 + + self.assertEqual(day.eto(), 5.3) + + +if __name__ == "__main__": + unittest.main()