diff --git a/.github/workflows/AutoPushToPypi.yml b/.github/workflows/AutoPushToPypi.yml index fd7ea972..2e9a8ca4 100644 --- a/.github/workflows/AutoPushToPypi.yml +++ b/.github/workflows/AutoPushToPypi.yml @@ -1,16 +1,12 @@ -name: Push easyclimate to pypi +name: published easyclimate to pypi on: + release: + types: + - published push: - branches: [ main ] tags: - - v* - # - TEST* - # paths: - # - 'rapid_videocr/**' - # - 'docs/doc_whl.md' - # - 'setup.py' - # - '.github/workflows/AutoPushToPypi.yml' + - 'v*' jobs: UnitTesting: diff --git a/.github/workflows/TEST_AutoPushToPypi.yml b/.github/workflows/TEST_AutoPushToPypi.yml deleted file mode 100644 index a037112f..00000000 --- a/.github/workflows/TEST_AutoPushToPypi.yml +++ /dev/null @@ -1,52 +0,0 @@ -name: Push easyclimate to TestPyPI - -on: - push: - branches: [ dev ] - tags: - # - v* - - TEST* - # paths: - # - 'rapid_videocr/**' - # - 'docs/doc_whl.md' - # - 'setup.py' - # - '.github/workflows/AutoPushToPypi.yml' - -jobs: - UnitTesting: - runs-on: ubuntu-latest - steps: - - name: Pull latest code - uses: actions/checkout@v3 - - - name: Set up Python 3.10 - uses: actions/setup-python@v4 - with: - python-version: '3.10' - architecture: 'x64' - - - name: Unit testings - run: | - pip install -r test_requirements.txt - pytest --cov src - - GenerateWHL_PushPyPi: - needs: UnitTesting - runs-on: ubuntu-latest - - steps: - - uses: actions/checkout@v3 - - - name: Run setup.py - run: | - pip install --upgrade build - python -m pip install --upgrade pip - python -m build - - - name: Publish distribution 📦 to TestPyPI - uses: pypa/gh-action-pypi-publish@release/v1 - with: - # For test - password: ${{ secrets.TEST_PYPI_API_TOKEN }} - repository-url: https://test.pypi.org/legacy/ - packages-dir: dist/ \ No newline at end of file diff --git a/.github/workflows/ci.yaml b/.github/workflows/ci.yaml index 200cd06d..97666c3b 100644 --- a/.github/workflows/ci.yaml +++ b/.github/workflows/ci.yaml @@ -1,4 +1,4 @@ -name: Workflow for Codecov example-python +name: Workflow for Codecov easyclimate-python on: [push, pull_request] jobs: run: diff --git a/.gitignore b/.gitignore index c1065b80..189d8840 100644 --- a/.gitignore +++ b/.gitignore @@ -2,4 +2,5 @@ **/.vscode dist **/.pytest_cache -.coverage \ No newline at end of file +.coverage +.coverage* \ No newline at end of file diff --git a/docs/requirements.txt b/docs/requirements.txt index f74fa8f0..d280d605 100644 --- a/docs/requirements.txt +++ b/docs/requirements.txt @@ -4,6 +4,7 @@ geocat.viz <= 2023.10.0 cartopy >= 0.20 xeofs >= 1.1.0 matplotlib +pandas fast-barnes-py oceans intel-fortran-rt @@ -18,13 +19,15 @@ flox xarray-datatree xarray_regrid gsw_xarray +pooch +tqdm -sphinx == 7.2.6 -recommonmark -sphinx-markdown-tables -sphinxcontrib-jupyter -sphinx-inline-tabs -sphinx-gallery -sphinx-autoapi -sphinx-copybutton -sphinx-book-theme \ No newline at end of file +sphinx == 6.2.1 +recommonmark == 0.7.1 +sphinx-markdown-tables == 0.0.17 +sphinxcontrib-jupyter == 0.5.10 +sphinx-inline-tabs == 2023.4.21 +sphinx-gallery == 0.15.0 +sphinx-autoapi == 3.0.0 +sphinx-copybutton == 0.5.2 +sphinx-book-theme == 1.0.1 \ No newline at end of file diff --git a/release_requirements.txt b/release_requirements.txt index 6d674cb8..32a23bed 100644 --- a/release_requirements.txt +++ b/release_requirements.txt @@ -4,6 +4,7 @@ geocat.viz <= 2023.10.0 cartopy >= 0.20 xeofs >= 1.1.0 matplotlib +pandas fast-barnes-py oceans intel-fortran-rt @@ -17,4 +18,6 @@ pymannkendall flox xarray-datatree xarray_regrid -gsw_xarray \ No newline at end of file +gsw_xarray +pooch +tqdm \ No newline at end of file diff --git a/src/easyclimate/core/diff.py b/src/easyclimate/core/diff.py index c828b58d..b03ca00c 100644 --- a/src/easyclimate/core/diff.py +++ b/src/easyclimate/core/diff.py @@ -889,7 +889,7 @@ def calc_vertical_water_flux( specific_humidity_data: xr.DataArray, omega_data: xr.DataArray, g = 9.8 -) -> xr.Dataset: +) -> xr.DataArray: """ Calculate vertical water vapor flux. @@ -903,7 +903,11 @@ def calc_vertical_water_flux( omega_data: :py:class:`xarray.DataArray`. The vertical velocity data (:math:`\\frac{\\mathrm{d} p}{\\mathrm{d} t}`). g: :py:class:`float`, default: `9.8`. - The acceleration of gravity. + The acceleration of gravity. + + Returns + ------- + The vertical water flux. (:py:class:`xarray.DataArray `). """ water_flux = -omega_data *specific_humidity_data /g return water_flux diff --git a/src/easyclimate/core/tutorial.py b/src/easyclimate/core/tutorial.py index 80f62dbb..9cf3f524 100644 --- a/src/easyclimate/core/tutorial.py +++ b/src/easyclimate/core/tutorial.py @@ -43,6 +43,7 @@ def _construct_cache_dir(path): "shum_202201_mon_mean": 4, "uwnd_202201_mon_mean": 4, "vwnd_202201_mon_mean": 4, + "omega_202201_mon_mean": 4, "mini_HadISST_ice": 4, "PressQFF_202007271200_872": 'csv', } @@ -102,6 +103,7 @@ def open_tutorial_dataset( * ``"shum_202201_mon_mean"``: Absolute humidity of the NCEP reanalysis subset * ``"uwnd_202201_mon_mean"``: Zonal wind of the NCEP reanalysis subset * ``"vwnd_202201_mon_mean"``: Meridional wind of the NCEP reanalysis subset + * ``"omega_202201_mon_mean"``: Vertical velocity of the NCEP reanalysis subset * ``"mini_HadISST_ice"``: Hadley Centre Sea Ice and Sea Surface Temperature data set (HadISST) subset * ``"PressQFF_202007271200_872"``: Observational data from European stations (from https://github.com/EXCITED-CO2/xarray-regrid) diff --git a/src/easyclimate/core/variability.py b/src/easyclimate/core/variability.py index 21c1492b..41d1c806 100644 --- a/src/easyclimate/core/variability.py +++ b/src/easyclimate/core/variability.py @@ -22,7 +22,7 @@ def calc_climatological_mean(data_input: xr.DataArray, dim = 'time', **kwargs) - ------- :py:class:`xarray.DataArray`. """ - return data_input.mean(dim = dim) + return data_input.mean(dim = dim, **kwargs) def calc_climatological_seasonal_mean(data_input: xr.DataArray, dim = 'time', **kwargs) -> xr.DataArray: """ diff --git a/test/const_define.py b/test/const_define.py new file mode 100644 index 00000000..503d6100 --- /dev/null +++ b/test/const_define.py @@ -0,0 +1 @@ +TEST_DATA_PATH = 'test/data' \ No newline at end of file diff --git a/test/data/test_input_waveletpytest.nc b/test/data/test_input_waveletpytest.nc new file mode 100644 index 00000000..6d6f3cee Binary files /dev/null and b/test/data/test_input_waveletpytest.nc differ diff --git a/test/data/test_output_calc_absolute_vorticity.nc b/test/data/test_output_calc_absolute_vorticity.nc new file mode 100644 index 00000000..4c0bfff8 Binary files /dev/null and b/test/data/test_output_calc_absolute_vorticity.nc differ diff --git a/test/data/test_output_calc_gradient.nc b/test/data/test_output_calc_gradient.nc new file mode 100644 index 00000000..98ac750c Binary files /dev/null and b/test/data/test_output_calc_gradient.nc differ diff --git a/test/data/test_output_calc_helmholtz.nc b/test/data/test_output_calc_helmholtz.nc new file mode 100644 index 00000000..af91671c Binary files /dev/null and b/test/data/test_output_calc_helmholtz.nc differ diff --git a/test/data/test_output_calc_irrotational_component.nc b/test/data/test_output_calc_irrotational_component.nc new file mode 100644 index 00000000..e6ff99ed Binary files /dev/null and b/test/data/test_output_calc_irrotational_component.nc differ diff --git a/test/data/test_output_calc_nondivergent_component.nc b/test/data/test_output_calc_nondivergent_component.nc new file mode 100644 index 00000000..e4eccce8 Binary files /dev/null and b/test/data/test_output_calc_nondivergent_component.nc differ diff --git a/test/data/test_output_calc_planetary_vorticity.nc b/test/data/test_output_calc_planetary_vorticity.nc new file mode 100644 index 00000000..e6670b87 Binary files /dev/null and b/test/data/test_output_calc_planetary_vorticity.nc differ diff --git a/test/data/test_output_calc_relative_vorticity_and_horizontal_divergence.nc b/test/data/test_output_calc_relative_vorticity_and_horizontal_divergence.nc new file mode 100644 index 00000000..41401dbe Binary files /dev/null and b/test/data/test_output_calc_relative_vorticity_and_horizontal_divergence.nc differ diff --git a/test/data/test_output_calc_rossby_wave_source.nc b/test/data/test_output_calc_rossby_wave_source.nc new file mode 100644 index 00000000..29db015b Binary files /dev/null and b/test/data/test_output_calc_rossby_wave_source.nc differ diff --git a/test/data/test_output_calc_streamfunction_and_velocity_potential.nc b/test/data/test_output_calc_streamfunction_and_velocity_potential.nc new file mode 100644 index 00000000..99c9f5e2 Binary files /dev/null and b/test/data/test_output_calc_streamfunction_and_velocity_potential.nc differ diff --git a/test/data/test_output_calc_wind_speed.nc b/test/data/test_output_calc_wind_speed.nc new file mode 100644 index 00000000..066fe520 Binary files /dev/null and b/test/data/test_output_calc_wind_speed.nc differ diff --git a/test/data/test_output_waveletpytest.nc b/test/data/test_output_waveletpytest.nc new file mode 100644 index 00000000..28360373 Binary files /dev/null and b/test/data/test_output_waveletpytest.nc differ diff --git a/test/test_core_diff.py b/test/test_core_diff.py index 1f03e475..37c568cc 100644 --- a/test/test_core_diff.py +++ b/test/test_core_diff.py @@ -7,7 +7,6 @@ import numpy as np import xarray as xr import pandas as pd -from .util import round_sf_np # Sample Data Declaration t_data = xr.DataArray( @@ -92,6 +91,27 @@ coords = {'lat': np.array([5.0, 2.5, 0.0]), 'lon': np.array([0.0, 2.5, 5.0])} ) +omega_data_sample1 = xr.DataArray( + np.array([ 0.02437097, 0.03512904, 0.04098387, 0.02112902, 0.01658062, 0.02814513, 0.03577418, 0.04062096, 0.03781452, 0.02902419, 0.00769354, -0.01137904]), + dims = "level", + coords = {'level': np.array([1000., 925., 850., 700., 600., 500., 400., 300., 250., 200., 150., 100.])} +) + +q_data_sample1 = xr.DataArray( + np.array([17.74024, 13.399196, 9.72387, 5.151291, 3.0942986, 1.7616936 , 1.09875, 0.32160482]), + dims = "level", + coords = {'level': np.array([1000., 925., 850., 700., 600., 500., 400., 300.])} +) + +t_data_sample1 = xr.DataArray( + np.array([299.68384, 294.31772, 290.3024 , 283.88385, 277.28384, 269.40726, 259.0016 , 245.04675, 234.3129 , 221.59273, 206.91853, 191.03226, 194.09918, 203.64514, 211.79112, 221.36288, 229.66693]), + dims = "level", + coords = {'level': np.array([1000., 925., 850., 700., 600., 500., 400., 300., 250., 200., 150., 100., 70., 50., 30., 20., 10.])} +) + +# ----------------------------------------------------------- +# TEST + def test_calc_gradient(): result_data = ecl.calc_gradient(t_data, dim = 'lon').data refer_data = np.array([[ 0.12 , 0.08 , 0.04 ], [ 0.1 , 0. , -0.1 ], [ 0.025, -0.105, -0.235]]) @@ -222,6 +242,63 @@ def test_calc_horizontal_water_flux(): assert np.isclose(result_data1, refer_data1).all() assert np.isclose(result_data2, refer_data2).all() +def test_calc_vertical_water_flux(): + result_data = ecl.calc_vertical_water_flux(q_data_sample1, omega_data_sample1).data + refer_data = np.array([-0.04411703, -0.0480307 , -0.04066549, -0.0111063 , -0.00523524, -0.0050595 , -0.00401091, -0.00133305]) + assert np.isclose(result_data, refer_data).all() + +def test_calc_water_flux_top2surface_integral(): + specific_humidity_data_sample = ecl.open_tutorial_dataset('shum_202201_mon_mean').shum.isel(time = 0).sel(lat = slice(5, 0), lon = slice(0, 5)) + u_data_sample = ecl.open_tutorial_dataset('uwnd_202201_mon_mean').uwnd.isel(time = 0).sel(lat = slice(5, 0), lon = slice(0, 5)) + v_data_sample = ecl.open_tutorial_dataset('vwnd_202201_mon_mean').vwnd.isel(time = 0).sel(lat = slice(5, 0), lon = slice(0, 5)) + surface_pressure_data_sample = ecl.open_tutorial_dataset('pressfc_202201_mon_mean').pres.isel(time = 0).sel(lat = slice(5, 0), lon = slice(0, 5)) + + result_data = ecl.calc_water_flux_top2surface_integral( + specific_humidity_data_sample, + u_data_sample, + v_data_sample, + surface_pressure_data_sample, + surface_pressure_data_units = 'hPa', + vertical_dim = 'level', + vertical_dim_units = 'hPa', + ) + + result_data1 = result_data['qu'].data.flatten() + result_data2 = result_data['qv'].data.flatten() + refer_data1 = np.array([ -61719.65 , -53384.867, -46633.44 , -94774.73 , -84097.77 , -68822.93 , -127009.766, -119459.65 , -109226.625]) + refer_data2 = np.array([ 1478.0461, 4159.2134, -2113.73 , 19404.598 , 21102.201 ,16411.918 , 31748.395 , 31434.625 , 23815.54 ]) + + assert np.isclose(result_data1, refer_data1).all() + assert np.isclose(result_data2, refer_data2).all() + +def test_calc_divergence_watervaporflux(): + result_data = ecl.calc_divergence_watervaporflux( + q_data_500hpa, + u_data_500hpa, + v_data_500hpa, + specific_humidity_units = 'g/kg', + ).data.flatten() + refer_data = np.array([ 4.63115041e-10, 3.28054031e-10, 6.85267982e-11, -3.62682915e-10, -2.54289074e-10, -2.22359484e-10, -6.92384309e-10, -5.19564327e-10, -3.74825155e-10]) + assert np.isclose(result_data, refer_data).all() + +def test_calc_divergence_watervaporflux_top2surface_integral(): + specific_humidity_data_sample = ecl.open_tutorial_dataset('shum_202201_mon_mean').shum.isel(time = 0).sel(lat = slice(5, 0), lon = slice(0, 5)) + u_data_sample = ecl.open_tutorial_dataset('uwnd_202201_mon_mean').uwnd.isel(time = 0).sel(lat = slice(5, 0), lon = slice(0, 5)) + v_data_sample = ecl.open_tutorial_dataset('vwnd_202201_mon_mean').vwnd.isel(time = 0).sel(lat = slice(5, 0), lon = slice(0, 5)) + surface_pressure_data_sample = ecl.open_tutorial_dataset('pressfc_202201_mon_mean').pres.isel(time = 0).sel(lat = slice(5, 0), lon = slice(0, 5)) + + result_data = ecl.calc_divergence_watervaporflux_top2surface_integral( + specific_humidity_data_sample, + u_data_sample, + v_data_sample, + surface_pressure_data_sample, + vertical_dim = 'level', + specific_humidity_units = 'g/kg', + surface_pressure_data_units = 'hPa', + vertical_dim_units = 'hPa', + ).data.flatten() + refer_data = np.array([-3.16235016e-08, -5.25733773e-08, -6.11471233e-08, -1.84145484e-08, 3.67623974e-10, 1.46561626e-08, -1.93079642e-08, 1.67279636e-08, 3.28154620e-08]) + assert np.isclose(result_data, refer_data).all() def test_calc_u_advection(): result_data = ecl.calc_u_advection(u_data_500hpa, t_data).data.flatten() @@ -235,4 +312,16 @@ def test_calc_v_advection(): refer_data = np.array([ 1.33300385e-07, -9.17452042e-08, -9.71419810e-08, -5.18270458e-07, -4.12133849e-07, -6.69200313e-08, -3.96483196e-07, -1.24665542e-07, 2.81711745e-07]) + assert np.isclose(result_data, refer_data).all() + +def test_calc_p_advection(): + result_data = ecl.calc_p_advection( + omega_data_sample1, + t_data_sample1, + vertical_dim = 'level', + vertical_dim_units = 'hPa' + ) + refer_data = np.array([-1.96316937e-05, -2.19707321e-05, -1.90053499e-05, -1.10027766e-05, + -1.20015419e-05, -2.57278011e-05, -4.35738635e-05, -6.68585797e-05, + -8.86902508e-05, -7.95099400e-05, -2.35118198e-05, 1.82339871e-05]) assert np.isclose(result_data, refer_data).all() \ No newline at end of file diff --git a/test/test_core_stat.py b/test/test_core_stat.py new file mode 100644 index 00000000..347198d8 --- /dev/null +++ b/test/test_core_stat.py @@ -0,0 +1,163 @@ +""" +pytest for stat.py +""" +import pytest + +import easyclimate as ecl +import numpy as np +from .util import round_sf_np_new + +sst_data = ecl.tutorial.open_tutorial_dataset('mini_HadISST_sst').sst +sic_data_Barents_Sea = ecl.tutorial.open_tutorial_dataset('mini_HadISST_ice').sic +sic_data_Barents_Sea_12 = ecl.get_specific_months_data(sic_data_Barents_Sea, 12) + +def test_calc_linregress_spatial(): + result_data = ecl.calc_linregress_spatial(sic_data_Barents_Sea_12.sel(lon = slice(34.5, 36.5), lat = slice(78.5, 80.5)), dim = 'time').compute() + result_data1 = result_data['slope'].data + result_data2 = result_data['intercept'].data + result_data3 = result_data['rvalue'].data + result_data4 = result_data['pvalue'].data + result_data5 = result_data['stderr'].data + result_data6 = result_data['intercept_stderr'].data + + refer_data1 = np.array([[-0.01689814, -0.01618345, -0.01640629], + [-0.01011993, -0.00922373, -0.0091192 ], + [-0.00641115, -0.0054169 , -0.00600519]]) + refer_data2 = np.array([[1.05593576, 1.04461794, 1.04132889], + [1.01817285, 1.01218153, 1.01741975], + [0.89857147, 0.91509416, 0.94763013]]) + refer_data3 = np.array([[-0.58339207, -0.57217978, -0.57376992], + [-0.47457306, -0.4485609 , -0.45343254], + [-0.32601495, -0.28470031, -0.33127693]]) + refer_data4 = np.array([[5.01586941e-05, 7.52887062e-05, 7.11385575e-05], + [1.49647207e-03, 2.88846392e-03, 2.56361219e-03], + [3.51172733e-02, 6.76380713e-02, 3.21088821e-02]]) + refer_data5 = np.array([[0.00371969, 0.00366767, 0.00370284], + [0.00296779, 0.00290584, 0.00283422], + [0.00293946, 0.00288389, 0.00270435]]) + refer_data6 = np.array([[0.08858534, 0.08734656, 0.08818416], + [0.07067876, 0.06920341, 0.06749764], + [0.07000405, 0.06868049, 0.06440477]]) + + assert np.isclose(result_data1, refer_data1).all() + assert np.isclose(result_data2, refer_data2).all() + assert np.isclose(result_data3, refer_data3).all() + assert np.isclose(result_data4, refer_data4).all() + assert np.isclose(result_data5, refer_data5).all() + assert np.isclose(result_data6, refer_data6).all() + +def test_calc_linregress_spatial(): + result_data = ecl.calc_linregress_spatial(sic_data_Barents_Sea_12.sel(lon = slice(34.5, 36.5), lat = slice(78.5, 80.5)), dim = 'time').compute() + result_data1 = result_data['slope'].data + result_data2 = result_data['intercept'].data + result_data3 = result_data['rvalue'].data + result_data4 = result_data['pvalue'].data + result_data5 = result_data['stderr'].data + result_data6 = result_data['intercept_stderr'].data + + refer_data1 = np.array([[-0.01689814, -0.01618345, -0.01640629], + [-0.01011993, -0.00922373, -0.0091192 ], + [-0.00641115, -0.0054169 , -0.00600519]]) + refer_data2 = np.array([[1.05593576, 1.04461794, 1.04132889], + [1.01817285, 1.01218153, 1.01741975], + [0.89857147, 0.91509416, 0.94763013]]) + refer_data3 = np.array([[-0.58339207, -0.57217978, -0.57376992], + [-0.47457306, -0.4485609 , -0.45343254], + [-0.32601495, -0.28470031, -0.33127693]]) + refer_data4 = np.array([[5.01586941e-05, 7.52887062e-05, 7.11385575e-05], + [1.49647207e-03, 2.88846392e-03, 2.56361219e-03], + [3.51172733e-02, 6.76380713e-02, 3.21088821e-02]]) + refer_data5 = np.array([[0.00371969, 0.00366767, 0.00370284], + [0.00296779, 0.00290584, 0.00283422], + [0.00293946, 0.00288389, 0.00270435]]) + refer_data6 = np.array([[0.08858534, 0.08734656, 0.08818416], + [0.07067876, 0.06920341, 0.06749764], + [0.07000405, 0.06868049, 0.06440477]]) + + assert np.isclose(result_data1, refer_data1).all() + assert np.isclose(result_data2, refer_data2).all() + assert np.isclose(result_data3, refer_data3).all() + assert np.isclose(result_data4, refer_data4).all() + assert np.isclose(result_data5, refer_data5).all() + assert np.isclose(result_data6, refer_data6).all() + +def test_calc_linregress_spatial2(): + result_data = ecl.calc_linregress_spatial(sic_data_Barents_Sea_12.sel(lon = slice(34.5, 36.5), lat = slice(78.5, 80.5)), dim = 'time', engine = 'xarray').compute() + result_data1 = result_data['slope'].data + result_data2 = result_data['intercept'].data + result_data3 = result_data['rvalue'].data + result_data4 = result_data['pvalue'].data + result_data5 = result_data['stderr'].data + result_data6 = result_data['intercept_stderr'].data + + refer_data1 = np.array([[-0.01689814, -0.01618345, -0.01640629], + [-0.01011993, -0.00922373, -0.0091192 ], + [-0.00641115, -0.0054169 , -0.00600519]]) + refer_data2 = np.array([[1.0559357 , 1.04461782, 1.04132889], + [1.01817267, 1.01218159, 1.01741975], + [0.89857135, 0.91509404, 0.94763007]]) + refer_data3 = np.array([[-0.5833921 , -0.57217977, -0.57376993], + [-0.47457306, -0.4485609 , -0.45343254], + [-0.32601492, -0.2847003 , -0.33127695]]) + refer_data4 = np.array([[0.54670677, 0.55950058, 0.55297443], + [0.70697605, 0.73026712, 0.73443278], + [0.78751658, 0.82292746, 0.81043421]]) + refer_data5 = np.array([[0.02779874, 0.02749916, 0.02741885], + [0.02672881, 0.02656663, 0.02669478], + [0.02362677, 0.02404782, 0.02487059]]) + refer_data6 = np.array([[0.66203424, 0.6548997 , 0.65298707], + [0.63655362, 0.63269127, 0.63574306], + [0.5626778 , 0.57270522, 0.59229961]]) + + assert np.isclose(result_data1, refer_data1).all() + assert np.isclose(result_data2, refer_data2).all() + assert np.isclose(result_data3, refer_data3).all() + assert np.isclose(result_data4, refer_data4).all() + assert np.isclose(result_data5, refer_data5).all() + assert np.isclose(result_data6, refer_data6).all() + +def test_calc_detrend_data(): + result_data = ecl.calc_detrend_data(sic_data_Barents_Sea_12.sel(lon = slice(34.5, 36.5), lat = slice(78.5, 80.5)), time_dim = 'time').mean(dim = ('lat', 'lon')).data + result_data = round_sf_np_new(result_data) + refer_data = np.array([-5.9e-02, -4.9e-02, -4.5e-02, -4.9e-01, -2.2e-02, -4.7e-03, + 2.0e-02, 3.3e-02, 1.3e-02, 2.9e-02, -3.5e-04, 2.5e-02, + 5.3e-02, 3.0e-02, 5.6e-02, 1.5e-01, 5.8e-02, 1.5e-01, + 1.5e-01, 1.5e-01, 2.2e-02, 1.2e-01, 1.8e-01, 1.5e-01, + 1.2e-01, -3.3e-01, -9.8e-02, 1.6e-01, -1.2e-01, 2.9e-01, + 1.3e-01, -4.2e-01, 1.8e-01, 2.7e-01, -2.0e-01, -6.2e-01, + -2.3e-01, -6.0e-01, 3.9e-01, 5.8e-02, 3.9e-01, -7.8e-02]) + assert np.isclose(result_data, refer_data).all() + +def test_calc_ttestSpatialPattern_spatial(): + sic_data_Barents_Sea_12_spatial_mean = sic_data_Barents_Sea_12.sel(lon = slice(34.5, 36.5), lat = slice(78.5, 80.5)).mean(dim = ('lat', 'lon')) + test1 = sic_data_Barents_Sea_12_spatial_mean.isel(time = slice(0, 20)) + test2 = sic_data_Barents_Sea_12_spatial_mean.isel(time = slice(21, None)) + result_data = ecl.calc_ttestSpatialPattern_spatial(test1, test2) + result_data1 = result_data['statistic'].data + result_data2 = result_data['pvalue'].data + refer_data1 = 3.43382768 + refer_data2 = 0.00142438 + assert np.isclose(result_data1, refer_data1).all() + assert np.isclose(result_data2, refer_data2).all() + +def test_calc_skewness_spatial(): + sic_data_Barents_Sea_12_detrend = ecl.calc_detrend_data(sic_data_Barents_Sea_12.sel(lon = slice(34.5, 36.5), lat = slice(78.5, 80.5)), time_dim = 'time') + result_data = ecl.calc_skewness_spatial(sic_data_Barents_Sea_12_detrend, dim = 'time') + result_data1 = result_data['skewness'].data + result_data2 = result_data['pvalue'].data + refer_data1 = np.array([[-0.24333526, -0.32189173, -0.27430525], + [-1.3397094 , -1.5588326 , -1.6165946 ], + [-1.8677251 , -2.209491 , -2.330299 ]]) + refer_data2 = np.array([[7.70400089e-01, 6.26686378e-01, 7.18524740e-01], + [4.16622922e-04, 3.66448314e-05, 1.56112432e-05], + [1.88511919e-06, 6.86471937e-08, 1.30767304e-08]]) + assert np.isclose(result_data1, refer_data1).all() + assert np.isclose(result_data2, refer_data2).all() + +def test_calc_kurtosis_spatial(): + sic_data_Barents_Sea_12_detrend = ecl.calc_detrend_data(sic_data_Barents_Sea_12.sel(lon = slice(34.5, 36.5), lat = slice(78.5, 80.5)), time_dim = 'time') + result_data = ecl.calc_kurtosis_spatial(sic_data_Barents_Sea_12_detrend, dim = 'time') + refer_data = np.array([[2.7300231, 2.8368442, 2.7555804], + [4.667509 , 5.464381 , 5.85138 ], + [6.32586 , 7.463421 , 8.428185 ]]) + assert np.isclose(result_data, refer_data).all() diff --git a/test/test_core_variability.py b/test/test_core_variability.py new file mode 100644 index 00000000..dd732175 --- /dev/null +++ b/test/test_core_variability.py @@ -0,0 +1,224 @@ +""" +pytest for variability.py +""" +import pytest + +import easyclimate as ecl +import numpy as np +import xarray as xr +import pandas as pd + +sst_timeseries = xr.DataArray( + np.array([29.347528, 29.097271, 29.414158, 29.682327, 29.765087, 29.57674 , + 29.308327, 28.854443, 28.85187 , 28.844013, 28.935112, 28.792313, + 28.6857 , 28.806227, 28.847729, 29.080112, 29.265385, 29.086498, + 28.835426, 28.771742, 28.81437 , 29.067184, 29.301987, 28.995514, + 28.773684, 28.641384, 28.770582, 29.13947 , 29.389273, 29.19037 , + 29.083426, 29.017258, 29.187487, 29.556942, 29.831827, 29.698898, + 29.366127, 29.203444, 29.0869 , 29.385927, 29.388098, 29.417614, + 29.29614 , 29.305815, 29.38174 , 29.594328, 29.757042, 29.799171, + 29.330843, 29.243286, 29.418528, 29.68187 , 30.018873, 29.844812, + 29.370714, 29.427555, 29.571928, 29.725399, 29.908756, 29.902645, + 29.699043, 29.540981, 29.444414, 29.913399, 29.845999, 29.67817 , + 29.325956, 28.979944, 29.214754, 29.23727 , 29.879799, 29.828169, + 29.456327, 29.282442, 29.508482, 29.6117 , 29.852713, 29.888113, + 29.680742, 29.479841, 29.620369, 29.676613, 29.685156, 29.587656, + 29.343454, 29.086227, 29.113672, 29.399113, 29.718683, 29.696186, + 29.512985, 29.451885, 29.49677 , 29.541456, 29.792229, 29.554129, + 29.260298, 29.066557, 29.263298, 29.4407 , 29.649672, 29.6781 , + 29.418974, 29.404299, 29.438728, 29.844227, 30.02553 , 29.862041, + 29.43343 , 29.241543, 29.350756, 29.498257, 29.80587 , 29.7013 , + 29.594614, 29.44226 , 29.465628, 29.968756, 29.956043, 30.124586]), + dims = 'time', coords = {'time': pd.date_range('2010-01-01', '2019-12-31', freq='M')} +) + +u_data_sample = xr.DataArray( + np.array([-5.2400055, -1.949997 , -3.3600006, -5.220001 , -5.6100006, + -6.0099945, -4.25 , -2.75 , -3.9199982, -3.4199982, + -6.0099945, -4.199997 , -7.8600006, -2.800003 , -6.3600006, + -6.779999 , -5.3899994, -5.0099945, -4.419998 , -4.199997 , + -2.2700043, -3.600006 , -6.630005 , -7.5899963, -6.3600006, + -8.25 , -3.4100037, -7.220001 , -5.5700073, -4.300003 , + -3.0200043, -4.0400085, -2.8899994, -4.930008 , -4.149994 , + -7.2100067, -2.7599945, -5.9400024, -5.1600037, -2.369995 , + -6.5400085, -4.6900024, -4.1600037, -2.5200043, -3.369995 , + -4.9799957, -7.5899963, -6.169998 , -3.100006 , -2.4700012, + -6.180008 , -6.6600037, -5.5 , -3.680008 , -6.300003 , + -4.0700073, -2.3163245, -4.443382 , -6.664166 , -7.143548 , + -8.664516 , -5.034821 , -5.090321 , -6.2341657, -4.68629 , + -7.069165 , -3.5564506, -3.5879028, -3.556666 , -4.7411284, + -6.964999 , -6.316935 , -4.7403226, -9.636208 , -6.0008063, + -5.1908326, -5.1298375, -5.497499 , -3.5637095, -4.6508055, + -2.9158323, -2.7991927, -5.6016665, -5.955645 , -6.0774198, + -4.3589277, -5.5903215, -6.139166 , -5.458063 , -5.650831 , + -5.0185485, -5.5419345, -3.5658329, -6.1096764, -5.8241663, + -5.004838 , -6.2935476, -4.885714 , -5.116128 , -7.0016665, + -7.820967 , -6.6649985, -4.5129037, -2.7483861, -3.4741662, + -5.074192 , -5.454166 , -1.6782249, -5.7241926, -1.1714278, + -5.7370963, -6.636666 , -6.5782256, -4.6766653, -5.293547 , + -4.511289 , -3.5266657, -3.0112894, -5.5858326, -4.0266113]), + dims = 'time', coords = {'time': pd.date_range('2010-01-01', '2019-12-31', freq='M')} +) + +v_data_sample = xr.DataArray( + np.array([ 0.1499939 , -1.6399994 , -0.47999573, -0.47000122, -0.66999817, + -0.9900055 , -0.04000854, -0.33999634, -0.91000366, -0.2400055 , + 0.08000183, -0.33000183, 1.3199921 , 1.2299957 , -0.1000061 , + 0.30000305, -0.15000916, -0.32000732, -0.2400055 , -0.93000793, + -0.80999756, -0.9900055 , 0.11000061, 0.5399933 , 0.84999084, + 0.5399933 , 0.3899994 , 0.02000427, -0.58000183, -0.9900055 , + -0.30000305, -0.58000183, -0.66999817, -0.3600006 , -0.7200012 , + 1.1100006 , 0.7599945 , 0.00999451, -0.7700043 , -0.58999634, + -0.40000916, 0.09999084, -0.94000244, -0.44999695, -0.29000854, + -0.3500061 , 0.22000122, -0.07000732, 0.6300049 , 0.6300049 , + -0.68000793, -1.4700012 , -0.2400055 , -0.44999695, -1.1900024 , + -1.0599976 , -1.2403259 , -0.79515547, -0.81999993, -0.5935484 , + 1.7290322 , 0.16428573, 0.03225811, -1.0258332 , -0.7814516 , + -1.1516665 , -0.60483867, -0.83387095, -0.87333333, -0.6064515 , + 0.22666669, 1.8475807 , 0.6564516 , 2.9586208 , 0.707258 , + -1.0899999 , -0.87822574, -0.6491667 , -1.1580644 , -1.3169358 , + -0.18749999, -0.22419356, -0.85833323, 0.04838711, 0.6233871 , + 0.03392854, -0.6959676 , -0.41 , -0.8241936 , -0.37499997, + -0.7637095 , -0.6129031 , -0.8325001 , -0.50000006, -0.6058334 , + 1.8274194 , 0.72016126, 1.2785712 , -0.17177422, -1.4308336 , + -0.62661296, -1.2241668 , -1.2975805 , -0.4685484 , -0.4208334 , + -0.2669355 , 0.4108334 , 1.8129032 , 1.1870967 , 0.9285714 , + 0.9766128 , -0.8391666 , -0.08145163, -1.1983331 , -1.3419353 , + -1.8306453 , -0.4808334 , -0.04838709, 0.48499998, 0.7798387 ]), + dims = 'time', coords = {'time': pd.date_range('2010-01-01', '2019-12-31', freq='M')} +) + +uv_data_sample = xr.Dataset( + data_vars = {'u': u_data_sample, 'v': v_data_sample} +) + +def test_calc_climatological_mean(): + result_data = ecl.calc_climatological_mean(sst_timeseries).data + refer_data = np.array(29.41912939) + assert np.isclose(result_data, refer_data).all() + +def test_calc_climatological_seasonal_mean(): + result_data = ecl.calc_climatological_seasonal_mean(sst_timeseries).data + refer_data = np.array([29.3350306 , 29.37734163, 29.45836823, 29.5057771 ]) + assert np.isclose(result_data, refer_data).all() + +def test_calc_seasonal_cycle_mean(): + result_data = ecl.calc_seasonal_cycle_mean(sst_timeseries).data + refer_data = np.array([29.2696434, 29.1209362, 29.2218519, 29.4832875, 29.6699653, + 29.5757903, 29.3427304, 29.2135042, 29.3043644, 29.5056188, + 29.7073481, 29.6145122]) + assert np.isclose(result_data, refer_data).all() + +def test_calc_seasonal_cycle_std(): + result_data = ecl.calc_seasonal_cycle_std(sst_timeseries).data + refer_data = np.array([0.29286213, 0.23928619, 0.24478873, 0.24081124, 0.2320231 , + 0.25181168, 0.23270414, 0.26229246, 0.26926783, 0.33472974, + 0.31933075, 0.39431879]) + assert np.isclose(result_data, refer_data).all() + +def test_calc_seasonal_cycle_var(): + result_data = ecl.calc_seasonal_cycle_var(sst_timeseries).data + refer_data = np.array([0.08576823, 0.05725788, 0.05992152, 0.05799006, 0.05383472, + 0.06340912, 0.05415122, 0.06879733, 0.07250517, 0.112044, + 0.10197213, 0.1554873]) + assert np.isclose(result_data, refer_data).all() + +def test_remove_seasonal_cycle_mean(): + result_data = ecl.remove_seasonal_cycle_mean(sst_timeseries).data + refer_data = np.array([ 0.0778846, -0.0236652, 0.1923061, 0.1990395, 0.0951217, + 0.0009497, -0.0344034, -0.3590612, -0.4524944, -0.6616058, + -0.7722361, -0.8221992, -0.5839434, -0.3147092, -0.3741229, + -0.4031755, -0.4045803, -0.4892923, -0.5073044, -0.4417622, + -0.4899944, -0.4384348, -0.4053611, -0.6189982, -0.4959594, + -0.4795522, -0.4512699, -0.3438175, -0.2806923, -0.3854203, + -0.2593044, -0.1962462, -0.1168774, 0.0513232, 0.1244789, + 0.0843858, 0.0964836, 0.0825078, -0.1349519, -0.0973605, + -0.2818673, -0.1581763, -0.0465904, 0.0923108, 0.0773756, + 0.0887092, 0.0496939, 0.1846588, 0.0611996, 0.1223498, + 0.1966761, 0.1985825, 0.3489077, 0.2690217, 0.0279836, + 0.2140508, 0.2675636, 0.2197802, 0.2014079, 0.2881328, + 0.4293996, 0.4200448, 0.2225621, 0.4301115, 0.1760337, + 0.1023797, -0.0167744, -0.2335602, -0.0896104, -0.2683488, + 0.1724509, 0.2136568, 0.1866836, 0.1615058, 0.2866301, + 0.1284125, 0.1827477, 0.3123227, 0.3380116, 0.2663368, + 0.3160046, 0.1709942, -0.0221921, -0.0268562, 0.0738106, + -0.0347092, -0.1081799, -0.0841745, 0.0487177, 0.1203957, + 0.1702546, 0.2383808, 0.1924056, 0.0358372, 0.0848809, + -0.0603832, -0.0093454, -0.0543792, 0.0414461, -0.0425875, + -0.0202933, 0.1023097, 0.0762436, 0.1907948, 0.1343636, + 0.3386082, 0.3181819, 0.2475288, 0.1637866, 0.1206068, + 0.1289041, 0.0149695, 0.1359047, 0.1255097, 0.2518836, + 0.2287558, 0.1612636, 0.4631372, 0.2486949, 0.5100738]) + assert np.isclose(result_data, refer_data).all() + +def test_calc_climate_monthly_std(): + result_data = ecl.calc_climate_monthly_std(sst_timeseries).data + refer_data = np.array(0.2803478) + assert np.isclose(result_data, refer_data).all() + +def test_calc_climate_monthly_var(): + result_data = ecl.calc_climate_monthly_var(sst_timeseries).data + refer_data = np.array(0.07859489) + assert np.isclose(result_data, refer_data).all() + +def test_calc_horizontal_wind_components_std(): + result_data = ecl.calc_horizontal_wind_components_std(uv_data_sample) + result_data1 = result_data['u'].data + result_data2 = result_data['v'].data + result_data3 = result_data['sigma_s'].data + result_data4 = result_data['sigma_d'].data + + refer_data1 = np.array([-5.2400055, -1.949997 , -3.3600006, -5.220001 , -5.6100006, + -6.0099945, -4.25 , -2.75 , -3.9199982, -3.4199982, + -6.0099945, -4.199997 , -7.8600006, -2.800003 , -6.3600006, + -6.779999 , -5.3899994, -5.0099945, -4.419998 , -4.199997 , + -2.2700043, -3.600006 , -6.630005 , -7.5899963, -6.3600006, + -8.25 , -3.4100037, -7.220001 , -5.5700073, -4.300003 , + -3.0200043, -4.0400085, -2.8899994, -4.930008 , -4.149994 , + -7.2100067, -2.7599945, -5.9400024, -5.1600037, -2.369995 , + -6.5400085, -4.6900024, -4.1600037, -2.5200043, -3.369995 , + -4.9799957, -7.5899963, -6.169998 , -3.100006 , -2.4700012, + -6.180008 , -6.6600037, -5.5 , -3.680008 , -6.300003 , + -4.0700073, -2.3163245, -4.443382 , -6.664166 , -7.143548 , + -8.664516 , -5.034821 , -5.090321 , -6.2341657, -4.68629 , + -7.069165 , -3.5564506, -3.5879028, -3.556666 , -4.7411284, + -6.964999 , -6.316935 , -4.7403226, -9.636208 , -6.0008063, + -5.1908326, -5.1298375, -5.497499 , -3.5637095, -4.6508055, + -2.9158323, -2.7991927, -5.6016665, -5.955645 , -6.0774198, + -4.3589277, -5.5903215, -6.139166 , -5.458063 , -5.650831 , + -5.0185485, -5.5419345, -3.5658329, -6.1096764, -5.8241663, + -5.004838 , -6.2935476, -4.885714 , -5.116128 , -7.0016665, + -7.820967 , -6.6649985, -4.5129037, -2.7483861, -3.4741662, + -5.074192 , -5.454166 , -1.6782249, -5.7241926, -1.1714278, + -5.7370963, -6.636666 , -6.5782256, -4.6766653, -5.293547 , + -4.511289 , -3.5266657, -3.0112894, -5.5858326, -4.0266113]) + refer_data2 = np.array([ 0.1499939 , -1.6399994 , -0.47999573, -0.47000122, -0.66999817, + -0.9900055 , -0.04000854, -0.33999634, -0.91000366, -0.2400055 , + 0.08000183, -0.33000183, 1.3199921 , 1.2299957 , -0.1000061 , + 0.30000305, -0.15000916, -0.32000732, -0.2400055 , -0.93000793, + -0.80999756, -0.9900055 , 0.11000061, 0.5399933 , 0.84999084, + 0.5399933 , 0.3899994 , 0.02000427, -0.58000183, -0.9900055 , + -0.30000305, -0.58000183, -0.66999817, -0.3600006 , -0.7200012 , + 1.1100006 , 0.7599945 , 0.00999451, -0.7700043 , -0.58999634, + -0.40000916, 0.09999084, -0.94000244, -0.44999695, -0.29000854, + -0.3500061 , 0.22000122, -0.07000732, 0.6300049 , 0.6300049 , + -0.68000793, -1.4700012 , -0.2400055 , -0.44999695, -1.1900024 , + -1.0599976 , -1.2403259 , -0.79515547, -0.81999993, -0.5935484 , + 1.7290322 , 0.16428573, 0.03225811, -1.0258332 , -0.7814516 , + -1.1516665 , -0.60483867, -0.83387095, -0.87333333, -0.6064515 , + 0.22666669, 1.8475807 , 0.6564516 , 2.9586208 , 0.707258 , + -1.0899999 , -0.87822574, -0.6491667 , -1.1580644 , -1.3169358 , + -0.18749999, -0.22419356, -0.85833323, 0.04838711, 0.6233871 , + 0.03392854, -0.6959676 , -0.41 , -0.8241936 , -0.37499997, + -0.7637095 , -0.6129031 , -0.8325001 , -0.50000006, -0.6058334 , + 1.8274194 , 0.72016126, 1.2785712 , -0.17177422, -1.4308336 , + -0.62661296, -1.2241668 , -1.2975805 , -0.4685484 , -0.4208334 , + -0.2669355 , 0.4108334 , 1.8129032 , 1.1870967 , 0.9285714 , + 0.9766128 , -0.8391666 , -0.08145163, -1.1983331 , -1.3419353 , + -1.8306453 , -0.4808334 , -0.04838709, 0.48499998, 0.7798387 ]) + refer_data3 = np.array(1.58087643) + refer_data4 = np.array(0.16853612) + assert np.isclose(result_data1, refer_data1).all() + assert np.isclose(result_data2, refer_data2).all() + assert np.isclose(result_data3, refer_data3).all() + assert np.isclose(result_data4, refer_data4).all() \ No newline at end of file diff --git a/test/test_wavelet_wavelet_top.py b/test/test_wavelet_wavelet_top.py new file mode 100644 index 00000000..2313a086 --- /dev/null +++ b/test/test_wavelet_wavelet_top.py @@ -0,0 +1,36 @@ +""" +pytest for wavelet_top.py +""" +import pytest + +import easyclimate as ecl +import xarray as xr +import numpy as np +import os +from .const_define import TEST_DATA_PATH + +def test_timeseries_wavelet_transform(): + inputdata = xr.open_dataset(os.path.join(TEST_DATA_PATH, 'test_input_waveletpytest.nc')).sst + refer_data = xr.open_dataset(os.path.join(TEST_DATA_PATH, 'test_output_waveletpytest.nc')) + + result_data = ecl.wavelet.timeseries_wavelet_transform(inputdata, dt = 0.25) + result_data1 = result_data.power.data.flatten() + result_data2 = result_data.sig.data.flatten() + result_data3 = result_data.global_ws.data.flatten() + result_data4 = result_data.global_signif.data.flatten() + result_data5 = result_data.coi.data.flatten() + result_data6 = result_data.coi_bottom.data.flatten() + + refer_data1 = refer_data.power.data.flatten() + refer_data2 = refer_data.sig.data.flatten() + refer_data3 = refer_data.global_ws.data.flatten() + refer_data4 = refer_data.global_signif.data.flatten() + refer_data5 = refer_data.coi.data.flatten() + refer_data6 = refer_data.coi_bottom.data.flatten() + + assert np.isclose(result_data1, refer_data1).all() + assert np.isclose(result_data2, refer_data2).all() + assert np.isclose(result_data3, refer_data3).all() + assert np.isclose(result_data4, refer_data4).all() + assert np.isclose(result_data5, refer_data5).all() + assert np.isclose(result_data6, refer_data6).all() \ No newline at end of file diff --git a/test/test_windspharm.py b/test/test_windspharm.py new file mode 100644 index 00000000..ff3dd6a1 --- /dev/null +++ b/test/test_windspharm.py @@ -0,0 +1,184 @@ +""" +pytest for windspharm.top.py +""" +import pytest + +import easyclimate as ecl +import numpy as np +import xarray as xr +import pandas as pd +import os +from .const_define import TEST_DATA_PATH +from .util import round_sf_np_new # Intel fortran outputs for Windows and linux are quit different + +u_data_sample = ecl.open_tutorial_dataset('uwnd_202201_mon_mean').uwnd.isel(time = 0).sel(level = 500) +v_data_sample = ecl.open_tutorial_dataset('vwnd_202201_mon_mean').vwnd.isel(time = 0).sel(level = 500) + +lon_start, lon_end, lat_start, lat_end = 20, 30, -10, 10 + +def test_calc_wind_speed(): + result_data = ecl.windspharm.calc_wind_speed( + u_data = u_data_sample, + v_data = v_data_sample, + ) + result_data = result_data.sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + refer_data = xr.open_dataset(os.path.join(TEST_DATA_PATH, 'test_output_calc_wind_speed.nc'))['result'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + assert np.isclose(round_sf_np_new(result_data.flatten()), round_sf_np_new(refer_data.flatten())).all() + +def test_calc_relative_vorticity_and_horizontal_divergence(): + result_data = ecl.windspharm.calc_relative_vorticity_and_horizontal_divergence( + u_data = u_data_sample, + v_data = v_data_sample, + ) + result_data1 = result_data['vrt'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + result_data2 = result_data['div'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + refer_data = xr.open_dataset(os.path.join(TEST_DATA_PATH, 'test_output_calc_relative_vorticity_and_horizontal_divergence.nc')) + refer_data1 = refer_data['vrt'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + refer_data2 = refer_data['div'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + assert np.isclose(round_sf_np_new(result_data1.flatten()), round_sf_np_new(refer_data1.flatten())).all() + assert np.isclose(round_sf_np_new(result_data2.flatten()), round_sf_np_new(refer_data2.flatten())).all() + +def test_calc_relative_vorticity(): + result_data = ecl.windspharm.calc_relative_vorticity( + u_data = u_data_sample, + v_data = v_data_sample, + ) + result_data1 = result_data.sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + refer_data = xr.open_dataset(os.path.join(TEST_DATA_PATH, 'test_output_calc_relative_vorticity_and_horizontal_divergence.nc')) + refer_data1 = refer_data['vrt'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + assert np.isclose(round_sf_np_new(result_data1.flatten()), round_sf_np_new(refer_data1.flatten())).all() + +def test_calc_divergence(): + result_data = ecl.windspharm.calc_divergence( + u_data = u_data_sample, + v_data = v_data_sample, + ) + result_data1 = result_data.sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + refer_data = xr.open_dataset(os.path.join(TEST_DATA_PATH, 'test_output_calc_relative_vorticity_and_horizontal_divergence.nc')) + refer_data1 = refer_data['div'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + assert np.isclose(round_sf_np_new(result_data1.flatten()), round_sf_np_new(refer_data1.flatten())).all() + +def test_calc_planetary_vorticity(): + result_data = ecl.windspharm.calc_planetary_vorticity( + u_data = u_data_sample, + v_data = v_data_sample, + ) + result_data = result_data.sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + refer_data = xr.open_dataset(os.path.join(TEST_DATA_PATH, 'test_output_calc_planetary_vorticity.nc'))['result'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + assert np.isclose(round_sf_np_new(result_data.flatten()), round_sf_np_new(refer_data.flatten())).all() + +def test_calc_absolute_vorticity(): + result_data = ecl.windspharm.calc_absolute_vorticity( + u_data = u_data_sample, + v_data = v_data_sample, + ) + result_data = result_data.sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + refer_data = xr.open_dataset(os.path.join(TEST_DATA_PATH, 'test_output_calc_absolute_vorticity.nc'))['result'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + assert np.isclose(round_sf_np_new(result_data.flatten()), round_sf_np_new(refer_data.flatten())).all() + +def test_calc_streamfunction_and_velocity_potential(): + result_data = ecl.windspharm.calc_streamfunction_and_velocity_potential( + u_data = u_data_sample, + v_data = v_data_sample, + ) + result_data1 = result_data['stream'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + result_data2 = result_data['pv'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + refer_data = xr.open_dataset(os.path.join(TEST_DATA_PATH, 'test_output_calc_streamfunction_and_velocity_potential.nc')) + refer_data1 = refer_data['stream'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + refer_data2 = refer_data['pv'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + assert np.isclose(round_sf_np_new(result_data1.flatten()), round_sf_np_new(refer_data1.flatten())).all() + assert np.isclose(round_sf_np_new(result_data2.flatten()), round_sf_np_new(refer_data2.flatten())).all() + +def test_calc_streamfunction(): + result_data = ecl.windspharm.calc_streamfunction( + u_data = u_data_sample, + v_data = v_data_sample, + ) + result_data1 = result_data.sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + refer_data = xr.open_dataset(os.path.join(TEST_DATA_PATH, 'test_output_calc_streamfunction_and_velocity_potential.nc')) + refer_data1 = refer_data['stream'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + assert np.isclose(round_sf_np_new(result_data1.flatten()), round_sf_np_new(refer_data1.flatten())).all() + +def test_calc_velocity_potential(): + result_data = ecl.windspharm.calc_velocity_potential( + u_data = u_data_sample, + v_data = v_data_sample, + ) + result_data1 = result_data.sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + refer_data = xr.open_dataset(os.path.join(TEST_DATA_PATH, 'test_output_calc_streamfunction_and_velocity_potential.nc')) + refer_data1 = refer_data['pv'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + assert np.isclose(round_sf_np_new(result_data1.flatten()), round_sf_np_new(refer_data1.flatten())).all() + +def test_calc_helmholtz(): + result_data = ecl.windspharm.calc_helmholtz( + u_data = u_data_sample, + v_data = v_data_sample, + ) + result_data1 = result_data['uchi'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data.flatten() + result_data1 = round_sf_np_new(result_data1) + result_data2 = result_data['vchi'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data.flatten() + result_data2 = round_sf_np_new(result_data2) + result_data3 = result_data['upsi'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data.flatten() + result_data3 = round_sf_np_new(result_data3) + result_data4 = result_data['vpsi'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data.flatten() + result_data4 = round_sf_np_new(result_data4) + refer_data = xr.open_dataset(os.path.join(TEST_DATA_PATH, 'test_output_calc_helmholtz.nc')) + refer_data1 = refer_data['uchi'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data.flatten() + refer_data1 = round_sf_np_new(refer_data1) + refer_data2 = refer_data['vchi'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data.flatten() + refer_data2 = round_sf_np_new(refer_data2) + refer_data3 = refer_data['upsi'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data.flatten() + refer_data3 = round_sf_np_new(refer_data3) + refer_data4 = refer_data['vpsi'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data .flatten() + refer_data4 = round_sf_np_new(refer_data4) + assert np.isclose(result_data1, refer_data1).all() + assert np.isclose(result_data2, refer_data2).all() + assert np.isclose(result_data3, refer_data3).all() + assert np.isclose(result_data4, refer_data4).all() + +def test_calc_irrotational_component(): + result_data = ecl.windspharm.calc_irrotational_component( + u_data = u_data_sample, + v_data = v_data_sample, + ) + result_data1 = result_data['uchi'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data.flatten() + result_data2 = result_data['vchi'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data.flatten() + refer_data = xr.open_dataset(os.path.join(TEST_DATA_PATH, 'test_output_calc_irrotational_component.nc')) + refer_data1 = refer_data['uchi'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data.flatten() + refer_data2 = refer_data['vchi'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data.flatten() + assert np.isclose(round_sf_np_new(result_data1.flatten()), round_sf_np_new(refer_data1.flatten())).all() + assert np.isclose(round_sf_np_new(result_data2.flatten()), round_sf_np_new(refer_data2.flatten())).all() + +def test_calc_nondivergent_component(): + result_data = ecl.windspharm.calc_nondivergent_component( + u_data = u_data_sample, + v_data = v_data_sample, + ) + result_data1 = result_data['upsi'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data.flatten() + result_data2 = result_data['vpsi'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data.flatten() + refer_data = xr.open_dataset(os.path.join(TEST_DATA_PATH, 'test_output_calc_nondivergent_component.nc')) + refer_data1 = refer_data['upsi'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data.flatten() + refer_data2 = refer_data['vpsi'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data.flatten() + assert np.isclose(round_sf_np_new(result_data1.flatten()), round_sf_np_new(refer_data1.flatten())).all() + assert np.isclose(round_sf_np_new(result_data2.flatten()), round_sf_np_new(refer_data2.flatten())).all() + +def test_calc_rossby_wave_source(): + result_data = ecl.windspharm.calc_rossby_wave_source( + u_data = u_data_sample, + v_data = v_data_sample, + ) + result_data = result_data.sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + refer_data = xr.open_dataset(os.path.join(TEST_DATA_PATH, 'test_output_calc_rossby_wave_source.nc'))['result'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + assert np.isclose(round_sf_np_new(result_data.flatten()), round_sf_np_new(refer_data.flatten())).all() + +def test_calc_gradient(): + result_data = ecl.windspharm.calc_gradient( + data_input = u_data_sample, + ) + result_data1 = result_data['zonal_gradient'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + result_data2 = result_data['meridional_gradient'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + refer_data = xr.open_dataset(os.path.join(TEST_DATA_PATH, 'test_output_calc_gradient.nc')) + refer_data1 = refer_data['zonal_gradient'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + refer_data2 = refer_data['meridional_gradient'].sel(lon = slice(lon_start, lon_end), lat = slice(lat_end, lat_start)).data + assert np.isclose(round_sf_np_new(result_data1.flatten()), round_sf_np_new(refer_data1.flatten())).all() + assert np.isclose(round_sf_np_new(result_data2.flatten()), round_sf_np_new(refer_data2.flatten())).all() \ No newline at end of file diff --git a/test/util.py b/test/util.py index 79cb816a..b6b18600 100644 --- a/test/util.py +++ b/test/util.py @@ -9,4 +9,12 @@ def round_sf_np( """ r=np.ceil(np.log(x)/np.log(10)) f=significant_figure - return np.round(x*(10**(f-r)),0)*(10**(r-f)) \ No newline at end of file + return np.round(x*(10**(f-r)),0)*(10**(r-f)) + +def round_sf_np_new( + arr: np.array, +) -> np.array: + """ + Take Two significant figures (NEW) + """ + return np.array(["{:.2g}".format(num) for num in arr], dtype=float) \ No newline at end of file diff --git a/test/windspharm/__init__.py b/test/windspharm/__init__.py new file mode 100644 index 00000000..e888c2a5 --- /dev/null +++ b/test/windspharm/__init__.py @@ -0,0 +1,50 @@ +"""Tests for the `windspharm` package.""" +# Copyright (c) 2012-2016 Andrew Dawson +# +# Permission is hereby granted, free of charge, to any person obtaining a copy +# of this software and associated documentation files (the "Software"), to deal +# in the Software without restriction, including without limitation the rights +# to use, copy, modify, merge, publish, distribute, sublicense, and/or sell +# copies of the Software, and to permit persons to whom the Software is +# furnished to do so, subject to the following conditions: +# +# The above copyright notice and this permission notice shall be included in +# all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR +# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, +# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE +# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER +# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, +# OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN +# THE SOFTWARE. +from __future__ import absolute_import + +import pytest + +# import windspharm +from easyclimate import windspharm +from .utils import error + + +# Create a mapping from interface name to VectorWind class. +solvers = {'standard': windspharm.standard.VectorWind} +try: + solvers['cdms'] = windspharm.cdms.VectorWind +except AttributeError: + pass +try: + solvers['iris'] = windspharm.iris.VectorWind +except AttributeError: + pass +try: + solvers['xarray'] = windspharm.xarray.VectorWind +except AttributeError: + pass + + +class VectorWindTest(object): + """Base class for vector wind tests.""" + + def assert_error_is_zero(self, f1, f2): + assert error(f1, f2) == pytest.approx(0., abs=1e-5) diff --git a/test/windspharm/data/gaussian/chi.ref.npy b/test/windspharm/data/gaussian/chi.ref.npy new file mode 100644 index 00000000..2399ae46 Binary files /dev/null and b/test/windspharm/data/gaussian/chi.ref.npy differ diff --git a/test/windspharm/data/gaussian/chigradu.ref.npy b/test/windspharm/data/gaussian/chigradu.ref.npy new file mode 100644 index 00000000..8deb65e4 Binary files /dev/null and b/test/windspharm/data/gaussian/chigradu.ref.npy differ diff --git a/test/windspharm/data/gaussian/chigradv.ref.npy b/test/windspharm/data/gaussian/chigradv.ref.npy new file mode 100644 index 00000000..9df2973f Binary files /dev/null and b/test/windspharm/data/gaussian/chigradv.ref.npy differ diff --git a/test/windspharm/data/gaussian/div.ref.npy b/test/windspharm/data/gaussian/div.ref.npy new file mode 100644 index 00000000..78233193 Binary files /dev/null and b/test/windspharm/data/gaussian/div.ref.npy differ diff --git a/test/windspharm/data/gaussian/psi.ref.npy b/test/windspharm/data/gaussian/psi.ref.npy new file mode 100644 index 00000000..6fce1513 Binary files /dev/null and b/test/windspharm/data/gaussian/psi.ref.npy differ diff --git a/test/windspharm/data/gaussian/uchi.ref.npy b/test/windspharm/data/gaussian/uchi.ref.npy new file mode 100644 index 00000000..cfb3608b Binary files /dev/null and b/test/windspharm/data/gaussian/uchi.ref.npy differ diff --git a/test/windspharm/data/gaussian/upsi.ref.npy b/test/windspharm/data/gaussian/upsi.ref.npy new file mode 100644 index 00000000..6f385143 Binary files /dev/null and b/test/windspharm/data/gaussian/upsi.ref.npy differ diff --git a/test/windspharm/data/gaussian/uwnd.ref.npy b/test/windspharm/data/gaussian/uwnd.ref.npy new file mode 100644 index 00000000..903ae9f1 Binary files /dev/null and b/test/windspharm/data/gaussian/uwnd.ref.npy differ diff --git a/test/windspharm/data/gaussian/vchi.ref.npy b/test/windspharm/data/gaussian/vchi.ref.npy new file mode 100644 index 00000000..a5fcefa8 Binary files /dev/null and b/test/windspharm/data/gaussian/vchi.ref.npy differ diff --git a/test/windspharm/data/gaussian/vpsi.ref.npy b/test/windspharm/data/gaussian/vpsi.ref.npy new file mode 100644 index 00000000..652ee128 Binary files /dev/null and b/test/windspharm/data/gaussian/vpsi.ref.npy differ diff --git a/test/windspharm/data/gaussian/vrt.ref.npy b/test/windspharm/data/gaussian/vrt.ref.npy new file mode 100644 index 00000000..2ee56598 Binary files /dev/null and b/test/windspharm/data/gaussian/vrt.ref.npy differ diff --git a/test/windspharm/data/gaussian/vrt_trunc.ref.npy b/test/windspharm/data/gaussian/vrt_trunc.ref.npy new file mode 100644 index 00000000..cc8ccbe4 Binary files /dev/null and b/test/windspharm/data/gaussian/vrt_trunc.ref.npy differ diff --git a/test/windspharm/data/gaussian/vwnd.ref.npy b/test/windspharm/data/gaussian/vwnd.ref.npy new file mode 100644 index 00000000..fd9d2c35 Binary files /dev/null and b/test/windspharm/data/gaussian/vwnd.ref.npy differ diff --git a/test/windspharm/data/regular/chi.ref.npy b/test/windspharm/data/regular/chi.ref.npy new file mode 100644 index 00000000..1450f8ce Binary files /dev/null and b/test/windspharm/data/regular/chi.ref.npy differ diff --git a/test/windspharm/data/regular/chigradu.ref.npy b/test/windspharm/data/regular/chigradu.ref.npy new file mode 100644 index 00000000..6a3f1c56 Binary files /dev/null and b/test/windspharm/data/regular/chigradu.ref.npy differ diff --git a/test/windspharm/data/regular/chigradv.ref.npy b/test/windspharm/data/regular/chigradv.ref.npy new file mode 100644 index 00000000..54ce661a Binary files /dev/null and b/test/windspharm/data/regular/chigradv.ref.npy differ diff --git a/test/windspharm/data/regular/div.ref.npy b/test/windspharm/data/regular/div.ref.npy new file mode 100644 index 00000000..218c6c4c Binary files /dev/null and b/test/windspharm/data/regular/div.ref.npy differ diff --git a/test/windspharm/data/regular/psi.ref.npy b/test/windspharm/data/regular/psi.ref.npy new file mode 100644 index 00000000..9a5aadf7 Binary files /dev/null and b/test/windspharm/data/regular/psi.ref.npy differ diff --git a/test/windspharm/data/regular/uchi.ref.npy b/test/windspharm/data/regular/uchi.ref.npy new file mode 100644 index 00000000..8ad04278 Binary files /dev/null and b/test/windspharm/data/regular/uchi.ref.npy differ diff --git a/test/windspharm/data/regular/upsi.ref.npy b/test/windspharm/data/regular/upsi.ref.npy new file mode 100644 index 00000000..2ca12223 Binary files /dev/null and b/test/windspharm/data/regular/upsi.ref.npy differ diff --git a/test/windspharm/data/regular/uwnd.ref.npy b/test/windspharm/data/regular/uwnd.ref.npy new file mode 100644 index 00000000..589c2794 Binary files /dev/null and b/test/windspharm/data/regular/uwnd.ref.npy differ diff --git a/test/windspharm/data/regular/vchi.ref.npy b/test/windspharm/data/regular/vchi.ref.npy new file mode 100644 index 00000000..358bd803 Binary files /dev/null and b/test/windspharm/data/regular/vchi.ref.npy differ diff --git a/test/windspharm/data/regular/vpsi.ref.npy b/test/windspharm/data/regular/vpsi.ref.npy new file mode 100644 index 00000000..963d1dbd Binary files /dev/null and b/test/windspharm/data/regular/vpsi.ref.npy differ diff --git a/test/windspharm/data/regular/vrt.ref.npy b/test/windspharm/data/regular/vrt.ref.npy new file mode 100644 index 00000000..b4ae5635 Binary files /dev/null and b/test/windspharm/data/regular/vrt.ref.npy differ diff --git a/test/windspharm/data/regular/vrt_trunc.ref.npy b/test/windspharm/data/regular/vrt_trunc.ref.npy new file mode 100644 index 00000000..bf8f0447 Binary files /dev/null and b/test/windspharm/data/regular/vrt_trunc.ref.npy differ diff --git a/test/windspharm/data/regular/vwnd.ref.npy b/test/windspharm/data/regular/vwnd.ref.npy new file mode 100644 index 00000000..60839e04 Binary files /dev/null and b/test/windspharm/data/regular/vwnd.ref.npy differ diff --git a/test/windspharm/reference.py b/test/windspharm/reference.py new file mode 100644 index 00000000..31da24b0 --- /dev/null +++ b/test/windspharm/reference.py @@ -0,0 +1,134 @@ +"""Reference solutions for testing the `windspharm` package.""" +# Copyright (c) 2012-2016 Andrew Dawson +# +# Permission is hereby granted, free of charge, to any person obtaining a copy +# of this software and associated documentation files (the "Software"), to deal +# in the Software without restriction, including without limitation the rights +# to use, copy, modify, merge, publish, distribute, sublicense, and/or sell +# copies of the Software, and to permit persons to whom the Software is +# furnished to do so, subject to the following conditions: +# +# The above copyright notice and this permission notice shall be included in +# all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR +# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, +# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE +# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER +# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, +# OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN +# THE SOFTWARE. +from __future__ import absolute_import +import os + +import numpy as np +from spharm import gaussian_lats_wts + + +def test_data_path(): + return os.path.join(os.path.abspath(os.path.dirname(__file__)), 'data') + + +def __read_reference_solutions(gridtype): + """Read reference solutions from file.""" + exact = dict() + for varid in ('psi', 'chi', 'vrt', 'div', 'uchi', 'vchi', 'upsi', 'vpsi', + 'chigradu', 'chigradv', 'uwnd', 'vwnd', 'vrt_trunc'): + try: + filename = os.path.join(test_data_path(), gridtype, + '{!s}.ref.npy'.format(varid)) + exact[varid] = np.load(filename).squeeze() + except IOError: + msg = 'required data file not found: {!s}' + raise IOError(msg.format(filename)) + return exact + + +def _wrap_cdms(reference, lats, lons): + try: + import cdms2 + except ImportError: + raise ValueError("cannot use container 'cdms' without cdms2") + londim = cdms2.createAxis(lons, id='longitude') + londim.designateLongitude() + latdim = cdms2.createAxis(lats, id='latitude') + latdim.designateLatitude() + for name in reference.keys(): + reference[name] = cdms2.createVariable(reference[name], + axes=[latdim, londim], + id=name) + + +def _wrap_iris(reference, lats, lons): + try: + from iris.cube import Cube + from iris.coords import DimCoord + except ImportError: + raise ValueError("cannot use container 'iris' without iris") + + londim = DimCoord(lons, + standard_name='longitude', + units='degrees_east') + latdim = DimCoord(lats, + standard_name='latitude', + units='degrees_north') + coords = list(zip((latdim, londim), (0, 1))) + for name in reference.keys(): + reference[name] = Cube(reference[name], + dim_coords_and_dims=coords, + long_name=name) + + +def _wrap_xarray(reference, lats, lons): + try: + import xarray as xr + except ImportError: + try: + import xray as xr + except ImportError: + raise ValueError("cannot use container 'xarray' without xarray") + londim = xr.IndexVariable('longitude', lons, + attrs={'standard_name': 'longitude', + 'units': 'degrees_east'}) + latdim = xr.IndexVariable('latitude', lats, + attrs={'standard_name': 'latitude', + 'units': 'degrees_north'}) + for name in reference.keys(): + reference[name] = xr.DataArray(reference[name], + coords=[latdim, londim], + attrs={'long_name': name}) + + +def _get_wrapper(container_type): + if container_type == 'cdms': + return _wrap_cdms + elif container_type == 'iris': + return _wrap_iris + elif container_type == 'xarray': + return _wrap_xarray + else: + raise ValueError('invalid container type: {!s}'.format(container_type)) + + +def reference_solutions(container_type, gridtype): + """Generate reference solutions in the required container.""" + container_type = container_type.lower() + if container_type not in ('standard', 'iris', 'cdms', 'xarray'): + raise ValueError("unknown container type: " + "'{!s}'".format(container_type)) + reference = __read_reference_solutions(gridtype) + if container_type == 'standard': + # Reference solution already in numpy arrays. + return reference + # Generate coordinate dimensions for meta-data interfaces. + if gridtype == 'gaussian': + lats, _ = gaussian_lats_wts(72) + else: + lats = np.linspace(90, -90, 73) + lons = np.arange(0, 360, 2.5) + _get_wrapper(container_type)(reference, lats, lons) + return reference + + +if __name__ == '__main__': + pass diff --git a/test/windspharm/test_description.txt b/test/windspharm/test_description.txt new file mode 100644 index 00000000..80dba7bd --- /dev/null +++ b/test/windspharm/test_description.txt @@ -0,0 +1,6 @@ +The computational tests compare the solution for the various quantities +to a reference solution obtained using the NCL Spherepack interface, which +is very well tested and used. + +The reference solutions are stored in NumPy binary files in the subdirectory +'data'. There is one file per reference field. diff --git a/test/windspharm/test_error_handling.py b/test/windspharm/test_error_handling.py new file mode 100644 index 00000000..be1e8b83 --- /dev/null +++ b/test/windspharm/test_error_handling.py @@ -0,0 +1,398 @@ +"""Tests for error handling in `windspharm`.""" +# Copyright (c) 2012-2016 Andrew Dawson +# +# Permission is hereby granted, free of charge, to any person obtaining a copy +# of this software and associated documentation files (the "Software"), to deal +# in the Software without restriction, including without limitation the rights +# to use, copy, modify, merge, publish, distribute, sublicense, and/or sell +# copies of the Software, and to permit persons to whom the Software is +# furnished to do so, subject to the following conditions: +# +# The above copyright notice and this permission notice shall be included in +# all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR +# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, +# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE +# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER +# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, +# OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN +# THE SOFTWARE. +from __future__ import absolute_import + +import pytest +import numpy as np +import numpy.ma as ma + +# from windspharm.tests import VectorWindTest, solvers +from . import VectorWindTest, solvers +from .reference import reference_solutions + + +class ErrorHandlersTest(VectorWindTest): + """Base class for all error handler tests.""" + interface = None + gridtype = None + + @classmethod + def setup_class(cls): + msg = 'missing dependencies required to test the {!s} interface' + if cls.interface not in solvers: + pytest.skip(msg.format(cls.interface)) + + +# ---------------------------------------------------------------------------- +# Tests for the standard interface + + +class TestStandardErrorHandlers(ErrorHandlersTest): + """Standard interface error handler tests.""" + interface = 'standard' + gridtype = 'regular' + + def test_masked_values(self): + # masked values in inputs should raise an error + solution = reference_solutions(self.interface, self.gridtype) + mask = np.empty(solution['uwnd'].shape, dtype=np.bool_) + mask[:] = False + mask[1, 1] = True + u = ma.array(solution['uwnd'], mask=mask, fill_value=1.e20) + v = ma.array(solution['vwnd'], mask=mask, fill_value=1.e20) + with pytest.raises(ValueError): + solvers[self.interface](u, v, gridtype=self.gridtype) + + def test_nan_values(self): + # NaN values in inputs should raise an error + solution = reference_solutions(self.interface, self.gridtype) + solution['vwnd'][1, 1] = np.nan + with pytest.raises(ValueError): + solvers[self.interface](solution['uwnd'], + solution['vwnd'], + gridtype=self.gridtype) + + def test_invalid_shape_components(self): + # invalid shape inputs should raise an error + solution = reference_solutions(self.interface, self.gridtype) + with pytest.raises(ValueError): + solvers[self.interface]( + solution['uwnd'][np.newaxis].repeat(2, axis=0), + solution['vwnd'][np.newaxis].repeat(2, axis=0), + gridtype=self.gridtype) + + def test_different_shape_components(self): + # different shape inputs should raise an error + solution = reference_solutions(self.interface, self.gridtype) + with pytest.raises(ValueError): + solvers[self.interface](solution['uwnd'], + solution['vwnd'][:-1], + gridtype=self.gridtype) + + def test_invalid_rank_components(self): + # invalid rank inputs should raise an error + solution = reference_solutions(self.interface, self.gridtype) + with pytest.raises(ValueError): + solvers[self.interface]( + solution['uwnd'][..., np.newaxis, np.newaxis], + solution['vwnd'][..., np.newaxis, np.newaxis], + gridtype=self.gridtype) + + def test_different_rank_components(self): + # different rank inputs should raise an error + solution = reference_solutions(self.interface, self.gridtype) + with pytest.raises(ValueError): + solvers[self.interface](solution['uwnd'][..., np.newaxis], + solution['vwnd'], + gridtype=self.gridtype) + + def test_invalid_gridtype(self): + # invalid grid type specification should raise an error + solution = reference_solutions(self.interface, self.gridtype) + with pytest.raises(ValueError): + solvers[self.interface](solution['uwnd'], solution['vwnd'], + gridtype='curvilinear') + + def test_gradient_masked_values(self): + # masked values in gradient input should raise an error + solution = reference_solutions(self.interface, self.gridtype) + vw = solvers[self.interface](solution['uwnd'], solution['vwnd'], + gridtype=self.gridtype) + mask = np.empty(solution['uwnd'].shape, dtype=np.bool_) + mask[:] = False + mask[1, 1] = True + chi = ma.array(solution['chi'], mask=mask, fill_value=1.e20) + with pytest.raises(ValueError): + vw.gradient(chi) + + def test_gradient_nan_values(self): + # NaN values in gradient input should raise an error + solution = reference_solutions(self.interface, self.gridtype) + vw = solvers[self.interface](solution['uwnd'], solution['vwnd'], + gridtype=self.gridtype) + solution['chi'][1, 1] = np.nan + with pytest.raises(ValueError): + vw.gradient(solution['chi']) + + def test_gradient_invalid_shape(self): + # input to gradient of different shape should raise an error + solution = reference_solutions(self.interface, self.gridtype) + vw = solvers[self.interface](solution['uwnd'], solution['vwnd'], + gridtype=self.gridtype) + with pytest.raises(ValueError): + vw.gradient(solution['chi'][:-1]) + + +# ---------------------------------------------------------------------------- +# Tests for the cdms interface + + +class TestCDMSErrorHandlers(ErrorHandlersTest): + """cdms interface error handler tests.""" + interface = 'cdms' + gridtype = 'regular' + + def test_non_variable_input(self): + # input not a cdms2 variable should raise an error + solution = reference_solutions('standard', self.gridtype) + with pytest.raises(TypeError): + solvers[self.interface](solution['uwnd'], solution['vwnd']) + + def test_different_shape_components(self): + # inputs not the same shape should raise an error + solution = reference_solutions(self.interface, self.gridtype) + with pytest.raises(ValueError): + solvers[self.interface](solution['uwnd'], + solution['vwnd'].reorder('xy')) + + def test_unknown_grid(self): + # inputs where a lat-lon grid cannot be identified should raise an + # error + solution = reference_solutions(self.interface, self.gridtype) + lat = solution['vwnd'].getLatitude() + delattr(lat, 'axis') + lat.id = 'unknown' + with pytest.raises(ValueError): + solvers[self.interface](solution['uwnd'], solution['vwnd']) + + def test_non_variable_gradient_input(self): + # input to gradient not a cdms2 variable should raise an error + solution = reference_solutions(self.interface, self.gridtype) + vw = solvers[self.interface](solution['uwnd'], solution['vwnd']) + dummy_solution = reference_solutions('standard', self.gridtype) + with pytest.raises(TypeError): + vw.gradient(dummy_solution['chi']) + + def test_gradient_non_variable_input(self): + # input to gradient not a cdms2 variable should raise an error + solution = reference_solutions(self.interface, self.gridtype) + vw = solvers[self.interface](solution['uwnd'], solution['vwnd']) + dummy_solution = reference_solutions('standard', self.gridtype) + with pytest.raises(TypeError): + vw.gradient(dummy_solution['chi']) + + def test_gradient_different_shape(self): + # input to gradient of different shape should raise an error + solution = reference_solutions(self.interface, self.gridtype) + vw = solvers[self.interface](solution['uwnd'], solution['vwnd']) + with pytest.raises(ValueError): + vw.gradient(solution['chi'][:-1]) + + def test_gradient_unknown_grid(self): + # input to gradient with no identifiable grid should raise an error + solution = reference_solutions(self.interface, self.gridtype) + vw = solvers[self.interface](solution['uwnd'], solution['vwnd']) + lat = solution['chi'].getLatitude() + delattr(lat, 'axis') + lat.id = 'unknown' + with pytest.raises(ValueError): + vw.gradient(solution['chi']) + + def test_truncate_non_variable_input(self): + # input to truncate not a cdms2 variable should raise an error + solution = reference_solutions(self.interface, self.gridtype) + vw = solvers[self.interface](solution['uwnd'], solution['vwnd']) + dummy_solution = reference_solutions('standard', self.gridtype) + with pytest.raises(TypeError): + vw.truncate(dummy_solution['chi']) + + def test_truncate_different_shape(self): + # input to truncate of different shape should raise an error + solution = reference_solutions(self.interface, self.gridtype) + vw = solvers[self.interface](solution['uwnd'], solution['vwnd']) + with pytest.raises(ValueError): + vw.truncate(solution['chi'][:-1]) + + def test_truncate_unknown_grid(self): + # input to truncate with no identifiable grid should raise an error + solution = reference_solutions(self.interface, self.gridtype) + vw = solvers[self.interface](solution['uwnd'], solution['vwnd']) + lat = solution['chi'].getLatitude() + delattr(lat, 'axis') + lat.id = 'unknown' + with pytest.raises(ValueError): + vw.truncate(solution['chi']) + + +# ---------------------------------------------------------------------------- +# Tests for the iris interface + + +class TestIrisErrorHandlers(ErrorHandlersTest): + """Iris interface error handler tests.""" + interface = 'iris' + gridtype = 'regular' + + def test_non_cube_input(self): + # input not an iris cube should raise an error + solution = reference_solutions('standard', self.gridtype) + with pytest.raises(TypeError): + solvers[self.interface](solution['uwnd'], solution['vwnd']) + + def test_different_shape_components(self): + # inputs not the same shape should raise an error + solution = reference_solutions(self.interface, self.gridtype) + solution['vwnd'].transpose([1, 0]) + with pytest.raises(ValueError): + solvers[self.interface](solution['uwnd'], solution['vwnd']) + + def test_unknown_grid(self): + # inputs where a lat-lon grid cannot be identified should raise an + # error + solution = reference_solutions(self.interface, self.gridtype) + solution['vwnd'].coord('latitude').rename('unknown') + with pytest.raises(ValueError): + solvers[self.interface](solution['uwnd'], solution['vwnd']) + + def test_gradient_non_cube_input(self): + # input to gradient not an iris cube should raise an error + solution = reference_solutions(self.interface, self.gridtype) + vw = solvers[self.interface](solution['uwnd'], solution['vwnd']) + dummy_solution = reference_solutions('standard', self.gridtype) + with pytest.raises(TypeError): + vw.gradient(dummy_solution['chi']) + + def test_gradient_different_shape(self): + # input to gradient of different shape should raise an error + solution = reference_solutions(self.interface, self.gridtype) + vw = solvers[self.interface](solution['uwnd'], solution['vwnd']) + with pytest.raises(ValueError): + vw.gradient(solution['chi'][:-1]) + + def test_gradient_unknown_grid(self): + # input to gradient with no identifiable grid should raise an error + solution = reference_solutions(self.interface, self.gridtype) + vw = solvers[self.interface](solution['uwnd'], solution['vwnd']) + solution['chi'].coord('latitude').rename('unknown') + with pytest.raises(ValueError): + vw.gradient(solution['chi']) + + def test_truncate_non_cube_input(self): + # input to truncate not an iris cube should raise an error + solution = reference_solutions(self.interface, self.gridtype) + vw = solvers[self.interface](solution['uwnd'], solution['vwnd']) + dummy_solution = reference_solutions('standard', self.gridtype) + with pytest.raises(TypeError): + vw.truncate(dummy_solution['chi']) + + def test_truncate_different_shape(self): + # input to truncate of different shape should raise an error + solution = reference_solutions(self.interface, self.gridtype) + vw = solvers[self.interface](solution['uwnd'], solution['vwnd']) + with pytest.raises(ValueError): + vw.truncate(solution['chi'][:-1]) + + def test_truncate_unknown_grid(self): + # input to truncate with no identifiable grid should raise an error + solution = reference_solutions(self.interface, self.gridtype) + vw = solvers[self.interface](solution['uwnd'], solution['vwnd']) + solution['chi'].coord('latitude').rename('unknown') + with pytest.raises(ValueError): + vw.truncate(solution['chi']) + + +# ---------------------------------------------------------------------------- +# Tests for the xarray interface + + +class TestXarrayErrorHandlers(ErrorHandlersTest): + """xarray interface error handler tests.""" + interface = 'xarray' + gridtype = 'regular' + + def test_non_dataarray_input(self): + # input not an xarray.DataArray should raise an error + solution = reference_solutions('standard', self.gridtype) + with pytest.raises(TypeError): + solvers[self.interface](solution['uwnd'], solution['vwnd']) + + def test_different_shape_components(self): + # inputs not the same shape should raise an error + solution = reference_solutions(self.interface, self.gridtype) + solution['vwnd'] = solution['vwnd'].transpose() + with pytest.raises(ValueError): + solvers[self.interface](solution['uwnd'], solution['vwnd']) + + def test_unknown_grid(self): + # inputs where a lat-lon grid cannot be identified should raise an + # error + solution = reference_solutions(self.interface, self.gridtype) + solution['vwnd'].coords.update( + {'unknown': ('latitude', + solution['vwnd'].coords['latitude'].values)}) + solution['vwnd'] = solution['vwnd'].swap_dims({'latitude': 'unknown'}) + del solution['vwnd'].coords['latitude'] + with pytest.raises(ValueError): + solvers[self.interface](solution['uwnd'], solution['vwnd']) + + def test_gradient_non_dataarray_input(self): + # input to gradient not an xarray.DataArray should raise an error + solution = reference_solutions(self.interface, self.gridtype) + vw = solvers[self.interface](solution['uwnd'], solution['vwnd']) + dummy_solution = reference_solutions('standard', self.gridtype) + with pytest.raises(TypeError): + vw.gradient(dummy_solution['chi']) + + def test_gradient_different_shape(self): + # input to gradient of different shape should raise an error + solution = reference_solutions(self.interface, self.gridtype) + vw = solvers[self.interface](solution['uwnd'], solution['vwnd']) + with pytest.raises(ValueError): + vw.gradient(solution['chi'][:-1]) + + def test_gradient_unknown_grid(self): + # input to gradient with no identifiable grid should raise an error + solution = reference_solutions(self.interface, self.gridtype) + vw = solvers[self.interface](solution['uwnd'], solution['vwnd']) + solution['chi'].coords.update( + {'unknown': ('latitude', + solution['chi'].coords['latitude'].values)}) + solution['chi'] = solution['chi'].swap_dims({'latitude': 'unknown'}) + del solution['chi'].coords['latitude'] + with pytest.raises(ValueError): + vw.gradient(solution['chi']) + + def test_truncate_non_dataarray_input(self): + # input to truncate not an xarray.DataArray should raise an error + solution = reference_solutions(self.interface, self.gridtype) + vw = solvers[self.interface](solution['uwnd'], solution['vwnd']) + dummy_solution = reference_solutions('standard', self.gridtype) + with pytest.raises(TypeError): + vw.truncate(dummy_solution['chi']) + + def test_truncate_different_shape(self): + # input to truncate of different shape should raise an error + solution = reference_solutions(self.interface, self.gridtype) + vw = solvers[self.interface](solution['uwnd'], solution['vwnd']) + with pytest.raises(ValueError): + vw.truncate(solution['chi'][:-1]) + + def test_truncate_unknown_grid(self): + # input to truncate with no identifiable grid should raise an error + solution = reference_solutions(self.interface, self.gridtype) + vw = solvers[self.interface](solution['uwnd'], solution['vwnd']) + solution['chi'].coords.update( + {'unknown': ('latitude', + solution['chi'].coords['latitude'].values)}) + solution['chi'] = solution['chi'].swap_dims({'latitude': 'unknown'}) + del solution['chi'].coords['latitude'] + with pytest.raises(ValueError): + vw.truncate(solution['chi']) diff --git a/test/windspharm/test_solution.py b/test/windspharm/test_solution.py new file mode 100644 index 00000000..ef817929 --- /dev/null +++ b/test/windspharm/test_solution.py @@ -0,0 +1,470 @@ +"""Test windspharm computations against reference solutions.""" +# Copyright (c) 2012-2018 Andrew Dawson +# +# Permission is hereby granted, free of charge, to any person obtaining a copy +# of this software and associated documentation files (the "Software"), to deal +# in the Software without restriction, including without limitation the rights +# to use, copy, modify, merge, publish, distribute, sublicense, and/or sell +# copies of the Software, and to permit persons to whom the Software is +# furnished to do so, subject to the following conditions: +# +# The above copyright notice and this permission notice shall be included in +# all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR +# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, +# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE +# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER +# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, +# OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN +# THE SOFTWARE. +from __future__ import absolute_import + +import numpy as np +import pytest + +# from windspharm.tests import VectorWindTest, solvers +from . import VectorWindTest, solvers +from .reference import reference_solutions + + +class SolutionTest(VectorWindTest): + """Base class for all solution test classes.""" + interface = None + gridtype = None + radius = None + legfunc = None + + @classmethod + def setup_class(cls): + msg = 'missing dependencies required to test the {!s} interface' + try: + cls.solution = reference_solutions(cls.interface, cls.gridtype) + except ValueError: + pytest.skip(msg.format(cls.interface)) + cls.pre_modify_solution() + try: + # gridtype argument only available for the standard interface + kwargs = {} + if cls.interface == 'standard': + kwargs['gridtype'] = cls.gridtype + if cls.radius is not None: + kwargs['rsphere'] = cls.radius + if cls.legfunc is not None: + kwargs['legfunc'] = cls.legfunc + cls.vw = solvers[cls.interface](cls.solution['uwnd'], + cls.solution['vwnd'], **kwargs) + except KeyError: + pytest.skip(msg.format(cls.interface)) + cls.post_modify_solution() + + @classmethod + def pre_modify_solution(cls): + pass + + @classmethod + def post_modify_solution(cls): + pass + + def test_magnitude(self): + # computed magnitude matches magnitude of reference solution? + mag1 = self.vw.magnitude() + mag2 = (self.solution['uwnd'] ** 2 + self.solution['vwnd'] ** 2) ** 0.5 + self.assert_error_is_zero(mag1, mag2) + + def test_vorticity(self): + # computed vorticity matches reference solution? + vrt1 = self.vw.vorticity() + vrt2 = self.solution['vrt'] + self.assert_error_is_zero(vrt1, vrt2) + + def test_divergence(self): + # computed divergence matches reference solution? + div1 = self.vw.divergence() + div2 = self.solution['div'] + self.assert_error_is_zero(div1, div2) + + def test_streamfunction(self): + # computed streamfunction matches reference solution? + sf1 = self.vw.streamfunction() + sf2 = self.solution['psi'].copy() + self.assert_error_is_zero(sf1, sf2) + + def test_velocitypotential(self): + # computed velocity potential matches reference solution? + vp1 = self.vw.velocitypotential() + vp2 = self.solution['chi'].copy() + self.assert_error_is_zero(vp1, vp2) + + def test_nondivergent(self): + # computed non-divergent vector wind matches reference solution? + upsi1, vpsi1 = self.vw.nondivergentcomponent() + upsi2, vpsi2 = self.solution['upsi'], self.solution['vpsi'] + self.assert_error_is_zero(upsi1, upsi2) + self.assert_error_is_zero(vpsi1, vpsi2) + + def test_irrotational(self): + # computed irrotational vector wind matches reference solution? + uchi1, vchi1 = self.vw.irrotationalcomponent() + uchi2, vchi2 = self.solution['uchi'], self.solution['vchi'] + self.assert_error_is_zero(uchi1, uchi2) + self.assert_error_is_zero(vchi1, vchi2) + + def test_gradient(self): + # computed gradient matches reference solution? + uchi1, vchi1 = self.vw.gradient(self.solution['chi']) + uchi2, vchi2 = self.solution['chigradu'], self.solution['chigradv'] + self.assert_error_is_zero(uchi1, uchi2) + self.assert_error_is_zero(vchi1, vchi2) + + def test_vrtdiv(self): + # vrtdiv() matches vorticity()/divergence()? + vrt1, div1 = self.vw.vrtdiv() + vrt2 = self.vw.vorticity() + div2 = self.vw.divergence() + self.assert_error_is_zero(vrt1, vrt2) + self.assert_error_is_zero(div1, div2) + + def test_sfvp(self): + # sfvp() matches streamfunction()/velocitypotential()? + sf1, vp1 = self.vw.sfvp() + sf2 = self.vw.streamfunction() + vp2 = self.vw.velocitypotential() + self.assert_error_is_zero(sf1, sf2) + self.assert_error_is_zero(vp1, vp2) + + def test_helmholtz(self): + # helmholtz() matches irrotationalcomponent()/nondivergentcomponent()? + uchi1, vchi1, upsi1, vpsi1 = self.vw.helmholtz() + uchi2, vchi2 = self.vw.irrotationalcomponent() + upsi2, vpsi2 = self.vw.nondivergentcomponent() + uchi, vchi, upsi2, vpsi2 = self.vw.helmholtz() + self.assert_error_is_zero(uchi1, uchi2) + self.assert_error_is_zero(vchi1, vchi2) + self.assert_error_is_zero(upsi1, upsi2) + self.assert_error_is_zero(vpsi1, vpsi2) + + def test_truncate(self): + # vorticity truncated to T21 matches reference? + vrt_trunc = self.vw.truncate(self.solution['vrt'], truncation=21) + self.assert_error_is_zero(vrt_trunc, self.solution['vrt_trunc']) + + +# ---------------------------------------------------------------------------- +# Tests for the standard interface + + +class StandardSolutionTest(SolutionTest): + """Base class for all standard interface solution test classes.""" + interface = 'standard' + + +class TestStandardRegular(StandardSolutionTest): + """Regular grid.""" + gridtype = 'regular' + + +class TestStandardGaussian(StandardSolutionTest): + """Gaussian grid.""" + gridtype = 'gaussian' + + +class TestStandardRegularSingleton(StandardSolutionTest): + """Singleton right-most dimension.""" + gridtype = 'regular' + + @classmethod + def pre_modify_solution(cls): + for field_name in cls.solution: + cls.solution[field_name] = \ + cls.solution[field_name][..., np.newaxis] + + +class TestStandardGaussianSingleton(StandardSolutionTest): + """Singleton right-most dimension.""" + gridtype = 'gaussian' + + @classmethod + def pre_modify_solution(cls): + for field_name in cls.solution: + cls.solution[field_name] = \ + cls.solution[field_name][..., np.newaxis] + + +class TestStandardMultiTime(StandardSolutionTest): + gridtype = 'regular' + + @classmethod + def pre_modify_solution(cls): + for field_name in cls.solution: + cls.solution[field_name] = \ + cls.solution[field_name][..., np.newaxis].repeat(5, axis=-1) + + +class TestStandardRadiusDefaultExplicit(StandardSolutionTest): + gridtype = 'regular' + radius = 6.3712e6 + + +class TestStandardRadius(StandardSolutionTest): + gridtype = 'regular' + radius = 6.3712e6 / 16. + + @classmethod + def post_modify_solution(cls): + # Divergence and vorticity should be scaled by the inverse of the + # radius factor. + cls.solution['vrt'] = cls.solution['vrt'] * 16. + cls.solution['div'] = cls.solution['div'] * 16. + cls.solution['vrt_trunc'] = cls.solution['vrt_trunc'] * 16 + # Stream function and velocity potential should be scaled by the + # radius factor. + cls.solution['psi'] = cls.solution['psi'] / 16 + cls.solution['chi'] = cls.solution['chi'] / 16 + + +class TestStandardLegfuncComputed(StandardSolutionTest): + """Computed Legendre functions.""" + gridtype = 'regular' + legfunc = 'computed' + + +# ---------------------------------------------------------------------------- +# Tests for the CDMS interface + + +class CDMSSolutionTest(SolutionTest): + """Base class for all CDMS interface solution test classes.""" + interface = 'cdms' + + def test_truncate_reversed(self): + # vorticity truncated to T21 matches reference? + vrt_trunc = self.vw.truncate(self.solution['vrt'][::-1], truncation=21) + self.assert_error_is_zero(vrt_trunc, self.solution['vrt_trunc']) + + +class TestCDMSRegular(CDMSSolutionTest): + """Regular grid.""" + gridtype = 'regular' + + +class TestCDMSGaussian(CDMSSolutionTest): + """Gaussian grid.""" + gridtype = 'gaussian' + + +class TestCDMSGridTranspose(CDMSSolutionTest): + gridtype = 'regular' + + @classmethod + def pre_modify_solution(cls): + for field_name in cls.solution: + cls.solution[field_name] = cls.solution[field_name].reorder('xy') + + def test_truncate_reversed(self): + # vorticity truncated to T21 matches reference? + vrt_trunc = self.vw.truncate(self.solution['vrt'][:, ::-1], + truncation=21) + self.assert_error_is_zero(vrt_trunc, self.solution['vrt_trunc']) + + +class TestCDMSInvertedLatitude(CDMSSolutionTest): + gridtype = 'regular' + + @classmethod + def pre_modify_solution(cls): + for field_name in ('uwnd', 'vwnd'): + cls.solution[field_name] = \ + cls.solution[field_name](latitude=(-90, 90)) + + @classmethod + def post_modify_solution(cls): + for field_name in ('uwnd', 'vwnd'): + cls.solution[field_name] = \ + cls.solution[field_name](latitude=(90, -90)) + + +class TestCDMSRadiusDefaultExplicit(CDMSSolutionTest): + gridtype = 'regular' + radius = 6.3712e6 + + +class TestCDMSRadius(CDMSSolutionTest): + gridtype = 'regular' + radius = 6.3712e6 / 16. + + @classmethod + def post_modify_solution(cls): + # Divergence and vorticity should be scaled by the inverse of the + # radius factor. + cls.solution['vrt'] = cls.solution['vrt'] * 16. + cls.solution['div'] = cls.solution['div'] * 16. + cls.solution['vrt_trunc'] = cls.solution['vrt_trunc'] * 16 + # Stream function and velocity potential should be scaled by the + # radius factor. + cls.solution['psi'] = cls.solution['psi'] / 16 + cls.solution['chi'] = cls.solution['chi'] / 16 + + +class TestCDMSLegfuncComputed(CDMSSolutionTest): + """Computed Legendre functions.""" + gridtype = 'regular' + legfunc = 'computed' + + +# ---------------------------------------------------------------------------- +# Tests for the Iris interface + + +class IrisSolutionTest(SolutionTest): + """Base class for all Iris interface solution test classes.""" + interface = 'iris' + + def test_truncate_reversed(self): + vrt_trunc = self.vw.truncate(self.solution['vrt'][::-1], truncation=21) + self.assert_error_is_zero(vrt_trunc, self.solution['vrt_trunc']) + + +class TestIrisRegular(IrisSolutionTest): + """Regular grid.""" + gridtype = 'regular' + + +class TestIrisGaussian(IrisSolutionTest): + """Gaussian grid.""" + gridtype = 'gaussian' + + +class TestIrisGridTranspose(IrisSolutionTest): + gridtype = 'regular' + + @classmethod + def pre_modify_solution(cls): + for field_name in cls.solution.keys(): + cls.solution[field_name].transpose([1, 0]) + + def test_truncate_reversed(self): + vrt_trunc = self.vw.truncate(self.solution['vrt'][:, ::-1], + truncation=21) + self.assert_error_is_zero(vrt_trunc, self.solution['vrt_trunc']) + + +class TestIrisInvertedLatitude(IrisSolutionTest): + gridtype = 'regular' + + @classmethod + def pre_modify_solution(cls): + for field_name in ('uwnd', 'vwnd'): + cls.solution[field_name] = cls.solution[field_name][::-1] + + @classmethod + def post_modify_solution(cls): + for field_name in ('uwnd', 'vwnd'): + cls.solution[field_name] = cls.solution[field_name][::-1] + + +class TestIrisRadiusDefaultExplicit(IrisSolutionTest): + gridtype = 'regular' + radius = 6.3712e6 + + +class TestIrisRadius(IrisSolutionTest): + gridtype = 'regular' + radius = 6.3712e6 / 16. + + @classmethod + def post_modify_solution(cls): + # Divergence and vorticity should be scaled by the inverse of the + # radius factor. + cls.solution['vrt'] = cls.solution['vrt'] * 16. + cls.solution['div'] = cls.solution['div'] * 16. + cls.solution['vrt_trunc'] = cls.solution['vrt_trunc'] * 16 + # Stream function and velocity potential should be scaled by the + # radius factor. + cls.solution['psi'] = cls.solution['psi'] / 16 + cls.solution['chi'] = cls.solution['chi'] / 16 + + +class TestIrisLegfuncComputed(IrisSolutionTest): + """Computed Legendre functions.""" + gridtype = 'regular' + legfunc = 'computed' + + +# ---------------------------------------------------------------------------- +# Tests for the Xarray interface + + +class XarraySolutionTest(SolutionTest): + """Base class for all Xarray interface solution test classes.""" + interface = 'xarray' + + def test_truncate_reversed(self): + vrt_trunc = self.vw.truncate(self.solution['vrt'][::-1], truncation=21) + self.assert_error_is_zero(vrt_trunc, self.solution['vrt_trunc']) + + +class TestXarrayRegular(XarraySolutionTest): + """Regular grid.""" + gridtype = 'regular' + + +class TestXarrayGaussian(XarraySolutionTest): + """Gaussian grid.""" + gridtype = 'gaussian' + + +class TestXarrayGridTranspose(XarraySolutionTest): + gridtype = 'regular' + + @classmethod + def pre_modify_solution(cls): + for field_name in cls.solution.keys(): + cls.solution[field_name] = cls.solution[field_name].transpose() + + def test_truncate_reversed(self): + vrt_trunc = self.vw.truncate(self.solution['vrt'][:, ::-1], + truncation=21) + self.assert_error_is_zero(vrt_trunc, self.solution['vrt_trunc']) + + +class TestXarrayInvertedLatitude(XarraySolutionTest): + gridtype = 'regular' + + @classmethod + def pre_modify_solution(cls): + for field_name in ('uwnd', 'vwnd'): + cls.solution[field_name] = cls.solution[field_name][::-1] + + @classmethod + def post_modify_solution(cls): + for field_name in ('uwnd', 'vwnd'): + cls.solution[field_name] = cls.solution[field_name][::-1] + + +class TestXarrayRadiusDefaultExplicit(XarraySolutionTest): + gridtype = 'regular' + radius = 6.3712e6 + + +class TestXarrayRadius(XarraySolutionTest): + gridtype = 'regular' + radius = 6.3712e6 / 16. + + @classmethod + def post_modify_solution(cls): + # Divergence and vorticity should be scaled by the inverse of the + # radius factor. + cls.solution['vrt'] = cls.solution['vrt'] * 16. + cls.solution['div'] = cls.solution['div'] * 16. + cls.solution['vrt_trunc'] = cls.solution['vrt_trunc'] * 16 + # Stream function and velocity potential should be scaled by the + # radius factor. + cls.solution['psi'] = cls.solution['psi'] / 16 + cls.solution['chi'] = cls.solution['chi'] / 16 + + +class TestXarrayLegfuncComputed(XarraySolutionTest): + """Computed Legendre functions.""" + gridtype = 'regular' + legfunc = 'computed' diff --git a/test/windspharm/test_tools.py b/test/windspharm/test_tools.py new file mode 100644 index 00000000..7a54af60 --- /dev/null +++ b/test/windspharm/test_tools.py @@ -0,0 +1,81 @@ +"""Tests for the `windspharm.tools` module.""" +# Copyright (c) 2012-2013 Andrew Dawson +# +# Permission is hereby granted, free of charge, to any person obtaining a copy +# of this software and associated documentation files (the "Software"), to deal +# in the Software without restriction, including without limitation the rights +# to use, copy, modify, merge, publish, distribute, sublicense, and/or sell +# copies of the Software, and to permit persons to whom the Software is +# furnished to do so, subject to the following conditions: +# +# The above copyright notice and this permission notice shall be included in +# all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR +# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, +# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE +# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER +# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, +# OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN +# THE SOFTWARE. +from __future__ import absolute_import + +import numpy as np +from numpy.testing import assert_array_equal + +# from windspharm.tools import (prep_data, recover_data, get_recovery, +# reverse_latdim, order_latdim) +from easyclimate.windspharm.tools import (prep_data, recover_data, get_recovery, + reverse_latdim, order_latdim) +# from windspharm.tests import VectorWindTest +from . import VectorWindTest + + +class TestTools(VectorWindTest): + """Tests for extra tools.""" + + def test_prep_recover_data(self): + # applying preparation and recovery should yield an identical data set + u = np.random.rand(12, 17, 73, 144) + up, uinfo = prep_data(u, 'tzyx') + ur = recover_data(up, uinfo) + assert_array_equal(u, ur) + + def test_get_recovery(self): + # recovery helper should produce the same result as the manual method + u = np.random.rand(12, 17, 73, 144) + up, uinfo = prep_data(u, 'tzyx') + ur1 = recover_data(up, uinfo) + recover = get_recovery(uinfo) + ur2, = recover(up) + assert_array_equal(ur1, ur2) + + def test_reverse_latdim(self): + # applying reversal to the latitude dimension twice should return it to + # its original + u = np.random.rand(12, 17, 73, 144) + v = np.random.rand(12, 17, 73, 144) + ur, vr = reverse_latdim(u, v, axis=2) + urr, vrr = reverse_latdim(ur, vr, axis=2) + assert_array_equal(u, urr) + assert_array_equal(v, vrr) + + def test_order_latdim(self): + # order_latdim should reverse a south-north latitude dimension + u = np.random.rand(12, 17, 73, 144) + v = np.random.rand(12, 17, 73, 144) + lat = np.arange(-90, 92.5, 2.5) + latr, ur, vr = order_latdim(lat, u, v, axis=2) + assert_array_equal(lat[::-1], latr) + assert_array_equal(u[:, :, ::-1], ur) + assert_array_equal(v[:, :, ::-1], vr) + + def test_order_latdim_null(self): + # order_latdim should not reverse a north-south latitude dimension + u = np.random.rand(12, 17, 73, 144) + v = np.random.rand(12, 17, 73, 144) + lat = np.arange(90, -92.5, -2.5) + latr, ur, vr = order_latdim(lat, u, v, axis=2) + assert_array_equal(lat, latr) + assert_array_equal(u, ur) + assert_array_equal(v, vr) diff --git a/test/windspharm/utils.py b/test/windspharm/utils.py new file mode 100644 index 00000000..627a8967 --- /dev/null +++ b/test/windspharm/utils.py @@ -0,0 +1,79 @@ +"""Utilities for constructing tests.""" +# Copyright (c) 2012-2016 Andrew Dawson +# +# Permission is hereby granted, free of charge, to any person obtaining a copy +# of this software and associated documentation files (the "Software"), to deal +# in the Software without restriction, including without limitation the rights +# to use, copy, modify, merge, publish, distribute, sublicense, and/or sell +# copies of the Software, and to permit persons to whom the Software is +# furnished to do so, subject to the following conditions: +# +# The above copyright notice and this permission notice shall be included in +# all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR +# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, +# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE +# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER +# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, +# OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN +# THE SOFTWARE. +from __future__ import absolute_import + +import numpy as np +try: + from iris.cube import Cube +except ImportError: + pass +try: + import xarray as xr +except ImportError: + try: + import xray as xr + except ImportError: + pass + + +def __tomasked(*args): + """Convert cdms2 variables or iris cubes to masked arrays. + + The conversion is safe, so if non-variables/cubes are passed they + are just returned. + + """ + def __asma(a): + try: + if isinstance(a, Cube): + # Retrieve the data from the cube. + a = a.data + except NameError: + pass + try: + # Retrieve data from cdms variable. + a = a.asma() + except AttributeError: + # The input is already an array or masked array, either extracted + # from an iris cube, or was like that to begin with. + pass + try: + if isinstance(a, xr.DataArray): + a = a.values + except NameError: + pass + return a + return [__asma(a) for a in args] + + +def error(a, b): + """Compute the error between two arrays. + + Computes RMSD normalized by the range of the second input. + + """ + a, b = __tomasked(a, b) + + return np.sqrt(((a - b)**2).mean()) / (np.max(b) - np.min(b)) + + +if __name__ == '__main__': + pass diff --git a/test_requirements.txt b/test_requirements.txt index 83cedabc..16b4fa12 100644 --- a/test_requirements.txt +++ b/test_requirements.txt @@ -4,6 +4,7 @@ geocat.viz <= 2023.10.0 cartopy >= 0.20 xeofs >= 1.1.0 matplotlib +pandas fast-barnes-py oceans intel-fortran-rt @@ -18,6 +19,8 @@ flox xarray-datatree xarray_regrid gsw_xarray +pooch +tqdm -pytest -pytest-cov \ No newline at end of file +pytest == 7.4.3 +pytest-cov == 4.1.0 \ No newline at end of file