Meteorology#

Warning

This is not meant to be a standalone notebook. This notebook is part of the process we have for adding entries to the NCL Index and is not meant to be used as tutorial or example code.

Functions covered#

NCL code#

; daylight_fao56
; Adapted from https://www.ncl.ucar.edu/Document/Functions/Crop/daylight_fao56.shtml

; ncl -n daylight_fao56.ncl >> daylight_fao56_output.txt

print("DOY, Latitude (Degrees), Daylight Hours")
do doy=0,365
	do lat=-66,66
		begin
		daylight_hours = daylight_fao56(doy, lat)
		print (doy +","+ lat +","+ daylight_hours)
		end
	end do
end do

; dewtemp_trh
; Adapted from https://www.ncl.ucar.edu/Document/Functions/Built-in/dewtemp_trh.shtml

; ncl -n dewtemp_trh.ncl >> dewtemp_trh_output.txt

print("Temperature (K), Relative Humidity (%), Dew Temperature (C)")
do tk=273,374
	do rh=1,100
		begin
		dewtemp = dewtemp_trh(tk,rh)-273.15
		print (tk +","+ rh +","+ dewtemp)
		end
	end do
end do

; satvpr_temp_fao56
; Adapted from https://www.ncl.ucar.edu/Document/Functions/Crop/satvpr_temp_fao56.shtml

; ncl -n satvpr_temp_fao56.ncl >> satvpr_temp_fao56_output.txt

print("Temperature (F), Saturation Vapor Pressure (kPa)")
do temp=33,212
	begin
		sat_vpr_pressure = satvpr_temp_fao56(temp, (/2, 2/))
		print (temp + "," + sat_vpr_pressure)
	end
end do

; satvpr_tdew_fao56
; Adapted from https://www.ncl.ucar.edu/Document/Functions/Crop/satvpr_tdew_fao56.shtml

; ncl -n satvpr_tdew_fao56.ncl >> satvpr_tdew_fao56_output.txt

print("Temperature (F), Actual Saturation Vapor Pressure (kPa)")
do temp=33,212
	begin
		act_sat_vpr_pressure = satvpr_tdew_fao56(temp, (/2, 2/))
		print (temp + "," + act_sat_vpr_pressure)
	end
end do

Python Functionality#

dewtemp_trh#

#### Collect NCL values for dewtemp_trh from geocat-datafiles
import geocat.datafiles as gdf
import numpy as np

dewtemp_data = gdf.get('applications_files/ncl_outputs/dewtemp_trh_output.txt')
dewtemp_data = np.loadtxt(dewtemp_data, delimiter=',', skiprows=6)
Downloading file 'applications_files/ncl_outputs/dewtemp_trh_output.txt' from 'https://github.com/NCAR/GeoCAT-datafiles/raw/main/applications_files/ncl_outputs/dewtemp_trh_output.txt' to '/home/runner/.cache/geocat'.
### Collect NCL `dewtemp` value and associated (temperature_kelvin, relative humidity) values
ncl_dewtemp = {}
tk_rh = tuple(map(tuple, dewtemp_data[::, 0:2]))
dewtemp_values = dewtemp_data[::, 2]
ncl_dewtemp = dict(zip(tk_rh, dewtemp_values))
### Collect Temperature (Kelvin) and Relative Humidity Pairs
tk_rh = []
for tk in range(273, 374 + 1):
    for rh in range(1, 100 + 1):
        tk_rh.append((tk, rh))

### Calculate GeoCAT-Comp `dewtemp` value and tk/rh
from geocat.comp import dewtemp

geocat_dewtemp = {}

for i, pair in enumerate(tk_rh):
    tk, rh = pair
    geocat_dewtemp[pair] = dewtemp(tk, rh) - 273.15

daylight_fao56#

#### Collect NCL values for daylight_fao56 from geocat-datafiles
import geocat.datafiles as gdf
import numpy as np

daylight_data = gdf.get('applications_files/ncl_outputs/daylight_fao56_output.txt')
daylight_data = np.loadtxt(daylight_data, delimiter=',', skiprows=6)
Downloading file 'applications_files/ncl_outputs/daylight_fao56_output.txt' from 'https://github.com/NCAR/GeoCAT-datafiles/raw/main/applications_files/ncl_outputs/daylight_fao56_output.txt' to '/home/runner/.cache/geocat'.
### Collect NCL `daylight_fao56` value and associated (doy, latitude) values
ncl_daylight = {}
doy_lat = tuple(map(tuple, daylight_data[::, 0:2]))
daylight_values = daylight_data[::, 2]
ncl_daylight = dict(zip(doy_lat, daylight_values))
### Collect DOY and Latitude Pairs
doy_lat = []
for doy in range(0, 365 + 1):
    for lat in range(-66, 66 + 1):
        doy_lat.append((doy, lat))

### Calculate GeoCAT-Comp `daylight_fao56` value and doy/lat
from geocat.comp import max_daylight

geocat_daylight = {}

for i, pair in enumerate(doy_lat):
    doy, lat = pair
    geocat_daylight[pair] = max_daylight(doy, lat)
/home/runner/micromamba/envs/geocat-applications/lib/python3.12/site-packages/geocat/comp/meteorology.py:1058: UserWarning: WARNING: max_daylight has limited validity for abs(lat) > 55 
  warnings.warn(

satvpr_temp_fao56#

#### Collect NCL values for satvpr_temp_fao56 from geocat-datafiles
import geocat.datafiles as gdf
import numpy as np

satvpr_temp_fao56_data = gdf.get(
    'applications_files/ncl_outputs/satvpr_temp_fao56_output.txt'
)
satvpr_temp_fao56_data = np.loadtxt(satvpr_temp_fao56_data, delimiter=',', skiprows=6)
Downloading file 'applications_files/ncl_outputs/satvpr_temp_fao56_output.txt' from 'https://github.com/NCAR/GeoCAT-datafiles/raw/main/applications_files/ncl_outputs/satvpr_temp_fao56_output.txt' to '/home/runner/.cache/geocat'.
### Collect NCL `satvpr_temp_fao56` value and associated (temp, satvpr_temp) values
ncl_satvpr_temp_fao56 = dict(
    zip(satvpr_temp_fao56_data[::, 0], satvpr_temp_fao56_data[::, 1])
)
### Calculate GeoCAT-Comp `saturation_vapor_pressure`
from geocat.comp import saturation_vapor_pressure

geocat_satvpr_temp_fao56 = {}

for temp in range(33, 212 + 1):
    geocat_satvpr_temp_fao56[temp] = saturation_vapor_pressure(temp)

satvpr_tdew_fao56#

#### Collect NCL values for satvpr_tdew_fao56 from geocat-datafiles
import geocat.datafiles as gdf
import numpy as np

satvpr_tdew_fao56_data = gdf.get(
    'applications_files/ncl_outputs/satvpr_tdew_fao56_output.txt'
)
satvpr_tdew_fao56_data = np.loadtxt(satvpr_tdew_fao56_data, delimiter=',', skiprows=6)
Downloading file 'applications_files/ncl_outputs/satvpr_tdew_fao56_output.txt' from 'https://github.com/NCAR/GeoCAT-datafiles/raw/main/applications_files/ncl_outputs/satvpr_tdew_fao56_output.txt' to '/home/runner/.cache/geocat'.
### Collect NCL `satvpr_tdew_fao56` value and associated (temp, act_sat_vapr_pressure) values
ncl_satvpr_tdew_fao56 = dict(
    zip(satvpr_tdew_fao56_data[::, 0], satvpr_tdew_fao56_data[::, 1])
)
### Calculate GeoCAT-Comp `actual_saturation_vapor_pressure`
from geocat.comp import actual_saturation_vapor_pressure

geocat_satvpr_tdew_fao56 = {}

for temp in range(33, 212 + 1):
    geocat_satvpr_tdew_fao56[temp] = actual_saturation_vapor_pressure(temp)

Comparison#

dewtemp_trh#

import math

for pair in ncl_dewtemp.keys():
    try:
        assert math.isclose(
            ncl_dewtemp[pair], geocat_dewtemp[pair], rel_tol=1e-04
        )  # within 4 decimal points
    except Exception:
        assert math.isclose(
            ncl_dewtemp[pair], geocat_dewtemp[pair], rel_tol=1e-02
        )  # within 2 decimal points
        print(f"{pair}:")
        print(f"\t{ncl_dewtemp[pair]}, {geocat_dewtemp[pair]}")
        print(f"\tDifference: {ncl_dewtemp[pair] - geocat_dewtemp[pair]}")
(274.0, 94.0):
	-0.0055542, -0.005570036309563875
	Difference: 1.5836309563875377e-05
(275.0, 87.0):
	-0.0839233, -0.08393959272956408
	Difference: 1.6292729564076902e-05
(275.0, 88.0):
	0.073761, 0.07374617196580857
	Difference: 1.482803419142198e-05
(277.0, 76.0):
	0.00259399, 0.0025903565796738803
	Difference: 3.633420326119678e-06
(278.0, 71.0):
	0.0261536, 0.026158530424595483
	Difference: -4.930424595483984e-06
(279.0, 66.0):
	-0.0278931, -0.027897415813697535
	Difference: 4.3158136975342265e-06
(281.0, 58.0):
	0.0712891, 0.07128017477202775
	Difference: 8.925227972245153e-06
(282.0, 54.0):
	0.0129395, 0.01291856405322278
	Difference: 2.0935946777218828e-05
(283.0, 50.0):
	-0.130432, -0.13044636060385528
	Difference: 1.4360603855290144e-05
(283.0, 51.0):
	0.144928, 0.14490829422311435
	Difference: 1.9705776885647897e-05
(284.0, 47.0):
	-0.0726318, -0.07262340593854333
	Difference: -8.39406145666799e-06
(285.0, 44.0):
	-0.0798035, -0.07982300528288988
	Difference: 1.950528288988118e-05
(286.0, 41.0):
	-0.160156, -0.1601734585760255
	Difference: 1.7458576025503048e-05
(287.0, 39.0):
	0.0386658, 0.0386467450997543
	Difference: 1.905490024570189e-05
(290.0, 32.0):
	-0.0795288, -0.07954541765832346
	Difference: 1.6617658323461737e-05
(291.0, 30.0):
	-0.115204, -0.11521555297190389
	Difference: 1.1552971903888709e-05
(294.0, 25.0):
	-0.110413, -0.11042499197856159
	Difference: 1.1991978561595729e-05
(295.0, 24.0):
	0.156677, 0.15665929836796977
	Difference: 1.7701632030242553e-05
(297.0, 21.0):
	-0.0614929, -0.06151243710161225
	Difference: 1.953710161224642e-05
(299.0, 19.0):
	0.162537, 0.16251872398265732
	Difference: 1.8276017342666595e-05
(302.0, 16.0):
	0.138275, 0.1382598702660971
	Difference: 1.5129733902913278e-05
(303.0, 15.0):
	0.0130615, 0.013067891707237322
	Difference: -6.391707237322214e-06

daylight_fao56#

import math

for pair in ncl_daylight.keys():
    assert math.isclose(
        ncl_daylight[pair], geocat_daylight[pair].flatten()[0], rel_tol=1e-05
    )  # within 5 decimal points

satvpr_temp_fao56#

import math

for key in ncl_satvpr_temp_fao56.keys():
    assert math.isclose(
        ncl_satvpr_temp_fao56[key], geocat_satvpr_temp_fao56[key], rel_tol=1e-05
    )  # within 5 decimal points

satvpr_tdew_fao56#

import math

for key in ncl_satvpr_tdew_fao56.keys():
    assert math.isclose(
        ncl_satvpr_tdew_fao56[key], geocat_satvpr_tdew_fao56[key], rel_tol=1e-05
    )  # within 5 decimal points