From f842287b6bee40e14ecef156cbebb0b1e2224d8e Mon Sep 17 00:00:00 2001 From: "Sadie L. Bartholomew" Date: Fri, 16 Jan 2026 17:47:45 +0000 Subject: [PATCH 01/16] Set NEXTVERSION tags ready for v3.19.0 release --- cf/functions.py | 6 +++--- cf/read_write/read.py | 4 ++-- 2 files changed, 5 insertions(+), 5 deletions(-) diff --git a/cf/functions.py b/cf/functions.py index 6b1850a388..c4f7abe5fb 100644 --- a/cf/functions.py +++ b/cf/functions.py @@ -251,7 +251,7 @@ def configuration( The new display data option. The default is to not change the current behaviour. - .. versionadded:: NEXTVERSION + .. versionadded:: 3.19.0 regrid_logging: `bool` or `Constant`, optional The new value (either True to enable logging or False to @@ -1504,7 +1504,7 @@ def total_memory(): def is_log_level_info(logger): """Return True if and only if log level is at least as verbose as INFO. - Deprecated at version NEXTVERSION and is no longer available. Use + Deprecated at version 3.19.0 and is no longer available. Use `cfdm.is_log_level_info` instead. .. versionadded:: 3.16.3 @@ -1525,7 +1525,7 @@ def is_log_level_info(logger): _DEPRECATION_ERROR_FUNCTION( "is_log_level_info", message="Use cfdm.is_log_level_info instead", - version="NEXTVERSION", + version="3.19.0", removed_at="5.0.0", ) # pragma: no cover diff --git a/cf/read_write/read.py b/cf/read_write/read.py index 52cc132e7e..c920b72b97 100644 --- a/cf/read_write/read.py +++ b/cf/read_write/read.py @@ -318,7 +318,7 @@ class read(cfdm.read): {{read store_dataset_shards: `bool`, optional}} - .. versionadded:: NEXTVERSION + .. versionadded:: 3.19.0 {{read cfa: `dict`, optional}} @@ -334,7 +334,7 @@ class read(cfdm.read): {{read group_dimension_search: `str`, optional}} - .. versionadded:: NEXTVERSION + .. versionadded:: 3.19.0 umversion: deprecated at version 3.0.0 Use the *um* parameter instead. From 590ae8eafe441d1bc363499d7060cd19ed9fbe04 Mon Sep 17 00:00:00 2001 From: "Sadie L. Bartholomew" Date: Fri, 16 Jan 2026 17:53:18 +0000 Subject: [PATCH 02/16] Update __init__.py with new version, date and pins --- cf/__init__.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/cf/__init__.py b/cf/__init__.py index deaac4a813..23c973a29e 100644 --- a/cf/__init__.py +++ b/cf/__init__.py @@ -85,16 +85,16 @@ from packaging.version import Version -__date__ = "2025-10-16" -__version__ = "3.18.2" +__date__ = "2026-01-16" +__version__ = "3.19.0" __cf_version__ = cfdm.__cf_version__ __Conventions__ = f"CF-{__cf_version__}" # Check the version of cfdm (this is worth doing because of the very # tight coupling between cf and cfdm, and the risk of bad things # happening at run time if the versions are mismatched). -_minimum_vn = "1.12.3.1" -_maximum_vn = "1.12.4.0" +_minimum_vn = "1.13.0.0" +_maximum_vn = "1.13.1.0" _cfdm_vn = Version(cfdm.__version__) if _cfdm_vn < Version(_minimum_vn) or _cfdm_vn >= Version(_maximum_vn): raise RuntimeError( From ad286e816d18fa783b35a2cd3a2cc8ee2b4004e0 Mon Sep 17 00:00:00 2001 From: "Sadie L. Bartholomew" Date: Fri, 16 Jan 2026 17:59:19 +0000 Subject: [PATCH 03/16] Update info & setup files for new cfdm version pin --- Changelog.rst | 2 +- docs/source/installation.rst | 4 ++-- requirements.txt | 2 +- 3 files changed, 4 insertions(+), 4 deletions(-) diff --git a/Changelog.rst b/Changelog.rst index a5eefa214e..380799e41c 100644 --- a/Changelog.rst +++ b/Changelog.rst @@ -23,7 +23,7 @@ Version NEXTVERSION (https://github.com/NCAS-CMS/cf-python/issues/902) * Reduce the time taken to import `cf` (https://github.com/NCAS-CMS/cf-python/issues/902) -* Changed dependency: ``cfdm>=1.12.4.0, <1.12.5.0`` +* Changed dependency: ``cfdm>=1.13.0.0, <1.13.1.0`` ---- diff --git a/docs/source/installation.rst b/docs/source/installation.rst index 4ec36e0d67..13c9f0ef88 100644 --- a/docs/source/installation.rst +++ b/docs/source/installation.rst @@ -225,8 +225,8 @@ Required * `scipy `_, version 1.10.0 or newer. -* `cfdm `_, version 1.12.3.1 or up to, - but not including, 1.12.4.0. +* `cfdm `_, version 1.13.0.0 or up to, + but not including, 1.13.1.0. * `cfunits `_, version 3.3.7 or newer. diff --git a/requirements.txt b/requirements.txt index 13ddf3cf40..66cd2d0e36 100644 --- a/requirements.txt +++ b/requirements.txt @@ -1,7 +1,7 @@ netCDF4==1.7.2 cftime>=1.6.4 numpy>=2.0.0 -cfdm>=1.12.3.1, <1.12.4.0 +cfdm>=1.13.0.0, <1.13.1.0 psutil>=0.6.0 cfunits>=3.3.7 dask>=2025.5.1 From 67e88ca1f841e635ec81f99e6edfdd3cbad0b52a Mon Sep 17 00:00:00 2001 From: "Sadie L. Bartholomew" Date: Fri, 16 Jan 2026 18:07:20 +0000 Subject: [PATCH 04/16] Remove effective duplicate entry from long description --- README.md | 1 - docs/source/introduction.rst | 2 -- setup.py | 2 -- 3 files changed, 5 deletions(-) diff --git a/README.md b/README.md index 140e6faeb9..44f2c3b51c 100644 --- a/README.md +++ b/README.md @@ -96,7 +96,6 @@ of its array manipulation and can: * test whether two field constructs are the same, * modify field construct metadata and data, * create subspaces of field constructs, -* write field constructs to netCDF datasets on disk, * incorporate, and create, metadata stored in external files, * read, write, and create data that have been compressed by convention (i.e. ragged or gathered arrays, or coordinate arrays compressed by diff --git a/docs/source/introduction.rst b/docs/source/introduction.rst index 9e9b3208a9..537ff8a088 100644 --- a/docs/source/introduction.rst +++ b/docs/source/introduction.rst @@ -93,8 +93,6 @@ The `cf` package can: * create subspaces of field constructs, -* write field constructs to netCDF datasets on disk, - * incorporate, and create, metadata stored in external files, * read, write, and create data that have been compressed by convention diff --git a/setup.py b/setup.py index 2ed416976d..4620e64a4a 100755 --- a/setup.py +++ b/setup.py @@ -198,8 +198,6 @@ def compile(): * create subspaces of field constructs, -* write field constructs to netCDF datasets on disk, - * incorporate, and create, metadata stored in external files, * read, write, and create data that have been compressed by convention From fc886eba08d5f380cf569a0f7ce4cd6fe273bcc6 Mon Sep 17 00:00:00 2001 From: "Sadie L. Bartholomew" Date: Fri, 16 Jan 2026 18:10:11 +0000 Subject: [PATCH 05/16] Update & improve order in v3.19.0 changelog entry --- Changelog.rst | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/Changelog.rst b/Changelog.rst index 380799e41c..32211c82d7 100644 --- a/Changelog.rst +++ b/Changelog.rst @@ -1,7 +1,7 @@ -Version NEXTVERSION +Version 3.19.0 -------------- -**2026-01-??** +**2026-01-16** * Write Zarr v3 datasets with `cf.write` (https://github.com/NCAS-CMS/cf-python/issues/895) @@ -9,7 +9,6 @@ Version NEXTVERSION `cf.read` (https://github.com/NCAS-CMS/cf-python/issues/894) * Reduce the time taken to import `cf` (https://github.com/NCAS-CMS/cf-python/issues/902) -* New optional dependency: ``zarr>=3.1.3`` * New function to control the creation of cached elements during data display: `cf.display_data` (https://github.com/NCAS-CMS/cf-python/issues/913) @@ -23,6 +22,7 @@ Version NEXTVERSION (https://github.com/NCAS-CMS/cf-python/issues/902) * Reduce the time taken to import `cf` (https://github.com/NCAS-CMS/cf-python/issues/902) +* New optional dependency: ``zarr>=3.1.3`` * Changed dependency: ``cfdm>=1.13.0.0, <1.13.1.0`` ---- From a2ac17ed036ee0d9052ff5567ce7fd7c37b36071 Mon Sep 17 00:00:00 2001 From: "Sadie L. Bartholomew" Date: Fri, 16 Jan 2026 18:23:05 +0000 Subject: [PATCH 06/16] Add all undoc'd new methods, all for cf.Data class --- docs/source/class/cf.Data.rst | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/docs/source/class/cf.Data.rst b/docs/source/class/cf.Data.rst index b89f262de9..41b4c787b7 100644 --- a/docs/source/class/cf.Data.rst +++ b/docs/source/class/cf.Data.rst @@ -659,6 +659,11 @@ Performance ~cf.Data.fits_in_memory ~cf.Data.section ~cf.Data.persist + ~cf.Data.cache_elements + ~cf.Data.get_cached_elements + ~cf.Data.nc_clear_dataset_shards + ~cf.Data.nc_dataset_shards + ~cf.Data.nc_set_dataset_shards .. rubric:: Attributes From 792e054bcc3836103ad9be813000c18009b85b43 Mon Sep 17 00:00:00 2001 From: "Sadie L. Bartholomew" Date: Fri, 16 Jan 2026 18:26:23 +0000 Subject: [PATCH 07/16] Create new docs release link on releases page + typo fix --- docs/source/releases.rst | 7 ++++++- 1 file changed, 6 insertions(+), 1 deletion(-) diff --git a/docs/source/releases.rst b/docs/source/releases.rst index ff5fe0923b..48ddbd85f2 100644 --- a/docs/source/releases.rst +++ b/docs/source/releases.rst @@ -13,10 +13,15 @@ Documentation for all versions of cf. :local: :backlinks: entry +**CF-1.13** +----------- + +* `Version 3.19.0 `_ (2026-01-16) + **CF-1.12** ----------- -* `Version 3.18.2 `_ (2025-10-16) +* `Version 3.18.2 `_ (2025-10-16) * `Version 3.18.1 `_ (2025-08-20) * `Version 3.18.0 `_ (2025-06-05) * `Version 3.17.0 `_ (2025-04-02) From 41eec1aacb99c31491ae71065e2b3bcd9c44e981 Mon Sep 17 00:00:00 2001 From: "Sadie L. Bartholomew" Date: Fri, 16 Jan 2026 18:30:41 +0000 Subject: [PATCH 08/16] Update RELEASE checklist to remove outdated step + typo fix --- RELEASE.md | 9 ++------- 1 file changed, 2 insertions(+), 7 deletions(-) diff --git a/RELEASE.md b/RELEASE.md index 0e6713120e..f76e96fcc4 100644 --- a/RELEASE.md +++ b/RELEASE.md @@ -26,8 +26,7 @@ - [ ] Ensure that the requirements on dependencies & their versions are up-to-date and consistent in both the `requirements.txt` and in `docs/source/installation.rst` (paying particular attention to - `cfdm`); and in the `_requires` list and `Version` checks in - `cf/__init__.py`. + `cfdm`); and in the cfdm-only `Version` checks in `cf/__init__.py`. - [ ] Make sure that `README.md` is up to date. @@ -40,7 +39,7 @@ deprecated methods and keyword arguments that can be completely removed, i.e. those with a ``removed_at`` version that is at or before the version being released. Remove any reference to them in - the method, class, or fucntion (including, if appropriate, the + the method, class, or function (including, if appropriate, the ``@_deprecated_kwarg_check`` decorator), and remove them from the relevant documentation ``.rst`` files. @@ -67,10 +66,6 @@ ./test_tutorial_code ``` -- [ ] **Follow all of the steps outlined externally in [`RECIPES.md`](./RECIPES.md)**, - notably so that the correct Sphinx-related environment is prepared for - documentation building. - - [ ] Ensure that the [PDF for Cheat Sheet](docs/_downloads/cheatsheet.pdf) is updated to include any API changes. The PDF is created using Canva keeping in mind the colours and fonts of the website. The same could From ecdc483bbc16f1ba87266b9b9355be9c87b5f7d6 Mon Sep 17 00:00:00 2001 From: "Sadie L. Bartholomew" Date: Fri, 16 Jan 2026 18:42:33 +0000 Subject: [PATCH 09/16] Add sharding to spelling checker false positives --- docs/source/spelling_false_positives.txt | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/docs/source/spelling_false_positives.txt b/docs/source/spelling_false_positives.txt index 656661c165..36d5ada565 100644 --- a/docs/source/spelling_false_positives.txt +++ b/docs/source/spelling_false_positives.txt @@ -449,8 +449,9 @@ seterr setitem setmask setprop -Signell shapefiles +sharding +Signell sinh spacings spatiotemporal From e1ef2d0ce8708ce5a37ff3d5ef5c6b36ac40d1d4 Mon Sep 17 00:00:00 2001 From: "Sadie L. Bartholomew" Date: Fri, 16 Jan 2026 22:39:52 +0000 Subject: [PATCH 10/16] Prevent Sphinx bulid errors on section can't 'begin w/ a transition' --- docs/source/aggregation_rules.rst | 2 -- docs/source/api_reference.rst | 2 -- 2 files changed, 4 deletions(-) diff --git a/docs/source/aggregation_rules.rst b/docs/source/aggregation_rules.rst index c9732ee35e..6015cc4952 100644 --- a/docs/source/aggregation_rules.rst +++ b/docs/source/aggregation_rules.rst @@ -6,8 +6,6 @@ **Aggregation rules** ===================== ----- - Version |release| for version |version| of the CF conventions. *David Hassell and Jonathan Gregory (2012)* diff --git a/docs/source/api_reference.rst b/docs/source/api_reference.rst index 24a60e03de..af73b1a377 100644 --- a/docs/source/api_reference.rst +++ b/docs/source/api_reference.rst @@ -4,8 +4,6 @@ **API reference** ================= ----- - Version |release| for version |version| of the CF conventions. * **Construct classes** From 27af3e50596554d24bc447ead434a24c56ba47aa Mon Sep 17 00:00:00 2001 From: "Sadie L. Bartholomew" Date: Fri, 16 Jan 2026 22:46:12 +0000 Subject: [PATCH 11/16] Prevent more bulid errors on section can't 'begin w/ a transition' --- Changelog.rst | 4 ---- docs/source/cheat_sheet.rst | 2 -- 2 files changed, 6 deletions(-) diff --git a/Changelog.rst b/Changelog.rst index 32211c82d7..80e7f1efd6 100644 --- a/Changelog.rst +++ b/Changelog.rst @@ -1214,7 +1214,6 @@ version 2.2.4 version 2.2.3 -------------- ----- **2018-06-21** @@ -1253,7 +1252,6 @@ version 2.2.1 version 2.2.0 ------------- ----- **2018-06-04** @@ -1320,7 +1318,6 @@ version 2.1.3 version 2.1.2 ------------- ----- **2017-11-28** @@ -1353,7 +1350,6 @@ version 2.1 version 2.0.6 ------------- ----- **2017-09-28** diff --git a/docs/source/cheat_sheet.rst b/docs/source/cheat_sheet.rst index 8a9b8dead7..dce63c93e0 100644 --- a/docs/source/cheat_sheet.rst +++ b/docs/source/cheat_sheet.rst @@ -7,8 +7,6 @@ **Cheat Sheet** =============== ----- - Version |release| for version |version| of the CF conventions. This cheat sheet provides a summary of some key functions and methods in From 966267af8e0659aca834e4bd7ddd5c763b0352c6 Mon Sep 17 00:00:00 2001 From: "Sadie L. Bartholomew" Date: Fri, 16 Jan 2026 23:53:27 +0000 Subject: [PATCH 12/16] Remove deprecated function from functions listing --- docs/source/function.rst | 2 -- 1 file changed, 2 deletions(-) diff --git a/docs/source/function.rst b/docs/source/function.rst index 8892f325c1..31de3bce8e 100644 --- a/docs/source/function.rst +++ b/docs/source/function.rst @@ -20,7 +20,6 @@ Reading and writing cf.read cf.write - cf.netcdf_lock Aggregation ----------- @@ -157,7 +156,6 @@ Active storage reductions cf.active_storage cf.active_storage_url cf.active_storage_max_requests - cf.netcdf_lock Miscellaneous ------------- From 6926312da2c73d5d06841e3ce141224175abc0ac Mon Sep 17 00:00:00 2001 From: "Sadie L. Bartholomew" Date: Mon, 19 Jan 2026 13:35:21 +0000 Subject: [PATCH 13/16] Update date of v3.19.0 release to 19th Jan --- Changelog.rst | 2 +- cf/__init__.py | 2 +- docs/source/releases.rst | 2 +- 3 files changed, 3 insertions(+), 3 deletions(-) diff --git a/Changelog.rst b/Changelog.rst index 80e7f1efd6..efef383fb9 100644 --- a/Changelog.rst +++ b/Changelog.rst @@ -1,7 +1,7 @@ Version 3.19.0 -------------- -**2026-01-16** +**2026-01-19** * Write Zarr v3 datasets with `cf.write` (https://github.com/NCAS-CMS/cf-python/issues/895) diff --git a/cf/__init__.py b/cf/__init__.py index 23c973a29e..125aa6286f 100644 --- a/cf/__init__.py +++ b/cf/__init__.py @@ -85,7 +85,7 @@ from packaging.version import Version -__date__ = "2026-01-16" +__date__ = "2026-01-19" __version__ = "3.19.0" __cf_version__ = cfdm.__cf_version__ __Conventions__ = f"CF-{__cf_version__}" diff --git a/docs/source/releases.rst b/docs/source/releases.rst index 48ddbd85f2..74986da574 100644 --- a/docs/source/releases.rst +++ b/docs/source/releases.rst @@ -16,7 +16,7 @@ Documentation for all versions of cf. **CF-1.13** ----------- -* `Version 3.19.0 `_ (2026-01-16) +* `Version 3.19.0 `_ (2026-01-19) **CF-1.12** ----------- From 8640aa425627fc08458239f0cb8352bfe534f910 Mon Sep 17 00:00:00 2001 From: "Sadie L. Bartholomew" Date: Mon, 19 Jan 2026 17:25:16 +0000 Subject: [PATCH 14/16] Update recipes-docs conf.py esp. to separate source & gallery --- recipes-docs/source/conf.py | 15 +++++++++++++-- 1 file changed, 13 insertions(+), 2 deletions(-) diff --git a/recipes-docs/source/conf.py b/recipes-docs/source/conf.py index e5bfa04418..cc6782ef0b 100644 --- a/recipes-docs/source/conf.py +++ b/recipes-docs/source/conf.py @@ -1,8 +1,11 @@ +import faulthandler import os import sys from sphinx_gallery.sorting import FileNameSortKey +faulthandler.enable() # seg fault detection - as recipes prone to seg faulting + # Make main 'docs' conf.py importable sys.path.insert(0, os.path.abspath("../../docs/source")) @@ -19,8 +22,9 @@ # sphinx-gallery configuration sphinx_gallery_conf = { - "examples_dirs": "recipes", # path to recipe files + "examples_dirs": "recipes-source", # path to recipe files "gallery_dirs": "recipes", # path to save gallery generated output + "ignore_pattern": "/exclusions/", "run_stale_examples": False, # Below setting can be buggy: see: # https://github.com/sphinx-gallery/sphinx-gallery/issues/967 @@ -29,15 +33,22 @@ "doc_module": ("cf",), "inspect_global_variables": True, "within_subsection_order": FileNameSortKey, - "default_thumb_file": "_static/cf-recipe-placeholder-squarecrop.png", + "default_thumb_file": "../../docs/source/_static/cf-recipe-placeholder-squarecrop.png", "image_scrapers": ( "matplotlib", ), # Ensures Matplotlib images are captured "plot_gallery": True, # Enables plot rendering "reset_modules": ("matplotlib",), # Helps with memory management "capture_repr": (), + # "filename_pattern": r"plot", } +exclude_patterns = [ + "exclusions/**", +] + html_static_path = ["../../docs/source/_static"] html_logo = "../../docs/source/images/logo.svg" html_favicon = "../../docs/source/_static/favicon.ico" + +templates_path = ["../../docs/_templates"] From c211313efe05497608c8c7fc31b2b9b05a3d2264 Mon Sep 17 00:00:00 2001 From: "Sadie L. Bartholomew" Date: Mon, 19 Jan 2026 17:34:38 +0000 Subject: [PATCH 15/16] Update docs resources: move logos from templates to static dir. --- docs/_templates/sponsors.html | 8 ++++---- docs/{_templates => source/_static}/logo_EC.png | Bin docs/{_templates => source/_static}/logo_ERC.png | Bin docs/{_templates => source/_static}/logo_NCAS.png | Bin docs/{_templates => source/_static}/logo_NERC.png | Bin 5 files changed, 4 insertions(+), 4 deletions(-) rename docs/{_templates => source/_static}/logo_EC.png (100%) rename docs/{_templates => source/_static}/logo_ERC.png (100%) rename docs/{_templates => source/_static}/logo_NCAS.png (100%) rename docs/{_templates => source/_static}/logo_NERC.png (100%) diff --git a/docs/_templates/sponsors.html b/docs/_templates/sponsors.html index e5e4d4f6be..0ca7485220 100644 --- a/docs/_templates/sponsors.html +++ b/docs/_templates/sponsors.html @@ -16,8 +16,8 @@ and by NCAS.

- - - + + +
- + diff --git a/docs/_templates/logo_EC.png b/docs/source/_static/logo_EC.png similarity index 100% rename from docs/_templates/logo_EC.png rename to docs/source/_static/logo_EC.png diff --git a/docs/_templates/logo_ERC.png b/docs/source/_static/logo_ERC.png similarity index 100% rename from docs/_templates/logo_ERC.png rename to docs/source/_static/logo_ERC.png diff --git a/docs/_templates/logo_NCAS.png b/docs/source/_static/logo_NCAS.png similarity index 100% rename from docs/_templates/logo_NCAS.png rename to docs/source/_static/logo_NCAS.png diff --git a/docs/_templates/logo_NERC.png b/docs/source/_static/logo_NERC.png similarity index 100% rename from docs/_templates/logo_NERC.png rename to docs/source/_static/logo_NERC.png From fc335a22cc0e31b0cebcdac4f4a6712b07b15ecb Mon Sep 17 00:00:00 2001 From: "Sadie L. Bartholomew" Date: Mon, 19 Jan 2026 17:45:08 +0000 Subject: [PATCH 16/16] Add new recipes-source for recipe scripts & context RST only --- recipes-docs/source/index.rst | 92 +----- recipes-docs/source/recipes-source/README.rst | 64 +++++ recipes-docs/source/recipes-source/cfplot.png | Bin 0 -> 83646 bytes .../source/recipes-source/gallery.rst | 121 ++++++++ recipes-docs/source/recipes-source/index.rst | 7 + .../source/recipes-source/plot_01_recipe.py | 77 +++++ .../source/recipes-source/plot_02_recipe.py | 96 +++++++ .../source/recipes-source/plot_03_recipe.py | 35 +++ .../source/recipes-source/plot_04_recipe.py | 40 +++ .../source/recipes-source/plot_05_recipe.py | 64 +++++ .../source/recipes-source/plot_06_recipe.py | 65 +++++ .../source/recipes-source/plot_07_recipe.py | 65 +++++ .../source/recipes-source/plot_08_recipe.py | 143 ++++++++++ .../source/recipes-source/plot_09_recipe.py | 82 ++++++ .../source/recipes-source/plot_10_recipe.py | 96 +++++++ .../source/recipes-source/plot_11_recipe.py | 92 ++++++ .../source/recipes-source/plot_12_recipe.py | 117 ++++++++ .../source/recipes-source/plot_13_recipe.py | 268 ++++++++++++++++++ .../source/recipes-source/plot_14_recipe.py | 181 ++++++++++++ .../source/recipes-source/plot_15_recipe.py | 163 +++++++++++ .../source/recipes-source/plot_16_recipe.py | 57 ++++ .../source/recipes-source/plot_17_recipe.py | 109 +++++++ .../source/recipes-source/plot_18_recipe.py | 143 ++++++++++ .../source/recipes-source/plot_19_recipe.py | 108 +++++++ .../source/recipes-source/plot_20_recipe.py | 97 +++++++ .../source/recipes-source/plot_22_recipe.py | 121 ++++++++ .../source/recipes-source/plot_23_recipe.py | 174 ++++++++++++ 27 files changed, 2588 insertions(+), 89 deletions(-) create mode 100644 recipes-docs/source/recipes-source/README.rst create mode 100644 recipes-docs/source/recipes-source/cfplot.png create mode 100644 recipes-docs/source/recipes-source/gallery.rst create mode 100644 recipes-docs/source/recipes-source/index.rst create mode 100644 recipes-docs/source/recipes-source/plot_01_recipe.py create mode 100644 recipes-docs/source/recipes-source/plot_02_recipe.py create mode 100644 recipes-docs/source/recipes-source/plot_03_recipe.py create mode 100644 recipes-docs/source/recipes-source/plot_04_recipe.py create mode 100644 recipes-docs/source/recipes-source/plot_05_recipe.py create mode 100644 recipes-docs/source/recipes-source/plot_06_recipe.py create mode 100644 recipes-docs/source/recipes-source/plot_07_recipe.py create mode 100644 recipes-docs/source/recipes-source/plot_08_recipe.py create mode 100644 recipes-docs/source/recipes-source/plot_09_recipe.py create mode 100644 recipes-docs/source/recipes-source/plot_10_recipe.py create mode 100644 recipes-docs/source/recipes-source/plot_11_recipe.py create mode 100644 recipes-docs/source/recipes-source/plot_12_recipe.py create mode 100644 recipes-docs/source/recipes-source/plot_13_recipe.py create mode 100644 recipes-docs/source/recipes-source/plot_14_recipe.py create mode 100644 recipes-docs/source/recipes-source/plot_15_recipe.py create mode 100644 recipes-docs/source/recipes-source/plot_16_recipe.py create mode 100644 recipes-docs/source/recipes-source/plot_17_recipe.py create mode 100644 recipes-docs/source/recipes-source/plot_18_recipe.py create mode 100644 recipes-docs/source/recipes-source/plot_19_recipe.py create mode 100644 recipes-docs/source/recipes-source/plot_20_recipe.py create mode 100644 recipes-docs/source/recipes-source/plot_22_recipe.py create mode 100644 recipes-docs/source/recipes-source/plot_23_recipe.py diff --git a/recipes-docs/source/index.rst b/recipes-docs/source/index.rst index ee2bddf946..3bbd15b42f 100644 --- a/recipes-docs/source/index.rst +++ b/recipes-docs/source/index.rst @@ -1,95 +1,9 @@ -:orphan: - **Recipes using cf** ==================== ---- -Version |release| for version |version| of the CF conventions. - -Click on the keywords below to filter the recipes according to their function: - -.. raw:: html - - - -.. raw:: html - - - -
- - - - - - - - - - - - - -
- -.. raw:: html - - -.. raw:: html - -
- - -.. raw:: html - -
- - -.. only:: html - - .. container:: sphx-glr-footer sphx-glr-footer-gallery - - .. container:: sphx-glr-download sphx-glr-download-python - - :download:`Download all examples in Python source code: recipes_python.zip ` - - .. container:: sphx-glr-download sphx-glr-download-jupyter - - :download:`Download all examples in Jupyter notebooks: recipes_jupyter.zip ` - - -.. only:: html - - .. rst-class:: sphx-glr-signature +.. toctree:: + :maxdepth: 2 - `Gallery generated by Sphinx-Gallery `_ + recipes-source/gallery diff --git a/recipes-docs/source/recipes-source/README.rst b/recipes-docs/source/recipes-source/README.rst new file mode 100644 index 0000000000..69ef485e26 --- /dev/null +++ b/recipes-docs/source/recipes-source/README.rst @@ -0,0 +1,64 @@ +**Recipes using cf** +==================== + +---- + +Version |release| for version |version| of the CF conventions. + +Click on the keywords below to filter the recipes according to their function: + +.. raw:: html + + + +.. raw:: html + + + +
+ + + + + + + + + + + + + +
+ +.. raw:: html + diff --git a/recipes-docs/source/recipes-source/cfplot.png b/recipes-docs/source/recipes-source/cfplot.png new file mode 100644 index 0000000000000000000000000000000000000000..9e4718c22df5bccd5df80342a563e6c53bd6b973 GIT binary patch literal 83646 zcmce;cUY5Y6F!Q)Vn;^XK8ZcGs1+JY}Ald+xdCeR^G0fpzPltt>1otXRdX z8Z0d9i7YJZq&IJbzmW~&Ul0EgcaXc`plM^~;C$2WHjBzl2irR~4tFeXaXH<#v$wRd z77-8^5I)Ui;ox9vFCi#M_;Z1Pjh(sRG1YaQ@Q^LGihA}eEIV(af7YbRq*$`7VPU~u zmDRczGuGpLPs_2Z@~2xfb~Dx+voGz{Gp@hZomN!;XDxSThq_M1-y`Yx9_`Fi&ncNj zBU^90;m_IXm-@ZSh3%^PF3YvVT^foOSj8!~l?!pLSmA-XYs3Z5tY@RaA}LG4_bsgB zhQtL-{CC4B;OCR~lk0;2{*gsxo$Q`}|9J6oFq$g-1nt7Wz~JZI@%7l%Kd}7R{QvJ| zyVf)2##>B+q>}flUXUH973f~O=D%m}-ZrwajN!d64_$|1owbaNj7Bd;)Fx#a7B70N zM2eW#4A@!rmbiEtXGl*E)`!*(G}lFncu2^Ih-j*-tJ@VC6wNd%#~K>E8115D7zjk@ zmb1`qM;M-jXMFsqNTB7PX|1)~C%w?Wx6nX9C00W0GoA0u8STm5QZa|`pK#<-#zb+^ z!eF$R>x{*T6DO+Sbz>dbtVYwLuFb7?x>1IF zJn+KhGy3^jqwORy>+TC43j-0>K|RIJWTSguHq^fmx}tF8?t@zqbUdjna&3RAaBwz>8CgNl8h(=g#S{v$ID@dla1`bZJd47<(=xRR$bRX(#Dt+^&5t z@%72^nQr^WLMf-JKn~ZrPIJe_;gk#*j#bm^hs4lO$L_mT!-Yd}p+c48E`1*K*ojiQ zkxI-3bq-0Vu-X9)0|SFtyl|z(r_1TiMCtU3rSHbsLNK;6w(wtT_=In&-iA zf=LcPIUacZ@eU3xEv-ajFQzsucuR@PwCG>&SV%YuUu!LE+QRfuMrP(U;!GEwJlv$f zseR|pohQ0&7q3N$V#US9+vYr0m~O%oU0P7W{cZosqQGnzEMhd!-fq;hBnJf z3q_xOwrT0;Bp>fxST->+kVDd{ z_24b;xIYN`#Rp?$4**5fIo|K_jzC)MIIr>5-C)mF58lJw&D<6utO89Joq zD?6pM#^+&1I@3dqx(*I`J7kt@S=oiJ(((x{L-Af&E6hnRR~YVW4L|vmuy9V7i*v6_ z%6~q7!8}Gu1y*1Gc^Y-zCQhK}J2z)9JIiL_bu61-g!Ai0y?uM(n_XdD!1-Q_ zt=qQYIdUyqPKuaRZg0uCU5g`Qg$;`we!6iK^f*O_m<~3?j?yV3BF+;A$E&-YXawn4 z8Q8t!)p<63y^6v%Wedr<+pdbu#Xs97@RqAp*;42E8U6OmUTQ(FYerXJ?Ufzkm{w)0 zeLk~}0(^YJmX?--pRVk5WGu=ZFf^SgU-6u-=81m+gRfE)=83y^E0#(ta?E)8?3oC$ zE=Jrwy~lj2zgnIcDed7(nnp_mZ?}h%H)~#zA8pP__+6h#nUJV=AIl4W}PcH~hh*MrzWu(2$s2IJ+a1NLwbBD=f0tr!u4#@1kuw=t1mUFq$I z;3OK^!+J-dNx-e_*puhZ<@Yvbl$3PO47XG{jCC{}xnTP)fxedI^3(QnESq`?*!s)f7F0=H*x5#Ncy?CNXnK<+J zHo*t?v{}Z=A3OJX6rjI{YlOx!czlOpFRwc}ImL!cfj7v^DxQ$waPGcq)tPhKml_R5 z;OFNjYF2j;c||sciHXVV{I?)}X=x*LwVl0vMlAia?@lG+@-n?V!@Mz$jRA)?eYhnh z<8oibie0wfQ>dZA))*$&MDjapWiaZvsz}k?o^5rpS)?u3x zb>hratsoWmjrDx@OPWo&VM>ghUvo=~)M5#E;Rf-cBg|#S~)%4m}~%Pu8|I=i@@ z9T03+lVK(b7#7`Y-6(AOaTl2^0~1Y3OgtqZpjzriZ-=ExH7s@lAA7#fJz=lfMfK|H zYOA*N@JfI71ayeo?^YfxU7D(m629%{Iy;i|;KBM6Cr^%m>-q#&f0;l17*3oK(E)~F z0lp1+$^C*xnI$D9ekM7&xv-5%aMs&hR;Wg9q_18zcu(m@DSAwx)ZC2|JUkh4&-NoX zdGy}rH5?agU$#4rMP9BvC^gq%atQwS>$h)?-(Nk}F=&jS!@^4yO%#o`rfKp}ESfGe z7iU`%f~5IOu9b>pORp^1$9gUfvx&%}oiAS+5bP?nH@{*-Qginx#|;0h;?!MPX12YS z+vYghu99osnDv4u=T^nmSi2fttDYhX1V*u&f3IhPqoZVf`o zu{bmQZElhV<91y7S*1(_t1h$>Ng3ikl9G#tPMx={`pQZnCbpG(t>9#@Xd8RXhF1o0l3QC_{m8U|2;(fk1yqUr zk&F4PhtGXi^I8n;$|r!2#SJHhsx~*utemcnI^zHvmI@Kj@7XiWprcYV|Lj$JTe|XE z?wLcA+<`QkzA|Jld-m)Rwf%gl2f_%SfWX{$X3mB9j3!W@vs1ikQ{qqDzPD)NrH!k;-u>I@}R%DP}o|>x2=`j<6 z$xz8Z>DDJXS@P}UA&;AlE<~--#=7FkGX832mzwl~DGH|_cFJf$2vVRU)1fiuKP(q( zQ{y>UEVM4`Y}>$6-_53La1#>=funs0pBmVn4;ulF6>yDNYpVJn3{-RzheL3PDcJe# zoSd9=4x56rhVDNrPM$njIWtOaycQ{3am1*kA3G$Ki^Yaokad!fNpPyfZL&D)K5!@`^U#N0vSj|hASO;!QZ@dEzP?MrC zr|0}D8XjB$L=4!>bSIDE+X=AQXi#QL^FNc}e99Bq&s)dv^7Cgn#lmR@=L?=jsHv&R zZ038@>c~^#J&PJA8<&}k`rQggV4AJ1S|VTOaSS2OlKChM1zR7Tp_-BG;Z@WBD=+4Dc0a*)y-%ETN5FmEt(RLwlilP8x z?Bt@U?>BMYahmKkyjyA>=aT%j(E{9RA;1$p&oejPygS%ps~lFE45w3@CRpt-(pnwk z(2PBK?AUw%!{-~W#a%3X&5D-oD|5IbJ3)LruchvLQ-aFn%RY_cSHW-+LrfbkPG5u9 zhqATcSakf5oos(OxC2l1nS&VW(= zaH1^b?R1Ap5Ic4d{!RFC!{t+#8F=GRbQT!&N&5GDESuSUFf5yc!7XOBLvLa(zV$bwY9bN z043nGVM9;0;o^Blwi4mwNWhsJXX0|r0+2ob{P{B(Pt6*PUs*D*z{7;lK2ppm2Q}G~ zTbEvy0}|@C{lE(IxOmT&V(;!kPVXHG*EPD!kUoRUOKl&Ec>S8&zWPK2Y&ZCR@*M7p zY>1335z~^Ourn7@F5FNo*t|TkaCvbwqy9&GweM>AExTY&#$zN~{9_ zHByqr(f&4igT*$};F+(_1a%%5%E`&GF}7{p%D0@U#Vpfrd$46MSajdYa>08>SXkJL zG7rWf>>%_ZZp=#ZwRLI9Kd3R{ zfkrMcyh$;wCgO&|A|j$?eK$rvd9t}6UjIt>2JCW~vDms}CjQ`LZpE&++N^lvrXSyx zyvBa+&JB{il9JN36xCQ*EC~2IWAuU2oJo6tHdKj}V?YszK>VYuC0*Kt1-R!*vLEvc zx%uId_BKC)(-V~U1WBJa$z-Y7p_$>?cMrhW$oBJNjUbu5Yh3n7-(iNEioCo$Y0}b3JGi)JpNT)asKs@XSnCjgTR3}wq=Xo@e-gf};nI2$vYutlt>qV6Od12yP=cQ^i*HF*VKt)XHcE#k z;zYbLyl2j+@Wj~q)YR7ESr3P=n%MXoO}Dv91c65>U|U9jEZYHQyxJyMavf5rHVE{* z)vkHDJ)q_P%Oms^7H(< zSz4;9fq~MC78E$zmbP4THIh*&4je9dcF{5=R%dci%EK%FPUma5Wca1;j%IZQ z4((BiP#(1`3b~}j0Kq4OZmn8VapYnH0ai|FgMonoHhZPmW~~xX(nGwoAdaVncJXYYMR*V zI(h~WXybU{V-BYqq+=CCKC=!jrYn^$bq3PJual2F4788S`=`!I)h|!)-@7*fKyp^Q z_$)J$0W8KWoJk8vRM?D%02$_f6;ng-N&`5JS>!Vql2dT9Gv2>d!=E$)f@-Qg^F?o9u}67 zR-m<*fQf3VtM4*g+g9$z!pcD+ARd%e3}g_1VaMU-it*98H)hUPE|rR>Bqm;kK|3CC=7wSQ!H6}1;mQ~v?JmrvnUs2XjnVVr zmDWLy&hg#d2IG)@5D&R`sT+DVlz^dVPn%$s6csIgP>LS0yDgfx_tDd{$hoFc^~1^0 z)*3((lu}tYx_RGInWK(BI&A$39Q*a5h2ACc>`hkjfVA7@3h-3b1cjh>;t6Gx#OQ9{ zk~IPbYv-oINV}z7-aAKqlO;Pn>!Fr0*{D?vY23qOV#eg6jVAsmWoJTuohhRAV}@a9 zU=Y(Fv*Iqb^y7o;w}A5w;CQgdP*Ph{e=3guA}xZkLg-dF5w6HBFBkyDF1q8mZMQ%8E;CP|c)%ELBV@zNev<1&zq zG+g+$-2;4)#Nzh^{$-(8>FJX_5OG@5wfW=Q%N&L^5IF_dz!hl8o@HO}@-a%Mu``QB z=4K)b*QJd+j-&z{y3DVWY0L1{c`3GUt=7QArj*A=T@{f5Ey8Jrvt{hNX{U~olHcou z1d4T!#4Px!5U>0B8Am+p?yVLQ;mHKsJ(+n7BLDk~o^1W%tY--$#7dOU3+9q1r-rlS2qi2*g zHbys=9uw=+_v*_iDr*2bl9E+AheHnqRNDSz_wiQXi2WK%uRPx2sH3@%UdF@SPEln{ z%M8ANAYj*hSNHZ{1Hgr(*-@%85$^56k!^a4w25rz3K`3TwbY_l_T8lut?>nLk@ilh zncy7~-?`8h9+GfK)zH*riF+XZb)98eL(=RDpmr8OlL8Rcab-(W%1k)T63g@7 zkPi$ZhlYh=uzaUa&&(fDJ7V|o7_QifWaPzk^k9q^0Oi~U2+|Sg+flF`B8J2L$9uLa zn(Qb5y>25LprgAe_z0kLnOjRC)9YSNsDb<(zX1!LQG(uW50Q2c>yZoBA-D%3r##vX zIWgq($hIIujs_}Ccl01*I-URe@?f-$18~Hn@Iu5FXGrljVV%59y?WNX%;r{H>&+q# z$uvruzXRFw3E+0ma_w6|35V}5Rie)g7xZ~*!f^1V)-oDJz;-*9>9hbz$d^SxEGS0% z!3k|b_uHm|3Tkj8UT%t9MwI-1O=>sS=fl*cGFA<4}$tO3gqxs1MzR zRQSHrSciH>)ThIc?S2BAehrh;&MVW_?xTyQbQljoHm3ocl@f*OzTioPBNuQ`=rV$5 zo^WFP6pD-?&o(mTfe=g7W)9^cN*T7Bx-b-POjvo$>3OQ!>Tn*M2>tBrY((wIJ+sz^ z;}nbmrsW-344(G)w|5hWBZK@t*{IZY*mtM&bvrw|*s2|{kt#P+lLJ`wu1o2~czg1eo9=HO!arpA^fa|Qd%W+^=h+$Dthu*lcc<)^W2f|Qwm>s9- zmNVy~CyFPftFqkcgsK$NrO$uNw;xhPc?r!$!P@eh!f{b zxvS)pjGT2Lz>u}M zo;iQM;Ok>M8~>wHI3jY_ZXVVuo4}9Fi#Zp;|4LbMDfV}SplG;>CibR+xnCQf%-%F59^VDpYE z1cvo5yM>pM%{Cu61zcoqoa+cL3eNx~NSnJy>;SPU+m+A3lITJza7fHbi{w{cAANol zRUV9t;>;Ri(ojCv;IYRVNMm`T18nEU8EKZ!f~!d>z~#~(e{b^(2xvP~(EU2vx~O{G zKcLT(Q9r&j%n@PpeYzYj?P%KZlf07CcXOT10#xS^#MT7n!slv z-dgwYOJle`+NMmed}X;9NS`0xb&1)X#zTm;z!{+d1<1KooY3N^PjLP5{q7|+*AC{nTQ-K(kYv5h4_=XxCLCYGC!zew6nw3W%T{rOx9q;ICStvTz7%_A&U_zgKWLy@31) zAky*bBS07ft-@LD0b4INU$z9547Q*xP|!)IcR+#RFcy)8nfWW=F0$@r+F!l;U;(nE zo;YJ29kT|B(M-RpvpYE?T4NwNLHuK~et|U*6EdFe6#KENPHF}4@IdpiRV5}TuaZ-t zsuI=zR8~gulS^{Gj_RU%v(8;ISi%dGfXBE`RVRC`EF5oz%P!xKQ1>SyM3nT+@)|CTxkx=aYI|k_CK`6YOBVN08>$P6?EjiMt3j+X~ z9vmGSWvIXZ1u*g27bh5qw`P}e&3hUccn;+~)6++xrVwI^$VCz0@SAxlCwkOkFTXQ>9#ZmaKt_g`0-`Waz7}~aP<-$p+{U*oO)W9H5}zSd zY2d}2W$Q}-2VgjKca~k~)`Jk!eCwW(xmn4Hw;ztU^oI;0ZUV3a6zrfzS5bC6=wkC7 zVE((rt2mrvuqrBICYAm|=CO{?*#T=I%Y>YNm2I>?qig6G9%2fGmU_s6>ajerrmtVW zz6N&T0?4{MPUPm_5BwZal}=F@hy@Qmh2i^m(%3+KA<_*gmwQ60C*6YLSS+Jrv9>H4E3 z03Ly$f)wXyIVa05K2*h8?>09kw6ZiLgSJBjhl@hH3?N9@US9q}aQC|< zaHe92yg);Bp6o^Zl_;D!;78!DfjJ|HMI`Z|7<%GKFuDpVq|Ta6>~(kYYq_*?aN@B6c2Bh^-ja#Y$3VCj0t;rG5eI^fc7`JNx_l zTMLuXWHR4Oy?T697}59GYu7{ppRLlJ5SGsd$t=I^1h7^*1V7Oo9XF0w-|5Or>oo}w zGV98-1j-V)-vFVqO)#iw2wX1^(~c_GR}~bZ8)7A+YXyqy?m&*yo~^oi4Mf&|Pi;~X z)L8C`fXu@0dIgosvYL1T{1otAu~nj`AYg!^xRjeyX}p?8xs~UccyMQaXL-S|f(Ye( z3;-9iftqmEEc1!(LgYu;yf_i}T<|`iepyn3N`AC_P4{;B5G&J;URrt31XUMN9VIj- z<_&3NuaQQkJNOXNR*^=7W{IgMGOjVZAr=p zK(;arVZheuzQsO&p%GdE8Ccr{fIj?uenZFBSac(JdTAGuf>4f5E<^m*sT;PXpK@eT zX|yS4CL@dvKT+5cU;D>sH|lJ z(b9T#$bPkM+(-qL!c}l{UBDldG3WHd!^71Bv^KC|Ez%KN$%GOFuE-$`@c$}$@>-}h zUMe!Pv`pI}J^vSQE#JYOoGg+ow)+bmx_IvmoK%Wc8=u_OACpTghh)k$*GNM#x*cRj z+6D$kMv_E%OKi%%cucF9c=Q1|^$O~SXlAZ3<~*ly9MmbVflvm&@0bQ$a7YsB)}mE( zLepzG*Pj)DeA$u^nGIzB4@MMtxkYfWe3ugk=gT%+F>&vMYK=~oK^iCnvYQ?$2K@BE zZ@t&&{(;dhdOc+Q{L43QuEXKyJ9X+hkp7oTV@~?g>(Bo?9M;Kbo7GDP55p^elr22^ z{VV*TL*`?wt0FIN`Nv#a{*`+j{gTZBH4>LX#KS=41Ioj&>xm#K$sl}T|9?DU`u})( zWi-A1U;@aSGTff(Q2RPY1?_O(uo= z-e%%wq2gFX>;9`(`_v82aHf0(OUuVzd9Urv5(iJN_6HKS&f7%zw~idxxN@a;`I%` z6mh9pFE7>9x}`ysvK8HdaWyS1)TL=^b~RIL&&<)Z<9qeVUU#ZMTBr zrX4p3P8Onno`i2#^k68lvuj+vYCDAKRJ$X&WQ$jJ(@*%V=S$Zw)){G0@Ov6ku9y5$Knl6`Uo3bV%?J%wvShs$dwjy?AMsX+cUHCo#ZsgM$MH0?Vm^gdc#lT z=1y$Y@=26`^vk}qaK$=SJzZPxVlpVCyga+jPF&&7m1y~TvP4r#!ETHgik;}^Wpk96 zj(&9QZIAAZ2nwsb8?BwBqi?%hzt;cSamupT| z(RjQ=!ysl>Xm#IP|H(8iH)=^yqe84mPfa~FJK0zM-X4yFnEl{24gfJqs;W8w{2_Q| zu<`4xzgdAs?|sq)>`R}*jooT3N9<`osVA#pohA0*aKKtwjlI)|vW1gB(vz#{E9VoVOzJYbMtbRB?{@h zcp!faRC?ads#3FBNS8wt1>}gxB-XeVlKy!@#~WFQsB8@Rpra(9WJj~Bj<&$#*X0|b zfnM61V;0a>xu0Fn!)vm)eU56JbnMDn=M|OcM3W%W(uFhiU1sLW!^8KB1ux3WY3$v* z`wKJT^s!^Fa8EGUyLay*uMVDKqO}=4jYY+o{OPz6 zE9Rvd5DS^0NUw-QBA=(@NGtD?j8zHj`f%ny z`9%`CuJj#u`}KCmGi?<_^h<~w{PAApIvpTQ5h6GiK{B^v6dzr4Sw>CadUx)X?kCY-bjbK z-_x|R@;EN!5FxeGuBg&FzTv2>#z8&E7*DFkq(UWqw2U!<|8Y-gEwThI|4coTdLFB$ zznP<@u2{D@MLj(Y##}?6S4;8;fspU&lL(Zb_BP_KdwKZ}pNiqzExO;{JtXb^7P1z8 zU9M&mw04(Tw5cx%`mp2Yu8$L<&+BsL5i)obFaJ;vY5ivm@>h^xJ;xwb5ipgl01l5~ zk%k9|rF6@dWK?F-gg}Ymm6X&+@;g9`vjyc$6_C>zR?I@A76X9u*^3*!%DeLRIL>lu zzQ`EF5Ptr78?A#5epzCkw|OofkX|ZtMrE;o@6hC?_+1z=q0Y|EsEaOk6py(ofEO*m z-Xbw45S@NhujK`E0^(&wU1z>q4LJ=r%Yg)L$Jp!D;}M{*9@2)>F_zP~=HKSJ|A}AX zYi>)16CW{(#V7G3!ga!xt!{NlrE1e#;uoagf3WGZ`0DV^&W(~aTna(L#H<4Sj5bD` zP`e*LYAvr}l`a--*(c2_-MirSfxqrwFUQqQ(la59W+)QZ9>8?Q3Oj0Ju+0sz4SffU zj{oNh?d3d6&6T^QZpw$xD{b0w!3~|LC!3+rWmRsr37c0Sw)K~zdr(ZlAGXGaX8I*I zYw*wC1*6kQQPps87^iOjpa1>4CrXU9fPw`d-M{p1ygs@I{zn6fYtP}daQ%95On8L; z`(Pe6EQf^RmMsm9Sir0kT)Zdcm9WPYtuVjh#%A+!=1Dd-7oCuPyA5W@WlqF*m?v@SPMGk+Jmo?Mc{M_D#UPJVrT-yI>Wc$9DQG0vKO=!f zVr>fiahCggkN&E^iY+YXz%+HWV}2>jcQv$U&EENxPtp~{Ge*e`-ComXynHo>`V9|2 zxbXege|y{Jn*|I_ha}c;Ri}w#;fxa;)zuw2O#oD1uveS+|9x3_=<$I;SKScc){{~s z&-0|c;>pMA1Z&~=#=8dUq`rv!a`%6?L1=G}_QO<1h4z~UYxTV(z!WiEv_otDdqwbg zoK&eKSy#|ppSGD8wf)cUx+og-35`3M&O@uFvLKmj(RehYv>~+p99o;-(>HLDw4k*Y z7qG3pZMNvV(&w3*&#H|%!K44-t@Fgpr00*rzv9;yt8v`qMgnb&r(i*n2Nb(_+V$Xl zyN@L$YlO6J{=HOzbK|*uKJICXBj^Op$S8J|==lXcmvFF{@5mE{<8gHf=V#y!+@sdm(DB{r%Qcf{k7ImfvU#rEZcfUG3NPqmk47J^CY^ zvQUdGi|xyZ5MjfzN@dr}2e+R3Gnge-qJ(2k*SDo?m)s#X>~@;YbZ&y;x$yeUC!Ftj5}~jVV2%yl_xPVn7Zx;>{QB({oIG8ND6W3mG}3P zL&^b=Xm43)=^9iC?w|bcJ3Y%*WHnriCP@VSe8pHvH|HI)8Xtclsa@9o*?dF|jius1 z@_SYF5H?fC@4I3c?}dIJ_UKZ#q;;DvQ@hM`g(xV2T*QxF#YYRiS2atf0A#fL zy;kQ{H()I-HBDHwz_t@oTNhBvj?!d zJDxD@^eJCUT=S`4V(_0=U0}l$8tAy`zroR7Cy*TXVjgb)ZJL#7T-cnh`jRw7KEVPX z9%CO2!*l(Ap1$RFa9HT^@bHE9F0(!mb^YhPjm`p!N@{cg#hRL&d@YO^_Rpf1T;x4b z5s6PqayfsTSEQwD&wqb;JwG1|-nGbvt_(nZ>#^S^fW4cY=Sle0uty!`R^=~9icCg5p0c7-K)LIQVM^{;$@E2e`D%gftfASm3t zHMGQ%FrC89E8bbebNSCXAx7{3_7M_lFcdjcpUY+0cyx5wdU-u|_ii2A|Cs@dN1WPj zrmRLk4LC9{e>pE73xUpS=!5y6IbB*0))rT#dt(}ye5at7KdB8JGYl+Ah0r_AgA%a$xcl{5@Ak`e14K0 z{CB&TfzI<40%4D=#^&GWNp)>o*GaMN-sh3SD|>1G9g9I&IJVenLo!9hYQCDFs*Pdg zoLCY2y!ZF^$;sYt^Vy%uRUPVw;qgiQ=jsZ14Sgo-Bhh?t>|G(uRec$XsoFG0s+@gR`cwc!r)R1Rcvi2cN^XYqstOi#H ziBS30OG+d(Ki|xoEA8n56vUxq#XbAI=?aL@r{_z~_ujjLKrPa5(tF+1666D#fow^G z`UsbI#s>5#QNH`ziOnrZOh!_vTTgjhxLk1S_P73N{CrnUP$}ifX8x6E;8zjB0{Wh! z_tT4W;~9u~f%*@wq6|c!>QH-L1S(Ui$ob~YR{*%Ru3vu%7H)l^fN3Ol2;rD_|kLwI5pa;P;ef<$xtVD0y|NGFc@W5h_T9-8Kih^-fe= zfr4oa?kF1ZOQG7@Pz{Z>_5`;lPKN6WED8Ain~8tFgJKFsPaCTZRIOU&L6TwdxqY>( zNM^nMOMM!K4#`wrzSv%WrJd^*78Zz!KonXrG?x`$K?_NDoUVaE=GU)Z{ere_+qN$E z5_-`k+x5bVRy{bsC$^lJfg$A+C!nqt*qbxhmOrQ}b z1u~e8V>Jm&OYTXuq(?+lKp3BmvMw4vX^8ZVoWP0fSt|~|7bIod9iR&V2<#@PBH{=uwSwiDsJ{(ZI@YEdG#eIKlg27m2RRKx z65;7njmwumPXxSt$({DJv%n@Er17V*i%^b4>>M&aq5}|-ykcV7P;+Ym1%3)DBu^tY z52?&*K7RZJ`rka;FPBk~8VQ`POu4MS#iwouapa)Q*C%3=t6yvixO>|IZu( zSA`NA5J?1~t_^ewkW5ka!j}gxNgU{@yDq;^g;PtZ|FWRCne|uRz^103QTr=zz%c}& z+GA+@e#06ZJG$egI@p?l0Bor(T~RsI&OWJGgEw8NKNua8YiEzI8o0M@Ft}zO<&cp3 z#R|A6Q~$6~O%E9op?>4<=#5RmYWPg@;hKkBu$H?fB>4o@>*qbQCTIkIU_n72V{{w(3PKCLqHv zb@SZ|n44-%NFe;d7k|0;sN#~6&vgnpNh$hl*pLt%Z#Mtu z66=3cg7t=H`H2|(XDH0hx^;C${pQc$u8Tk*eC`p1gqP_Y*U`RR?!ljFT-FW=i>%&D zi~1dYqKg08gq`zGth)wY@AqG>y-gM171vvwDE>rxRnm=90DcWjZk)E|zt>EgrMYH=(j}fV0lU)Je_@XY1>yWjQf_bxNY796+&0&z*42maIP!+i*XN9&ni6))vPY2k zBd-5lPxe`B`C;QIUD8=8R#MW9$sw7z%uHHg$c!#>NL+C9PymkgU8=skl z{5pzP9xFY4+MtU86q=&i{$}*LHiVPCnBPf*ck&q#O=0=P8O_NC$i#8Em0MqFXs&}~ zH%@?pKM>`IF-cYQN%*x%D*CE{6z11$lMVW%T~Rkh@4X3oKIk4C7JG(bag1`u@*(zB zs$y4x8Xs>q-)ojlzs^I{IY~tvMX4*dYI+*PN!-Ut?DA=8_bmy^Y7qtk8Qi@4WHsRQ zCUI}s@eA3C5}bmybWVCq2Mi6740vb9`70_mBqPLw!5{W(Kj(wF@GCq2lsC|?I8IZB zlz~z}x^a_2!QRh6PR&zIZco&fU(4E(Z=b@}vHkF%3(DG=^9f=W`dho-jHk09d%8msHjl#RkcsBn>>_v@`mJ!!a5g|SaOm7s-OS% z36G{bzrZM3B#9+kU@3Cv^1&MY{cDss+@6PM%4*zvrDb4H%+1?jPn?>v!2AB$jc5)G zB1v+kz}4oTOUvP|PHPvDJe#u)BbmBYPHvCvo+DTPbvqVC{%PBZn2b`=E3UMX?ybZz z@a+u+ked)Wxa1U+30>Di9-IvNWkL?)ohDLSNLIl!G@Itha;ir|vWm{AUYlQH$Bi^? z{vlsXF*i165CfHor@we9QA2U@_I#^gOq+3do?Q}!bO?jh(I$qH2u)2IN=gG|Gc%tH z@>eZS{a^3vot$BYlFimqH>PZ}9~&^^Z5aEc?`nTI{Z8%o zO|?6SW9$kVAMA-%JG%14^+Lj84PmxUa?Q?tSXzOy0w+g_+w0Ch!-9l1S!JN?%ZH@N z@`6F0HdBA?5G~!{e#onexp-f&VeJPMIV6=#D%VfOS zRhgUN+AeDYon>D{A+3%*4y5P$~khW2Q0H1~G{3(8Z zWe~&bMxP_Jg2IBsW2Q-oBr(5QkNdnGsC_7h^f<@2ebu2*YLDg_&RSZXV+`fTN+vy6 z-$#YSbhfs~AH(}&T>&tu(u<0&6NT@3H{E>7RjoEkB~ww4Gc@7z^YAEux);?bKzuwk zJx%di@dVX2hkcz;(frqEfT7KKxqU95Mefh`s%ju(+xlj$VThD_pg4~8DU`e6iGKhVlk3Wvu^QY zaiIBjgdU@8rMIid z?N+pw3y0U^n^ZPTCu6N%$jr~G=I3Xcx)#bvp94j-T{1Tqid9SU9DclOLq4}07M_VL2QCwE-%=Bhp>Sn|W?D130hF9^C_NuUbigC=ljKe(P> z-uDDOoTY7`{|0T5kC<2^j^)>F?*Rrt#Sy%aj7A(KbYkTwIga(?)f21xkqH|wI^VyO@2m~UmnQ1CpClhzx+?v=`yH0V%M zz(oOR=(qt3H6uJBT0}u&!E|(bc zyU4rHqoX1p0qu~%H~dU6k<;@{*XKi9i5i;xkL?(IlKRjo$%MtkK781`%!DT5<(`yu zPS0Ry#a*nDi77t$c_SNk2MBpmUpz*tKRE@3QP7e_U??_yM&fhlFQ6BAi;TtR;tH}D zFTQ*rJUoa>b? zY>P*Bb|W_U&_qvWV!tU=V=YHo9|3b>Mu{Eo;Xi@7SDyU!)~KbOrR9Wi!Mc+?LzvR( zz_4h3!+Wpnsq0739KL*bjECVue|{`{b|Wk%#^5K2Sz1sr?9-VY=p|9d)qJ2zVfp>pzNwz#Z8p;laT#^YTOu9ow|K%DudNhsJw~bs_mnfHCaj&@f)8sFIx*pOW z4BC~1UMk2bm%j|nEJ4e$C3Go{G$%dlb)QzbZEija#W;0q>ntM>D7rS5j^)%MQSd0z za6w}|#eKn!`Hgc032~2l%&3>QBi7^P)3?JXs|Sb^g*>hNcjZ+BTMKMHRH?+A@|b^R z&W3eTjv~uyTwz13kTynYr9_NATqhD+c=QvC59V<$w*59Iqy!SIW>iX540YRYQP$&~SD zf|;IWUIBrubNhduKN}~uwePJOXUp~$n2N4jWs zc2HIKNzBW8w-F2U&qmiWSEeJ>Cx?BuxtGg}e1i1SN~WM8N{k;q%`%s&Xr~CCn{;*v z`aAe|dDY$BOHkWd`D6%@_$9ari~L6pcDgSXdHVx6Ke!_X+L4N^*^ z%{9DR!KP%d>FV--7Pc+Y{rnP%dmzPfQtT;y407(RromxRe54y+Xoc{1n>Y3SNV7NZYp^*9(vfsjqISO~rtx_L}Jc6M1 z-yZc4VF>h*hu6;p(x*$HTUu{bA+JO8fIHN!qS@^^K7J$6T=xV_oCNg?Qu=y4f}t*Niq1~i7jZrgwQ zI;aic>q@#qF3?3AfqLlA!O)7U1s^Lp1cK}5r}bV#LyRoQ zzCM9F&?Hd+tXTT_rPpoVDLHnQVx;1b%2Kj(?I4d0S*2m9QVfU3iCkWRz49C6O#!N$Qx)_&_lomUOKlyLKzmpnDVWr zrUr`gIzR8CK1i1#4&**R!2+dZW}XHi1s2@|jS63>QOI~jcY`44)~)F@#jRHXYfzvn z>!x4(JTBh6oNO)sSPSSEj&HTpFTOhy$7xC11ZS2F-mL+&br7WJdzSB+4uWR!kh@*f zUgWb*C_n#S5#NddWZ?|-!*8&?jox%fte6`?d3 zPdUDj)RpMpW%+T}WaEY-5D*_1+7Ib0m6_SDLb(7X$P1?jg9R#fzfV`3dUReE=-ggi zC`T&Vl%D}8?wlLVvPtM)PVaA7+57tr`gM@|}l(m>~Ja)d@@` zal?n9GTE;U#X(UHr4C2zYR?pWhZT2ZbU~9zI((2q>?1#7puPFFZX@x;Q1Ip(4pu&J zm0R9oDhEW#td^)TAMs+EJkcHe8;pIn^&Mqq%~38zyN`%?MPDlz$5$P(_csFR2^V=4 zf}{iKgQMj51Oi(q<#xH0dj&{pBGf4^ZN&EKU?5ZXSSf`L*j~~W+k>*iD@f%BuBOmq z9GyMrGV>|nMyr|n+DeYtnn$BS>DhCGu%eR6Y#D`675E5A1=I?pH>s^FEwm?*dpfdz zKMH!`DA1o&Kt_GO&|tq~$`3v1r5~j5M4O;Eo6~moB!Xs^V$M1#wXz)*70pWv2nVTO z6xsGl4)6(RcP#)30G;nbVb(?nQsEN?uRtm#b(Md0hxFL6B;3O+?%^Kiq!!gsVG(Az z4v-xLe>Q9{l(FJGIF=d4yn<>)KQ2$UuUpRo(}M;T0o05GVnl5S_m1%0ic~Oo)Ve3( zIFbnBZOK+;B>dUcFxH+Qwh#Zp0Z26?4{g^W4uI)NoL)1Nxc6I$zsiv*m$VS_fm+1u z@V5n5Z{BE@YXZ}PNZ4FmuCy5{s0SiH=rfWCP z|AZQQ^N%mh#1m!4{IC^0U4z3d<5SZ~Il_&;P!`KjA-y7h$R{*!;yi15IA9~NU4BoW zHi6nW6PjBtv?`-)#(Mr>vYdb_PHP- zEGnd9zX(veFuj04#$uTB{Oifu~T-#?>8m2s1>X||1yLNye7zf zuavGT#Vl^zd_|!E&$6fT#k>5+ojRnPuC8u&8KcF%mCQ(RgNpFIo*%PNWXB(ye#oKa zf@pF#x@U5Va_#d%ukl#tL}Sm8!TRYX*Gs$39_MB87wUrV*nrIVO6h79FLgT+Y|?2w zb4|zlgXK5<*`Q_cI$h|>BgF$or$FUuPTek#*9SjgN`Xih*w`>9-E9Rzr%xHsl9ZGX z^gK|arWgpOM8?hF_y8lXK%5D@I5X)@6YG^d;R|S z{hmLb*Xwy|?z!i_uj{mD+u^L<&fBxJ1vAR^*MYkek=4JO652tFu{(E_KQ& z`%T@Fz5Vh}V;gq)nTfrG|B+z;_=aN}akt{OuKhRWViGAuqn!B4$TS?Ck zTA3{nQCh)D0wpU8)L1u|MkP3sE^pxBujp**7JblI6Tq_DZ`CsYXz0OLLmkfS+kWPNjc!k~?dVuqc|EYoxb)>|_gL`;2roq?!cs;^1GQ0W5sS0) zK9fy7vy?}x)MV&s*7SMvxI1_{Gi3H#t*o1K_{(R(gL+KM;aLFTSXFmdXnXH;+fM&E zz7S|cTY#rTBTSnfjDdDncUaOnVUNCsJ!iMAh#nfEZI1@J(UfaG!fPwkmls_V zMAg^%`F!RGydpYT29uQ{_C*IA_as}N;^FYcrvibLoFC3VaJBU%P5*EEa6E(8{sg{g zrlW&-0lHXbes^wZcjQ9PKWSXNI-@$+>v zRhNzhT}vNwenSD>Kle|*RycnM>;uD{la92xA3Pd34q`KtwdjM9^@R*ZS6iF-;Z9f2 zsRzQTwml4_Zu<|8E~MT%&QuF1KK}7S@vZ5C^9-5tV4PZ6TGqB7fu?rAT9s^70|S{^ zwH}PAU?Wb~-}xab=h@v)9|4ziF3eqk_M$2a&Fr}@JPi~CA%bx6@w=HJ7=k20p5Tz6 zJj6fwr1z&?08(j9(SLcFYATPw9Q zNHo7(U~muy=g~zZE!#@mqSf`w{A2A_)!IaW%y+!StLE&S%}OSMo?8t>-+6vH5BtV0 zbkMnECn7g;Sb|4`RE7#8TWpx`u*rB2*nkBc^&xl>Bsi|XcLp8?fX#5Trqa1_Nl0OO z+E$MJx_Eq`uK1Szp>Ru!y}CuL&RtD$^WAtmm6-?LKC=6{@FeAA*OU#HC}q>3MzJeA znTaCVS{KXC2VPiTM((=5pL-Way2FVUuYB-?WnjYcv2HRLm?D004uF-N8o%N4Yf;X&o=iTUJ_yENm^XXCVD)TV7fx~qNjDCCJyV>GD+C!XM z2T;GinR_=_ay<<;D>q1WU1UU_;fM5uuq23c#&6>oNidxlFTQo>#!@gcE9kN7*+&DA ze4yvqncPyK()R4{>5ed!q@SPe(IbO$X2*l0d6b4)9rE^tos;X?zzZEQG`~9EJ@ZkF z$q<%E)jh0@0xB6aKKnbG-KVxRHloKqP>;<+zp7vPo-5@&moo{}SJ(mtYHYNYi;eXE z^JIvlL3Y_tn~f8ia59R$GhfPNtf3PolMe)w{gFD@(hWtlU zJ5-6@XBn?J-PvtFmB6Gh%1m4@W zoigpm5_$40V;=|c4}pAO6+zXKqav0v##l7h_a>b+FMN*PR7Ky<@*4a4P+RoE(`NYt z)+s0Y+8<6iWc05U60fDQ{J8!kN*6ajcMK`q*!mNL6L#t%9_JM&%MQ)9Ym$095xIf_ z0^-3PK01k@F57W`sglx0&}8g>zS1ltC>SwZ4V{FRDa-57+YF~HTD;yVLTz{8+U;xQ zbQ;8&mooJ#%TnTOXKitDu~y+duG;+W-IEgi>QE#H```^Gi#Vj=*_w`ve+T9vlR`?WVFCXJs$5nEJ{Tys!qYilRT+(Y{pP2JB^|4b01UIf`%!zF+#%4s<{t&z z;R59wAuf3l@w_hs&d7GQ?2O5#chJ3JW`gwS zU5XA~RzLcts%^meVhel)$%_u1nPwt(@aWMjkd-!MPHGrFRQkKxJH(~aMmm6bSd zGyPO-Hs|GWi*X@7W&Akc7hIAW^6H6?3ZO-if2lpUqNz(Z8!(KV0}SU z7}EwU$a^kLHDStvxB%kbV4q*+{)=@c-@r$LaA+i2^>J#1t4_11!C_P+kD!hCes+#y zweC>!-E@6nE(=PZf)jm`KPQ|!)PtkQkzN= zM5$t(6g;6u+PF2?i<^6RqbDxsHwfgU(D*l838S7T7Op*Af9~f&+Dpk0Cb$2&@e!lU z)$>@Q;xaQsi@sKi(>gMA10|C7lhuKa=JC-c;;q`S=e*Fr_4(8h6inCIb8IG#TrN@s z!wojS_d4DS9P{H_c$ZzFF!AS=zMB8vlLRY={QWvi+li+(7ONYKAaLd!NsgYZoSSFO zwyUG9uo-)G?}2HW&FF~?tgZe-6A#VrnBV0;!NBDd%M_q0&T>5UMaB=NL$O%-msaP) z|95~$SRt&HeW9sKe#aAQv+5s!aB=cJU{}mdOI&6K(WtD1GIPEeKzIlV3Hn(_NQSU7 zf@@YXc3&q1V9C*>??u@qGuvl2Dml}4RQ;YWCRAM0E_Sn-k04j|)%!t%L*9c^)*vpZD4ZPL4stZWHE4VkXqnkhu{e z%zA))*REe*Y!kIaDR~+_ACuKvf%9o(2hmLc8srW05v>Dx{STV9 zAoUpxIXo%UqM$I_yBDzRSQ;%+R!np%iE12o+vI3+uY7efvrcL#fo>oNfE;v&5TOH; zJX-^2b35gV8wx$m9TV(JKMN-r3Bz^~*FcXG3=zZQItrnG?XVrNW<(TXeA0Qt`HI-b zRB~0dWsuHHj6|y)xOBY9W!xbxnYR9LX+ao(Lcao=1NV>r`A~doYCNPZM_(5f(#V7*0#M20byhYmfSVTHx`_`1tJ{So4hCr)0Wo1TSN@Gne&FGqK(AT~>L=X5 z@j%^ITuB1~XM{W(A@Xcq$k=?(?|S%NXwELP%#~VNw!}LrEG(?QzGU30?ey=<#Mrx@ zCv0=)|5w|qLww1hYnDxWA=i?w8zVbUD!oEg@6w*m(Q1kdoUHIUaS5?)+-6R9g zWHS?148&Ob`*nK{bJ%Wl=?oD|bQo`vejB+7Vw)^_hfBVlV!;0nA`1vGC7OGi8-I^B z5PF%x&Ui!!fbZ{!>2gkn?pXhO8HF7l9v;S`pII+1&LOp6U8O#J{$|O^rQv_n`}RWH zj-o3&Op4H~43s|agO~!?P{D0OeBq8X76YqKy=(=&P#lFQ29R zqf{V-^z#jzXO!>ZfhL}F4D8E$x$CM_?(nI(;6ZBPmypMROPi8T!+`_$n{qre(DU6z z_^6V){buwU&^LO2>w>Qlg)XgIz<;ePot0fzC#3%Ncm}+Oc^3QR&%MN2EFVU8mKGL? zuo=c^2eeUSwkJd;Z5Y*`%$50{Nixwi4Z9e*cdw4Np1KvJ1sOykks7gUMPI$>pOpXj z$^!>p&M72V5j;8}ymp-WU3Zs9v6h^3b8imzfGLv@yR|^YVZwz18Wzd#o&lpbd*RQ@ zKDKtQC!P$c0~Xk6;nxQW1M5U@^pxFYD&5@N`UXmw6GLThlUNc8?@6CZz70bMV$ps% zG#iEu5wrN<9Z`pxGo@(TUzH6qAk<<#Wv~Mf5Fsl_`CW+2$CPPuWxg$A0)`dX(~bx2 zo40DVz(Fo3%n6?!a2~)H^T~eh-@z+nXo9BPl>4sWm?%?f1OAjXc&rP5dsgQnbzi&# zADAE2Wjb$yd7&`@(1D-Npl zw0B}IjAO__QRh$AuP9}SJ+F{OV}j|EG|*L~4Jo~X!t9Rc|D@X&nEip$O}& z&$!d(Mvm|HCC=$wPN?WNJMzWlufsLmfAq$^eOqqW!WXNH0ImZcde9-p6EL-8JL&HZ z+QBDGB{lNz5p*S#5lTM2x>`^C43RexCraq*kEL2sVQI6~+`R6UQ;Sa*MtQhK5xc@k z?lnBTykO1Kfrc!xIr{HD9za=NV}Z6^oE)czVQgV?PF;uho|2aIJQd!*Vxn?**2Jd= zZhotBNItH;X5}(`JIQf~@%oXEMuff&|DMOLx?U_y41sxD-5}|E!M!d5_JId)`E;g` zz886fNTx=MNI}ad`ABAO!jL>jojlo^!9SG+0Z6_Hm=fQEZX8oG*DMq0TI#1vf3HduQE>MdH zJ~Bj3wW9dX-WD9c#4OKh=gy_{Nw;^H7@|ydLU+;SAI=X>cv&XLAs@h=)z^n3&!hM7 zpxr{tYkKp2JQ8h7K$#?}!z}MO%#kfa!1tLRO(pVSlCNyc(;+|bvID4mS!+^qcZOa3 zJ?0}M)SAP|+tn)3T>18zqSB(T5Bc{5Q_2j<0wJcRz?m%b^733k@+Z3!oQiWkF6270 z1w-Fk$`ZEd`~nrkPn<(xnz?okL#)$ z2`T}j>5sQ+fqEqm#ZD>>rFzws3bKhs*8=GH(>5{wt?nlbq}&6oNrn?&_NIazbT>q> z%Q@h7vo#cshUaxzpONDS?xV;X8se+oH-kF-VPd>*LfQY$ty?V9 z6?sO)$%@$g!RdvhKV%4m-e~t)a&jo_T81kJnp8nt1s*a}eMwm)1uhMJeN*BgvHJB4 zT}EYpeb~-#jag%iV=RE}OiWV=xmg7Ovcb8ZCp#X8ETE@-44B)6aq-f%{1cwg zF6HG}-ym5NxZ;GfWI>#7dG!P|ygnE(2^y6yC}6XUr$QS9LfPtfIltGgL5V~3TY@+# z6R_rOb{#oO5E}B!XD;4?R5lh*@M@3i?J}S@Y**+N(yvIgHf~lzyfsZz(T#|>oUO-B2ZbRBJH2+`B zY$hME{C3>UeosPxp+)^e_-Bj>vKZ#(Q6H|(l%YdT_u~1^C|Xy2y&u^2Qr8bLem^5c4cFeY zBVyCScjO7)xG|Sh1}P~i#MBS=q#$7}2308y;g?JtBBaH84<4NN=^&Svjcjf%Csbfd zs#oXdw%9dmI?4IXyrS1s;>&jkBGJDjj*?MMBrm~IOD2OUVq#*Hm6a{fiXrRp(4QY4 zWU=p%8*ePyMJ{_eVga?$j*~vWv}bC-d4?&FK_F@qv8yeqQ}%%RuR#2F|390J)wQ2L zTk2(Q6(2(>m`aP)Zxg@Vy?Gg)Vs?N&e(X8o=YxL|tRE|uVJ6NJW#(XaMAl=slI#D4 z(f>4}311(wD}WEDzg=H{0ik6qjtZ3EKbo~ca-BQFl=h>D%^31<>+&>aE+L)x@fHDt z$p3xb+Qyi8(b#mBJ$U}h{L@U5VO(sqHPL0LY7DpZObV*pUMYli-ITk>J4xlk>CtE|b|pWt|5 zk!C&QBaFz?HzW7iP=|TYo%Ar!Bm*~`kZui86J*0bKN+u+)G+UWb)XZB*h9>oLr!{; z3sP+rUY`*FKzQC~(QIv6s%->ycs5n&E(Db@ED_y-jdBVHEnLs0-9jsgpsrh@9KQm~ z{Z10He>j$91N%)A`pvXx^Cu4B&s_LYpup)q`KD7W>oL`KUo0T&BXEY ze^~&zmyn$hI5L(SLY~QHH8p+3z5%6ShWR{`;k%uYK*t_7Z%0g0MIrJF2b>#B#XOuY zUb%9>k#vWvcdBtH%QS6{=0PV{?&lp&l)!1$T=lWKyXV=g^-EkP{L@m&k21xN0yuZ^ z`bwwXLH2*gvG_k~K<<7KIyBVA#t9A=_X$DY3TpNjJ*-CU)$(+vT@Hu;Vd^pP;9LKD z$O>?T$|Nm)-O=vH%75=J83xWr5(VPt9X+*lvgj}6xJ6sJL`ENE z=L~esX&U(a)mlDV6_S^o_D*&@aq+l$*^D=@lN}@W5nGUd>*0SVvmah4=&t5qL&C~@ z=|Fem+FDbR#O;1l6+JQ-WzBM(e~bw;{J}@vFD)YaU(AV;KX?GpCcAzbWjWBlbTFs$ z(W~81G#-AES^M|V;H9&XW1s!#0ExeES9$Vn zm0XK=2pFUjPL(G{n*^ro#Diz3`;UCmZVlQ?OTJmgOFuwezQ)YRXPkTGOL7=Z`xtB2 zy0{5EB{UHN%M;^SKNjrc@}U(C{znt>`0~wuyQ16jNye|L{YY!m*Zxb?l%wq$hzDFG z7Fom(8WwClZe;FJ{>vXLC#81jq&9+X0g(Dx57qtqxz&BK(9@pwWE{&^MBgawUd6)P ztyq8RI1``0c;N!eZx5$8I5{k5eM-jD3p7b^b-QUf_HOVAUe3NtkZEt5&O=8v6$VXw zp_H=K#6pnVb`R_4N@3zW15`sD!+NThljmq~k13tTpW+%mt9+1yW30>4&7OV%N&$Db zn+OX{j|`Pnqp}{e_VDp3MPs78V{gJWqQC=47NrgGGVpDJRC- zCOQaf9$upyQn7M!|l&IK|)0V`t zXI&*(Do8Z(7;$z8a94sOLG1$bEm4@rZV!i*ZBPB_c^n5MWdi0AU~{NCdH8K18~h(s zBPw6uZOeNpjUci!Vf&UXlYe_ba#+%(X>Oy#->erwNB6j1q{$jCJ0kz%IHZ?b|Ie3Sj zlEUqq60&G~Q2R5nT_-p&2+CcEo#vMc0?f$|ma>~vjLnkw^mEtGV}SDvSDB0l#T2A2 zr0br_|H~-;9AlF)D0y~$?@P~AI?5SDZb1G5>4O`nMEN*(|6@~N-(<|sR<#nQmyFBC zQkE@wH3mI72+i4+BtLc5EdW|6bM3_otsE2ks|&5epP!~u(pM{$DG=9QOltY(T^6r5 zJRF5~iTY=tyVVtRSEz@6d9-`)NwB*7*+nQ)vs(*#ZiAgl0Hp}Q7vtlW)V?7<#tbM; zstDcUkcfw7ZPVgYFey1;_7z2AI)c%eJiWlhXFvs?4_4~r07;lD9zxe0N!UX#GIJuh zWEHteRI`w3!ZfC=D+M*qeJM6Nkf818Td<-%K^at~(jWYk|et~&z_{VL5ZVIu9Rh)5fJ8*qBGcOB^ zAKWb89Q7>)3G)b#%&YtRL{l@w=RW=B7=#JYfI%dIg28hL@ew0mC~IMZoy5+;e1sod zCia0GC;=SBj_{3t?QI-(t#EmHsdN|K){{q!^FS1Ym#=pAOU^UCeW~-)76}1d3r-h#Ad-nBjME1n=MV&Ht=cLbC!Oa!|7@B9E12z3Y}H0eZk} zT`^9|mK>?@j{I%5l6fRWdP}@AGI13- zii-24-r`XYjEW-Ml$!P@!T-t*N+Avkf%KU(-AcGf5F`%x-dqB}(~_DuwVbb10S_)7 zlOqT>o#0=7D7l~&DEWtfrT-Zb{x*rU+=kmmwRgzPfbHA2Tx(*((GFH11D#eAtFa^} z71sVMCHOkg01*I(WOB40j$I~vNbHiq01Tnw*BChO6`^~H4VE$<*z+*K*G)j`FiBB> z@B?_oUx=T939;jx40PYN<#3RBh=Toe9%IlM@Twywoe~8H8KMC%OjK`Hg}sv;+f=ECq;V#C3`z$fA#QhK zZEfj8l3{ApHhl$-;*%fH#UCh48k=zWL|gLNztf+hbNWw$y5-29PrCX`^qo3ei+byQ z2eUa{T%nb;qc}M|ohG1QkRYu)Xx9JjX3FWeam!x3^pU65X}T9|nojJNKs+2zf4B#N zKyC*mC829d;2na(y`0NvD|W!WUUBk8{!h9`fN8K2@25S?-%N&Ob;<&>49C1!z}pL) zNC#PJ7GQcf#oT&d{ix^9$rR^T?bDhE@59{>fz{3I(Lp z0YSzG1?2!v&%7o#_{cwkJ#x>GE<_jUAULHlqJGuGfO@hsj{aBhH&5_p2sdR?oKjE1 zSjgLR5>X{X*~-iG+nRjrA}qhipA%5j%=s`WFNN_HTok8FQqbI+bIn1=J2q;Ewgd=C ze~Vwi1fw?R<&L@FQba+RU~XJvW#e!U?zBz-bo+f(as5O5mm+l1cX?W*3PY^v>y4xX z>&8AS2vYxh>C5k2|C-snt4i1z>UOf&X(wH8ZX_n<=z3yQk2%|M6mPZ`_ODPH+_{sn zPjljfDV=V9drbby`a%Tx5%l-Sl#r`^|8GwXT|LrQzg&7=ezWtd&5ex}T`p*S+HK4X zm}vP;fi%voi;fn7)grcM-*0KTeK{QU=dc`56g!;ufBZHi-nl-798gjV)d6v5grsuQ zWZRo~&;vW{g3ik8v)&3#-C82WfZKp^K2FNX_P*Oix5;MrwQH6?)37}CScD6vkmZ%8 zKIiqrmo6(b+{3{bxSd6|=SD-gfN`4>4XWtBIPA+|)Kb>}b?ZoQw0OJ;a3P`2-SHox zBa1V?sXxB|LlP!a;$ZT^n6dAechEn&kigt0e(O&mGh^z0S=SM80CxYx*`#UURmpRV z2TR5zB>b@j)EfwCuD4G7f5|Z7mu1SYZZi{oSL@!q6l9<}4na0{#nt8=I~+vy1cO|R z2}byBL6MLMu(AnoEhE3<@~o4mqmR&lGZRH7Mo+n*Nc&%%f*IFn&o@Q*_TxA-H zi`YocaX4;%>l>U?tf4TqadmZN)q~yKj~O1Wpc=6q?@|%e;x7>RQ-W=*z6fdk4jy=z zc2$Ms7MO}ka(Y~lv1m9mA&TamD?XQ(XQ`yH9KO1#cHtqVN&;D4??_m^;6uCJp^#Mp}$+@|?m480*}gdkiI&6q#%DUw zW(QKxP?aak2hCJdrbo$`H}YN8_@M$wC(xQn$2>$gSUCHDf}lpKcFgOlv*RI=x}U7P zh(!{6co?RRftGuV;N56?d3~=|i0?vHV-TwLic;pk*M6wiI{aEpN8_miZ*?8{GVbFC z+5Sc9PRdfTiaBo046-bF!ESru_UmfX)!7#XN?p;W2;YpfY`ROXg(*FOY z+xS{vgmqm#P3;wzb6JkHr@^7X$0_2oOL%QjL3Th1T#&Ba?)`1iwJ*RyM}ld`=}%Wo z(bA2XZg3aH#d&Z}T_I2*PR~GaF+nTK?o}@4$395an$e4?lcwm3rVQP{@D?(RU1A~X z+`t15wIq41j|@8--pPiwBm~zfa(r_j8=g4TDoFRYzYG;L)TZAQ@X=u~&_9$CoW>h`<}@)vcoeqGestbc-oG~tAeV@j znME5sboT3Pl1}&V?CrOIrJlky=0dT6|Je82XyY?Z&x2ADnfU=Q0`M~dvEy1&nAT)Y za(B0l61`#xt%?06wJsA_8-Zc8R?u7Pd~MGZUI%Vhq3fPdhF~LFtNgXJQkIAV+=B{;|Ihi`ycv6MN#xxol8z3Ek zzmeao&%$rkhK|5ucW@F8hC2zX_nLh%x%FdT_s^GrhMZGi%YjPwY#A^B7C8d_&_ z-moBf^iYFM-?ILAIL2O1=xEb^9nWR^xl@9ZT8nh`Ti>HpNsZ$x-Sie&4YrnO1U_AT z^7CgHy;|p1?+Qx7LUv6`eN5MNVLTQ{0E+UN;byNlFNuc_ZwV92v*4^9e}{Pj2Px+)9K zm%VAJq)3?YjDxH@@M|&*fzU_QaDT*6S7DBSmr6d-KKRg~b$a5qZ5;lk&V||;3idg6 z_S$^ELmXh-a8Z#2EAq_YJaDv5y}C+b5@Vnof}8>0(70YilJ*Wn?x;Q-s)|bYDsQFr zUA>yK;m$}9)^}SFG=(X!oq6|^SI1@WG}0`VYjA&KzGG58|C7SPrYUlKf}&Wcm|o{F zg9@2yCpKlYqLPw2xJ1CsF}m0kxtt|D%+WNQG094$xX;}&T+2!Q?N!L=bKRohaSgJ* z#OA-E?-PQr#n=wz6~i#Zs^atU@XuReHd(~3fZZ4#0ltO^)5F@@FzFPW|BVAEkyO^K z$=W-CgK0bQTh>0&JM9G_u1UeBGcxR=iajsszC4xDVCJjGWciIkCeXP)XN47-T7H#c z_BlS9;7g9HA4i+mK#JfCMEr@-i{`3E%p}a-Ac9MnIegoYy@tHHTVUnducf(`Yj~s= z_`y~ETrkz?9K*D7db2kYek9pT;T9#Oh77wAB@+mc-GgituxaqB&Yg;c(cn*c1)HOe z9yQw~((={BC*i2o12$Cb+u`-@JiU54sxX|!_9qN2fCz)fBt<2HVay85|8P#O7L&)f z-dJV$&tAS7P3q%>8Hqnla8$(oJv~e4{Ys<+ZjGjNq^2d^n2WjfarQS?LL4D=F9DqE zvX*Cu|K5uDe66K`)v4eKQ*b4g*4EzQb1-&6Q_vtBbiON|)olkil4d}VP;K4J7Xy3b@{GUC^u`O68@|NQJ%+PYJUgh@o);$ybMeefeKE3NS4hSM^U`@RaWH(28E0 zh0{degFE~@W2C#);gtVUM)M&?*HpqHO*s(1H^(1pC-_)eSx}H@!S+d|0iJ&;Tjkp7 zh4>s%gYCp10{9BwFDRXws-${Gn^0#B~T)7nj{Px9E*(HVVfX5)48I-rPI*D;dca;uS6Zz?*Ko(c1SC52l07j6>Q~`N#&_jr(krK57vLZ zcv-H?1G*}$ZHy*bObpxUH=Dgxd<#08xHQnQhg^-XN>fA=L(Nv=0(^nb;N4~KpkU+| zc99^c;xFeOtpCQsT_yoDuB*Kfl*HhjzIl~e*V1(TCON5+Cz-;5U+)vR_JqBqpLZXJ zbMKyY@+{@iA^R!jerU(mwGuRCrV>%776G38>=usM&+y0%q$nchch;r-(_EEW;O2lM zR={Ma?e?Fq*I>`j(kI54))(*O@S8^%xyc1Cz5+;9KzLd;+HEP|#F!MMjc+2wg(xSM z#n`Gv123dW!;(160MYS8EKL91&0QA&osi5d!TX9So%OPO`~;UJqZbK>lgy1H1L;>3 z{@6gsg$a(cak-~kk7UZMTW1?Ft2|(7%cCc<_gfM+j=dt`p0{FUMseFi+OM;qd_KRF z%;2TNevL179e{+m!%xKsTzoAiCWTB!!oY?aShA;KlCT@bV}PtFK1qX>l7w24Rq~7@ zw}X&ZEsO;8^N)Y0ELdkYJj{ayFRZV=stv6@mU9e z)`IPdWA-rt$&lrro$TBmswBgAr}3_DfTixQ=ddN(n*u}z9kI`E;|Q7zp=a@}zFxPT z)(2RZNt`0>F?6>2GXlm;2@4IxvxyEO0wO6oVV{AG^?#epYBt1aC@-7CXl?-4GUcZT z9pe}}I=C}GpMIBh&uW2{g+*SKM*+#IhyaI-;UJ@ikQ*D4xmMujKg7d#Bpf%7!k;Dt z5Atb0BhpV=X+9}<_=g^Lg21_m78aJ}FAx6qucIo&m^8&f+wEqItzQYATaVixw=V0p zZ%TA-y;!F%TKTt-gklvBcxCI$W1lGyAQW@8JIE|o zM1W3R4}dtz8N*ZpjFy7`Vy$jb_o22}lIKxZ?)pI7%8iG430!zcgJ7E`izCg>ZAhJ# zRdp;?(d9F(ttuxtDK+8t0;&~DPau$wK+8bO=XNf_<@l_PEh=_Dso%r9aDlt!Hmt)K z4Z9fyUcbBlfUbokjSritHu0=gN9|m=5^*(75Oy735IQ6Zpk1H7m0@gTv~<}rI_iH* zwA6!hm2uIKTlqo6jeORx;|<65B*;)@=A2(aHcsjq9LtYCX&(>w|C8|R6LG(F3x|}Y zhcPZpM~NgPRLSNv!qydpw>EiTwCEjSCX90Qn?>|BFEK$qzv|0Ta(a-83GN^AOfMW! z@iA}8QpT|kN8%4`C*+Tmqoy8aTS3U=atRjk$JsplbHpgz36Iw9%pS-X@Ud z`ud+fFKHc$atq^0KZ!9p9&FA6?Gb;86HXhEc!Usadl*0oeDocc1c-TRMghuru@2+VH)x}Hydi*UpoPVK(IN(F{O|m4K1>UM{_b6 zfDNt#hF?A)SOYN@AR>T8bix87BlV;*NIb1zE;zYPO9=xlt1a<9n3^)g^V-O0zj4h& zhjEFf>u~w_1p%I_!kI8^vQ5}z#7TkBQo!D1VS1)DfJ@M0YrtS6^upIQIDwn$|q@fq?sTsIqE!c+reXyNg|wjz7wrH2ogSk!SjLHsIP z{+H!om~t_uW0OuuSM|o}R)j4Bju#pAoZu*nuX%a~f^)nMEvE;1c+1GkrO>lK$o-x& z&~b&6G|O8x7i?~R>Z8`MC2Us(Rg6v zcxut4svb_7H{FEf%{>BVn0w33pD1oYQjT(MCu1)@G+uJQ`zQk!^bl=f{ zhbr}F%C0>tTNQJx>bFi?YNZNiS%NS0?Z*BaO6(J?DR)XH(ihpX_Vwv7=y zylmOmL4nb39QfI+)Ls^r}S8xk2p)`s(^S4vsw!W{ydqYr&)@>H!L!vOF zVAo5}v>@t$LPjGc#dMCC`&&V}^|U78^>EDia@PA&Xy}&qOVo>hakUqaEdA7<_^qD{ z{Q{&0A=j_ZJ`{TnEtgM#Fpm4ofn8VVQm&QN6z5kXFc4FbXR)FBT0i-Ox zmI8qb7fN4V{&nK?yJ#L0mRF7C&QKx(j zqN=FC^7_pE#qoKK4iQF1vkpvSF1NO{eBO@x3=Msn!P@#%D9xGy-$x5(_KCAdlT*YR zGL^Cod7JxMBW6(Rxcqf>+zOgSLTiN5nE13q-*6|w(pq@F z*Q{{n-U|O@E~m*(m~QlWE*3to*iReep({vGvNM{j&CJ^SJtlu&`uiObwN1*;(s%3- zXp4T&VjO5{3VbcbOy89$aoC2aC-qIqMfh+L7Y5@6_}?L*>54HJjGpKvRBAg|xX*x} zT@_q=Cl$J#88d{k0!l)R$2TZiRUG{s`jfIOj?8%Ao8x$>SFc-~*>8F5l(NcW4?XGj z@p-+?%j_qwzF2NHra*1_fz)^^hVcm^^Efv$;)9WA5|I8D@_4j&ZYy#7H(hOW21dD zI2_y0fMZ^f7YiZY2YBlkV<9S0yw;Swt_**W^!)E>*uW;Vm;@y`6i~mL2bZuC{jZQ{b{^5J z8r+|Akj4rD4t)AwD09A|B`hY^O`$5TF=jtrv_znbt5vYl$S8vLzCVW}L?r z=Tmugn>Pm@t77+!Z>+mFH9Xd}P_ZNVZjx2i=4s5(*L$|CTmB1`&K^828jV)^mVU`gnp^%;DylF-pVvt{o97pGa^=Qy27`0>lqd3yQsCIUM9?ML*3Oi{O-KjF*5YP z?3{TlT9yxXl%$}*So(!f_ve%{MR4bq!r}hQZ(n%9Jpqxxb!0W>g@-`d8_(g)22DT` z;NN9JVmriZMvC(j9iw)LT1%0Q{J=7RF^df3fNQWl!X%Yx>px9vAMI;lco-Q`p`m6x zI;^0pdotW|hmyUMA@i0%&YSXmFLHN^G0V%-?_3ppqhP& zB3peBUs4ztEHeEg1op<`SNFdGD3%Li^H21u_Ttg$d1@YRH$b;M5kkJ6XfhRhuG%4~ zam;}I1jnZ}MV>}RoYeGdoS+@sA{M2XZblft)vBVH_}biR5zhd7rkGOJc8OZrN`aq1 zqsh_m<3HlZyZV}hFkg4fKRJv_0*;h0eQ_ZU(v#?tPvKh;Jr%SBn+K)I2?j z)Av_3K_QOA_xC7B^KH4Kw*aES1*FfeN!x5w;;wVNqeckA#Anb6Z0VabBcexHd0*x8 zWbSU(-5Ufpr_a5e`S%JRbTu^!)PYj@&66bxiB?yIS*C*WnU z*=+2AO;_f-Ocab9&D)jRnVJU1`N|vqg_xJHZI3z5WZ~%xFQJFHamefQ>lx_|1zgWR zI91K@r76hU3t~1k~&tfcER9bqbb&550wKSE)lm%;AQ<3yd~N#bcWrI0Wzzlo>$5d~o{^!?wgp#)rNI2?Cb=VicuR31?EU_fz-3$q z)Ze+%?a>2}P>0{T)iPdJb;T*f$jIlUF#WCxIVU?0P57snh!5UgRZF2)G-VgCRcnH& zMKW|{+FoPnh#oQ-(5rPl$OjKjEc=?NShfKt*UsF-3cz?TojWI-(Fesb<|iz;->a{w zIY-2A+8^!X-@|)nvq^bHT#%07h-gQpCwa036kd67(GHHbxw*Q|T|6UV?RsH)z#joK z2Z?U~bzzBhqyY$OZ!}n@@tujUjf|Fn%;}v_hwxI#K zLkq&THS{bw`epN9=Fr$iTpcp0^yEu)5b{yZMJ;CW8O^l3cnr6Gd z&@k8#y89xlmdn3z$^M=-YkQ_*00w_-Y~a%kN;%Vq|PLo?bNE z=rSl5y2&AF{sp0{FctL^Kjqj(m?mEp zS$L9f?ZL`v>uT*XfKJZ*=BpBa};Bw5+jFF}`sw<@= z88v+q7w^ssP)RoS?9a|&@D8+<+;h@zAE-a`d7p3ZozSJCIAfFY^JQ_yZQ8JbA8*^t z*0yoH2z~EQDGzS18Jmk6C2Nw{43HOU48`7RP#(51Tv;p@vQfswD1?>f*4Dy&>mb%F zd}xmh2G9RytEQqN@Axz(pNpq`No&j)au*{L{u1VRFuiD=0a>&+7OfY<^pk!zfPiaL`M!S7g@cVS8=tl?*s3y zC*byaFJ+SKqJy4LT)YtQ^PVaVvX+)QlHv*-{s8 zb2*o$&70fOAd0X)JQ};$YSXwK&$~|}C!;u}lLQVGqS0d3Z6q^bArO$ip^6>Tyy zy0SZbUu359tG$?6zO7GqhI#)EXpmvJwZWxjef!_vw;ER}CejjuNIVJ*F-AnlUoYnU zw6TeE@T#Siv5zs7CBN=1p~>p)I(L5Nj29Fn*=Z>$PXH2E0mBEjeE1`rFk>Q$7%p#YpvcPBVhFT5@`Ty z%tr3|`Af&>P)rJue%hL-s%FLBw@TdeLSxa+@aHnT(4~@8!4mO~)%#vcscC9z7s(kJ zW!@`x+~#^(2+Gsf^c_i6^pRfe#)a$JyZxnJ6y&cT4l}m08fj|c6AR{>h%$ZxW^jv! zW#(*HjX#cx3Ee9jOuhK_wyqqhZ%Bu6@?Nimi8p%hVoQWaCOiq=|Aj-lAC$EFneBzw zj+ygeF_ht1t+BgUQS6Zhtwv3ZnYgD~BYq`KSvjnw6FdJAyJ8L8Z`a@75)%Z%>Rsru zUze5r94UVD<>`seeb|uw=7!9Y?cSR%hyA6{d*oG^p>b9=a|)IWS&Puc#0~)&IsBy?Uv5ATYNbkbi{FWLLx3l`N&YWC z@kazk>#KFp3EtKxKEu2_IXOIXUch60lgrGR$zkI02djSc-dtE$S320UrV1Ly2SewU z7T49$?`>Va#v($HW@Ql{Ms;E^0(D7-^J~5S;drkhmj`k4Pcb~*y8fP4UUB`6)}foB zdv{$L*N`XyTj7$NW@SeUMKMe@DoR>$*nYvG+NF&qou!OAuEwdX(?KoVuitx(;_jKA z+5aZ4b!yn0xQNf!BW9dz@P0%oD<{7t9@aaibGmxpZ`0EHW!KyFC#B6;$t1)p8d%0< z97p~z(*t{|^{|g-Xy{&3t}rcJ@Q(+AyQBZgl+w$g zOHDG_6^cmC0KaPOjrk|LQ zlbnC}PDpt6qNDdIXE*v44c&3>Fv1jr#&eiQrMk7LI(-@wuvcPAYTH^`XItB-Zx=U{ zbvg7SKNf8>9On4~v$EEt@61emKQs@?-r4v>{_ic?KfJ=i0zHexs>bIIw45l-e8qq2 z)UP1{kgT&1COo|To5 zqsh(J$_s*%c3SLOUKexr?1<;g<%kx)0(ft?bSX2HRlXEo>MnTySBP$9WpGLhOMmH# zRmd0TIV01hr*1tG{wYqrCg1fIC6rKoGGZPN`Uq)<}-#f?fjk$ zECQork1c-)PG#2K9Q5kix6nN458}lY*+onCU$~G?`-Oa&_tHe{T1`z{i+mX4_J@dX zg-JUfxy)SUej}5Y^guxS#K`#V+my0YE~kgLUwe0;`p_K7l$2jSb*ly8vQs=~YaI-x z6$mOOwj?f)&`#gqI!%7t^w`h)^+pyS!VkN*rDpIl77IY#tUHGCR=>}KOGNM|sza-a5hZ0iB4xARDY29klrGE^a0l2d6-d25gIUkDc2 zA7FC%sMAa0acL+mppj*tE=+&%B1_#a;&VnStwK-yps__;F}K-f9exbQGkW(G8M@J| zB|9VbZ7?^FXv?E5Ks9mtfpCWO4D-WLQNIr8s;mF;vfBOFUw74!L{`DQH~C&!Jq=A9 z=ksevXI@B6ZT?`FBp5iWptCbDL@;Ik1vUi$L$k`p&ZyffUKzgUs~XDkL3}<4OA(vC&yBsi}efW1pW#8A>;@9&URY$>qWjY=dp4M7VQ;OU+9zJ8kqRG8OCeCT^4X zwdhTrZZ}hf@pQ^7gnWj)-ASaoV8|4MC^PGQfq~Ji0%hfXHhYS+;sW30kB;6e8cnJ{ z(KbcNfGl~q=vgU5B= zW__EDUqUXu zI@^ttGCR8g1JD>zdw3?_l%}*!uHpUrhf7$9F9qIZ8g^(b3(<{?j*x<3;prcS%?O)% zB_|L&KN6=$f}k1anGwFlzLRN@c6R#)ZePKWriInjQucLH2T*1L*(*q)q`1{Lpuy-V zA}hDi6FHFx9dRl4=9kP0vy8fT;Ey|~uFa8^TU6~Fj&BRWw~@&STYl`Y(}AqGaGuK> ziQIlzF-df-ua?>e*?Vq)Geb zWq%hWw->YWz1E)4-}bsJV%`om``58UjQ!+W&f;4@b_n(#8yWuV2#1(V)2KeO@bIM7 z=87?W1_gzudlwgNWmTtGGk46T=7Doc`*?{eetsuG2{yP%Otp-F)%>wB&VW1gN^zqgg>gUoRi1BMRBeF4pAmM(fR3`VJqlWoa% zpJBU6ID9p2Z@AHiZ}X@xo|gJ?k~H{+yl3&r@O}H;?6ZaGw`J}FA6zO) zs(Ej3jpC6`6(WXTI`T6bJ}xszZ)nepFqzktmEI?(1J9-<8`i9TMKg`raf}#_M&7;h zJA#k0jBY1CmwHNs?rOrkua$2)>yu^pwd=FxWatPi6N_uZ0553GEVlZu4M=fN#1+Oo)RHnIlT)9l4vl?c*;(=#mY^$K9KD>TsW41N*UefE| zE*%!b+QGmE*AnnM@%0W{mqJ;3%&IAWM%=XaOv~}=3S8j;aG|ZOs^ZIT7lR$m{UxKBZcWLb6)O16|6lVx3(u)xfa;QYe|;D_So*iWC&ZRRGu- z%{ry5d?};h{+1h*C7V*$FO4f6m7KeJ!?vu`dyJWAjvzW1N}t`7E`<5!^w4W@^GEzb z{vW2^JD%$Q{R6f|_yQBqLHscJ@Iwm6eg~LntGA zWY7D0_q~6Q`~K(i-$(EB99FbuA3lnO3Ug1;=sAx3}0YW;lW!`)(w|12|P4s&awz z-k;6!%dDZrNQd+v+ym~2Yqp5>vibbLGA%|#q;3=1sO~UU_^cC8CMjxWh7noKVBmvZ zrXp&vkH{hSt`h$gz=SVi|GNr!*8s*9F*^A(u1?C8=JXkM)XX?LDyuh|4HaW@`m}XU zDB%Mi`!a_M7TY(hfz6{515D_y5sAzR!WAp>7|oFr;YALNzzajA0-~Z{!C)AMJzP@K zb7)S+0c(Hg+YjC0^A8DUh>TE6olX{b3WMRD=kCd8i7l5MUmni64G=T4DZWO=C}ezS z6`x?hH;ZOdVDg4&$kx{vdz%`@U4tl70#0si82#N`Z)7SeeOHPdmERndY}dzGe<;QF zfBeHy!U#RqyK;D;2%s}R4!%8?HWo*Y<~8z<0?GBQ3AKgJ2Hs&{s5UGS9W<8N=o0;-rXRkeT03$ zB{6qv?`ahEnG(2tRG%NNs~Mds2&#Pu@;R2#rq7vQp1(1OyJQHs8=?j{DSA5!{iQSN z+|-iO4?c&ZKUlV?x(iFyblBiOsYWj3@(23zB`}a6CEOdRS;(#cjy}^mK!hS7p%%r* z$1Nsy9CBBr!5C}>uID#$D#Gw|JL6&eGl8qHKD+}a^U3ZiaML{JD9TF7fvlJ^jXHM~ zWuZx+#V}nX8xE@d;WwU_PF=ofYQAbcHvaU@&!j)c!Iav!3QPJn8qX?B1RTrfiK3)r zQ7PD?q$y@IzppM`$b1iQo4fBtfa#Df49IEYFSXR^e5+v%Y!2v@#Q?I_GS|Ufwg~li zH`B2RObyM-rr4OLGmZJ2zBr-~jb_516s7v<-8fV2Myfm|z3(wLlaY>Z{N|tN{r*Ps znunb*o@_1@RUAenrwhB0hZYB?G`nc$_i}KWL|7q>jg&JRlE?hE-w!xKV&5Xn0PtRc z!2^xCpFhQ#O~+!rI~hNdoW2AW?wcC7J-J$;e&?z&rCQc^jL5d58 zRCkhWKjsPiR03h6j-uuRiMyJX5K6D@l(eJ3x{wD=^8E`@`~iuZU^O!a(EDdVjsb&X z0BPn4fahRKCGEY;74iJ}3b>GrgTTBVh=DcBZ}~Fd_w;~JHrBkj)zinvdY!DmlG<(MOUL}v59C2xVI_lV$w#WT`MQev zFS&Arp?Iz#(uM5@0Ri5d_1RoNPBuK7=ZMY8iJkVl5eoa%o@{BI^v(zr_q)t5@NB>n zdmOS1>LB^&30%@CH8nN2eeYlI9ZrL;=n~keL1xn_>2SFoC8Zs0?bl~Y^2&`h<4#1$ zU`7M@8>@RU*$NI0MV2gKsE3&}N~TV*Dt0eAhE`5uQS45;Str>~wp2|!{E%~b>-S>Z z!PO38S{-yP(DM*J44O- zC-G#wKf}PIflFpEz`(@?Xl0I1R3Ru5^_Vo(-d&@s!wcp?5m$?aZNHphJq_`JlP{)i z#(&ZlD|UZSTRm)w=tr}n&dK~bhBp@7oD8wHM)PXj3|r~%)7LMAX}SZc(JCtR(_AoO zq`jkHlOImpO}wON8HsXSE$A>%N^tgi8NWZpEmD$+c9c~v%X{UW^a#~d062#Sh>(Jkipnc z$JgU{7G4@w)ebxB=8a{GESUe}+L%2e@9LE0?W9SVd zqcf1o^bat$2WeaFH9oX?<=K0B-hzb(z7k4WU@qu@_<9!)0f z_>PS{1RM80v-AH&abbQ;Sg+6Jhr|ho{7|P9EwB68)FXTs2i_4ARG8c6c*o%o66Y5* zR_r!%#h^VV;|}tF8(@CPI7RsX1vGp*lyICUm`Q;_4q!6}Y~`A!(?-i3+{E_T%Z{5a z=B+1wm%;hNZZ>gwcP4eR(3bP9u=J z$#vP>p$zQJ$D)HuHdekZv!b$=1vIBn3xkwviz1DAO)rivB@^^a$vfoF{ON5w^!%3m zxrrwuj6js~sZtChl#2Gq04=nwb2tLJO{Y-Qd*QOu(w4)EL?TQtlE@|#^u1ePX!so2 zhahl1IB{J9?Pk@%&YU+Ou`!8>5n$}afcy*&-#_R76_8b=m5$IMy`*D7*>4~C%`IyA z?0cHe>$A7!#`j9#54JTbIwd-(0uTM(79KBOxkXnYBpl+m6vx}hqx?34fKE7<4&*yF z1ydpX<@ST7kt%~5IXU^q{N-pcU$^M-vvZA`Ta1iF?aeqD4A8W6KMv4HL-dz{;L#HZ z4in)BXMARE?q@enN{ic(Y$0F2-U1w*k_?i;PmO`{``G3p{60=rlz~Bv$FFcr14~ao zQhiw6-gbx)&cE;zx5J__T`{OoKq07QAapSNAY^eQo!+uDWjRHfmN~angc;^H{?rSd zs>M0jnU&d`>t0^dIzv}A`B#L<-qt05L~AN!ta+G^)pez+CP}s2$qk{*8_bWgxT~gS zFM*F7s5T-R1D!`8TMWoJpPwD~2Gy~w$pN^)aLdTBQ&Uru1z)awmCVN-A?Za0)NJuN z+imx`MHaH3KDTxIAOXUv3>E^ctfE+~X1%y)wzlAwA=Y8peDDd+5xNiQ$_*c3*z_Vg z^N69Z&Y%W|Mt=77(d^Bci68O}?9%9>*nfyY#P9PCE^x4WqHR=Kl;)fASTlx@oe&9z z2}8U_fPupsaEv5{kYm0;_wWoh#K{30gJ=w<9?0t}u%G)4uGvac(Y28f@nw`6rwY6o zMdioI(C*tV?5n@#s*jE{4LS$)`K_IXRT)IFU!Pna?}{hW3*1t?k##-ntqx=f(@XDH zqIE0*^olT6+aMeHP|Cok&6s@wVC~LYGX{-J_Ih65q({CPffSKPDU=@xd05Kw< zgM_SqCkE%i^3#nuwZzMGuwASQ^HWhe!3*#nRimEG3-*j;iunO@;J8c1U5AlB{E$_T ztri^E+aflI7i(WyNxCE-&aKT-6FyKFl_dsUE6&cgCW)ab2-PCYwkY@p9w9 z2K1I+B2cmF@C5zA@|HYb!fSru(gGtCx17uZVxF4Jx;EX#uTj)HoAK!x5H>(e2?;TM zFx8x~a>~SdkLdWgD2z6Q)9CF@qQ4BB;Jw!lP^ji)0-CM-+y!*RBP1mcxnj7LG7h@u zqvDQn()c0&A!QIsM6C?N~S4TL(1D8xsD~ZgY8cq0~D| zOR&UD-TN#@YC4i54AiXx299rBn2TnLkh~<`Xg{gq$MIm%SSBG=Mr(APg&9R=jw3|> z00*~F6i=>cq1uW+2T)U5n~&sihldei)P?@$XcaggmBY#Oe`Fc3GNeL$0cmhr^Dj)lq1~I&sQcu(s_$AO_%~MmfK7V(>{#j{L!iQ! zFLgtuz*ZseL$F?P9>Feulug&$W(eL!yXALJ30}D(Zb$+bHJyWNjUO!}sg_Li8 zw3ZWmKbw?KnLL6Y!X)c8UpLjBRA|$7(i$?6!RxOP&;-N|L8dxK$vU0-D9pg-qw*P> zsEq-3lxNR-TDR@Vy3qr^-m6!B9vSwRSnkyeaozlS#+_RlYiJ*x^hSjLk!={NG~)Ho z2mq!C*gz(`CM5w^8V5Jp{P=( zwA^7`h_Wh3IXotu*(bGbEN&QMU4me&$h;{Q3@{#e^U?wdyAZiJN1I_^|0_~7Zu&02wa@62#(qtR&a+lT@aqmlA0|NZ_)w8o#NBU@#)NzC-#3oObAajmqum0LT^UFe$Wjb^Mwo(jexxD|N72hfNh`az%vpgO5)O>T zH*ba|+e@<|P)3~xJEQU9O&#UpD-``lD+HZ6=p(LqY>v_%yL=E-SU#`=1^a;_{(B!0 z=$&o=7cvb=ed zQhZ)t^!2mkv`~U~4%*Qk7b{-si|Yk)(ww$3UfkvE^6&$3SqK>$jxQ@AgTv35*GwWX z@)Cq)mYJ^j_@biVcifxKLc&$65hdGqQk?gQ z4EVKs#W-NzsHv?6T|EoKNPGK^LZp7RA(21De!nKGm6&*XyYfvQc6M}b?(h_p$T(ES zL#{-`#PszwA8~N;a0AK@B)d|B2g7afGl%Ry8`S>5(OT;k=KPH-+!2;9UaVro{|jQp z^!65w)KovndB(~Ovd-YSkFl|3Bp1qCQ7>ClA0JJ3Y;}gZ21{0%JH8jpM8~_p3#e8N zl{H~zl1${j6pDk(N-JO$4BQd+g?>MTTL}IQsvHkM9*a;5aVaTL;B7V*U%SUU>-q5- z@b2=>;m`c(PCp8|?(jpP#k`tlxo12`R5Psl3#K2OBa>^$r*rOnjekLmGbG6_WPyDQqJqU0v-MWvf(vR|GU zQ#ibzKS8|D;|Bg0PmE}-K9s8cAaKtlIz^7SEStmme#7sM@MWo1jNMyyOmOi7ymn=F zZ3L=UWox~^t|x1XV=hh)#lYYhEFOjm=(+9=uh94qZ6|CH^pFNi5ximsi($P4-A|$8 z`QDTLeIG7P3{@uyLo@&72o2`m1Ofs9Q*{LbintG(uKfD{6D0M6Z7hQJ=gfZNRYmK8 z6i4=tx(6hyWe69Cm8AGPfD0nNMwNd2!aHnebkiyQo}#AZJGlp7%1m+VMbSZEyY#F8 zdyU-vcU0#&dOv)44O7TdYkNm(RDGPltOI#tBn&x7?LA3%DG*xUydFuBtvz z)U@GUUTLb6=zJ8&3149NXIRe+k=Ce?$>eW$)HXVJUN4}!s>aTSy)rM`G#lUwt4oz$ zgPTjxY!p~J8L?e3D53({iCBQb-|5aKcDv;79Ngq}!2UtT6Y!ebtsa2_Q~)-v`ChVQ z<-a0?3xMFJt)r9h^UAyOdz##F=K^+BCY8$NzY|;Y%33eY+|9Ir&8^dZD=?($ifot% zloAWALLe>2k#PQIZbfI+HJJ)$=gWacrH%tzD_>o%co#UoNrd$8ekEUAPx|x2x9yYo zY*evuS-?&f)u(vV(>no|(A`aXaFkE&%T=C}>Gt`M3`Ejip~0wgPWyr6i;(HruvI1C zX;j^L6PdrhFwkHu5R#O?h2rFm`}gl>WWNV)BZ6^!AkQn*?b!1tg%CQ6n zGR-S1f(|Y&qQ|K0W2Co!hpRM^rVgZ|8%D8bfCKtx#lh!A_ARamJD`|SpDH&p^IGys zeR~@+r>&`rLMhCfqPm7i**TMh_nAkwc1LZt9doAGfoFu`ZE>6P9+P3pM8=K%Rn^U0 z4j#}tdh%)t^XO$Nm@Y;spiNJ6tHQ&*#VRS}Q4 zTFjhCDRP=|Z1$qxTjRQq*?w`BG6NWWb_Y=Qo~NS7QsRuqH4T8nl!f%3m9U1PSP=7Z z*o@mJo?QDR`ur!^)+bK2A08t?O>>$18xKha#kR;C#b|Lc+*132L_vVxnhqc_yu#5h zIGKh$Slqs_uF@5Kmyzl1z3wFiTY7D)moJ-B^9roP8A=zui&%4@)_s?Ol2Ny_UUxCb zVNc9N)K1SOCmUaIsS6PnH)o>ht-rBtKlIEGfYjets`)L4=$A+Jv~MX|dsPl@r;TQnIjL3D_?^n|~L!8la`y+pCLX{r(@&<1muwb)=Vh%%t@|m(!voa{9!@6 zI_W+iEyGvRbnY@*=xp-aP{Q?h=g}0XX1~MzE)J`8$3by4eKZjHqOMVX>(So!D?h2a zEa_Znrr>Z9h;k4POkz=o6G0oNZRF(YMEYp5s((xG2@iMoc`%YR@DHF9);@|@K+ z!D2yTUkB!4HajbDhX23hS>{%v7B-OsOKy&QuGf2a+1g(T=d@YdHrD20v ztsBDgj!Z?xnhz(Gy6}U)$cS{XdTK0A&pHL%_}vo{xG*A23+@;lybdcy^;L7gX_F55 z#PkoX;pyLft~<(WW-fGx(W=^*XldD(UAd`&gAl^uNCMC837L;TNK`PT=vsCE9bIhu z;+v4f)eF@;tYnwa^mnK+DaN7|7tSJ!LUGzGh|Jc+?}^f#M8Uu$>rD)^Ga`KFYH5CO zH*dhUKNrMx$o7tJ*#B?4exapu*BqjXxmdl1pJ}QDj(1FKY#3Mxks(}DV4y~T%#;x@ zs`UN{XD|lojCna>;tUGl=<&R=pADdWc?Jn3zk%pf2c$Oy(E7F(c>K!?JhQ{?yN9*h z72FOhmSa(}!LmzK(Bp4T!h%z{fISTqyK`X#t7u92V!@%wNp|Q@9ms2TT8+n~ybHdW z3>NaN?Hbt^5<&QU0S2!!IC3I0Ulfp1iD~p|lg$aY%-_g*#VWMLk!w@lc zI*^?au#idVSuEDd=nSa)%yB1_A(N~CM6-m9jNcgw`+`WW2iT8N1|9B>R)6dt8j9%a z)6H)Kk7QFTE7JGd{7tzibyVJQ1IlPI_ z)lpS=|Upzgcp6Q~CrXCJPEcMmuJv42eViNbHX#QC2hX*6Mku~!%4Q!^6k zEBAL=0;I*2 z9`s@DvMvTmzx%9}Orz3Q@Yi1SAnOn-ZY!X-E7*`eb zaK>_eQ#fd1M0X!~0D%6bwq1ZHWc!tuW}uf=#O1alp`*^-f!kI{esi-KD>|V8FeGl9@D$J<=`=IT=0lWG!CtMB?suL*q5>La!vmYR~@gC;O z(qP~L&pBSuLGkjnE~&y)NXLM0cw1i129%lCKPnZKY(0z%)n-Oy#X|}4X0H}{^=NBd z;GzzM@Yy*btL|2t=@$&R$BhpU(njQm0RC!um%^fzgVg1am2qd+=W7utwWA{Fs|%Mw zVQV6G5=>-@N+HHx!wz>M62z6DWm#vd*&xSl5Thj8bo21RwtwzgMY9S?DB@Gt5i5%` z_LlEfK(X=e*=H_61~4*s28T)y!o-zOQK`~%*4kOSUjtu-+R%__1!^n&>VJkff_e%J z_t*bWmi_jTXXb(!tK|*ir>h3%zC{qFt_$IVmyZp9%*LPO2;%<{0rFwB7iaC>E8dz-|PwBrR!GmEW(V+BI)sJD-;a`2k z!+;N^g@|Kv5s^O7BwC~29$y4=xxF!Lyo1APGW|oRB>J!yG#JtNYfslfuymVNf#w!e z(T6=WALOobx7qsF3u%SmkN@9k`~f9XP5J$Tf&!=LIoBhSYCk;XxFhbe zvbu|p<(Ia@ADB*!_fcly`%_ty7(;N_(Qpg^i*#P{+nnbNach(Dnix2jym2LNX^(<< zZX0e6jGI=dfi|h!xO?>!TqW>uuTP;y?ms=OyJ1j$8i0eYjZRJ)OuceL!NT><=OL6f z3gE;QTvl~fK=|mWKMMFeYaWL>S|NY1SS{9EVC1mCed~fm#AL9DK*nqeQ`g>X9dB^m zYADfDvaSIB_f8HvJCiIiDIL!}+}al6^EX)G;9igb7!E%$i8}lJkj}$X5|=1smeNy) zegEo0%&WHzH-4ABVD$|nbSlKhZRI*YxF=lp`1w!BK*(F~zQ77;fHr%?SiY9SLssEn zU_J-DsrH*}YXzJQuHD{n3@}jB!vFw{35IHASjea8725Wb4kOGvZ@%?=FD^bIO)P1C zdbFkhnR})Xku?-(dnc#HUd1OQTWf3BOhwqji{rW_x4tOCFFZMUy*6-s7|!se;qIx%Q?T$GlBQU&;^L^Y^gUlaw=OLM8gRIuS{f~9 z|KIYUUEfbTE&y@Q&W;0#D@i!d^A;oLF^@T#3#Ri5Xv0cfqAt`kHWc-b0fGEZ^$v7) z&e5Mx0dSXHk^yM()bb1j^Dc=AQ+>e74XJ7`%05~Nm(h8BIr=oZP4bcyqd1-#*5_K{+=x+ypT|+ zUaxm{etc}Yw+o2%l-k;k`p?6|(Ae3hK9P3%ntp>u>E()%NL`eS7?0pf_y&hU_|)jH z@k27{wu0f^|Awr#8o2}ut_y$gT%f-!37bhdP`$Ob#9Wkg^P5xwd7fQ!V9&j}t(CM9 zUB*8f>(_E~C(<~bbkhMP(!SSe0bhfLrP@P>O5D@O?zCKDG2x}i(DVto7NL7ZS>iotwvUSGKm0(}NXD=57C0BHszt)*&Xaj9DV ztm6_Z>E-+AR}l<6iV&35#}pl<>l+Fg(`E9(cd?1%Y#^QZ4t{mcKxIN%;4eD$=__#3 z;bOO@Y<-+(FX$`dr#(YuJkX|dc^pK>a>St|>+8{Je^fCX73}uJ(iZb=4F&3-0jB_m zv)h2xzcUam_+SYEI>Lm1XTHPf0-NWDTrHD$mM%oWKSb@n9HPy`&d8pdi=ORPRu6@? zkA+WP_S#j$$*CRJ-$P|gz^|Pn4UqvrB~!r&`)ea%i`y4yG0MuT#cSXzrb9hH?6B_- z_eq}7)7BH$lb?Dml{6H&5C$8HIN`~q@7ggvz-#&uOme0D=}qaqAUT#>e3m2$w~h&J z9Yn7NRDoaQLQwx@#$wy{8^DUBVLbf^8Z0h2;exi?%+3z0a*I#ak<;DVgLRyQVz>`kvBe${*eRc)SQVMgohOM?wqtD<2m3xyWupyC&_q~)h{kL4KDI*Y{orW_?9u5+ zl$r!5iWk+)0B#=Li9l=OL2HszlSuwh54tb=FwT#f9XO-Pe~^ZV4>EH1u2?_wQ*N^F zDlAA3m=+=Wl zZj|g&70dta80ZrW3q_!X$k6aMA_1GYc_U75$K}e=uan}Y>*KBiQe2c}>$<2s10H@C z2E5W+u9?4`>7HF5ZHj3b=`{d6;KImVk!=4cw@M7ErAh{u3BcHNb8&GE*Wmv-)Yxw~ zkZ9C1I$kkDhhO`KzUr@^jq6F~AU1=`erliphvV=eqhFXC1E z>{+ueSV>&JEG#1O1vYz-UP1lh&%#_$mow}dwC#2-|F z2B-P{&hk7HSD5O;SE#& zCc1G5JO%wG!q@r7EnWu+d2~ZVgWKZi;|9E4u!2pMU#0VuR-Za{mP~$|s{FwB!KM|V zG_O=gJ$X+D#duDLWWMk}gwxsSFH$}3o+CM(6nz<6*b=+oHxCq%H@5wvyBw4_E`Ni%%Y4)Q5abW9NAz58#@H?S?f8$+ANyn&GJUbm86cazd zW>-WlB`7>xv-%{~cax?=ffI6t=d-F1nhcQnSpqjJM=kcYR$gCrV0NjPJUI*&VdXHw z7Y6Kl-8ookKZH=PCjcRy%CFu?Z1ht-VFNT)mvY4$EN1-VZjR7XQ~cpHJ6Db_Sl zsj$1zGXL75@!Lpr44^+Vm9<=g!m14==9G{GtE=7&T=lH9PqfFp%Ub;pDw6LG$$03M zI=y;DK~23Orw)2E;20!CK-tC19PjV8S&2mE;<{@;f~*&qj1;VL?wXoL5s6Y0lam(UtsqeFSIhLBKJj6` z^TP#qQH0J-{}xiZ^kN;_&w(m`` z*A?+?yrVzPz^WD9t7e!86vlx9EE#NEM+3kZNg1)0X=!b3aL9NwwdV5Qm#9_u!6CU? z-v2y9fdHYR5sy=L#d)mvGILe!6ITK>`7hZD<0Ih$vqDcxF>KA5g?#wXA##J$pNmpI zH88OkeQui`FMnt}^p5XN!Ql!m(=iR4^~BA%mnURaj!j1b5>z@>gU1JJsgs|)NCVlP zaCGVYq(qg2MT{%71=nl(@5?7#1ZJ|+hc@2)FmNpZTN~x$QW8eiKz-c;H@xQ%;R>6g zZ6?SUjA-S+bprR`{C{LeljDA?Lg zF10T<1jm?BrO$H8)`v0%IMv-?|LyMzNL>azBwUR zv$JAiQ&F+h>4SsUZ-)6(zvoA5HrC2bBkSgP6TULP8D3Gk1yxP_qtI;btaRmUgLxM? zNEsul`lPu06~Tr_TfvYWrIv`QqaH8EHwBjvTQN5P)lR`4$=Z7m;z{0+fLoSZ~`@8Jc z8M(LK|Ir~*(0bC)%b(TjIw%n%KGNMV!mQ)pUUcv-3{SiFK%rww@doF)Bx(^d-HP2D zy_hHWV5&UcBco$=br@+C+^ZT{QH&fH_*+u)Jas{!ys;EL%#~H7*?98<_b7fVivl{^6}YDL?*T-XIl{Q&Z2Rc6|+>$FA)Pdp1~{{mn+K(+a>~G}~lwp_NwjL|)oS3#J`d)gifApWXkJDmj zD;QZlDFy?VXU|wMw$6h^Agy*&%xFWrE)lMX5z$%E(vmxH<+<4ycm=qdGhr?J`TLbx z#~W8A9bZym5~sceKGKZo>6<#dip7e{rayF*IZih@s-tD41_l59ns$J;-rUjrCycR3 z5_yQzM{?0mN47=}B@1KkH%MIG-g_|Cqs1F}^;FMH9v?>S9_qNRB)X)_-+0%%L_*=+_S&CQG zw1+^YK#sY7qrc>#ag9VdJJifMV!*_+DjKhucS}3q40aEMX>?acJXjDx^f1u3x8{pG z3xI`y2*sfoeTL@Y1sKdw0MxH?$-Iut?qMI84{^6)%gF6g=;nuK)>8Cux8FW3t4o9t zFN?)V?jNI{Jml$TQ4=AosH4pCT?+9WA$(*9y4xtGXW&S#$PNulp7Al#&_>{R}a204#B=+}Lr(1+q%F9vekWSt=cXr!f#|6{ zM`j&=i%9nNy8m+0wf&GhtUov;<@mD1T5i9D;VfT(Hg0Wunz{T0>*LLpO3&gYN{aQ} zhCh4hZXM{NPr%ncq%F2J2&%35tb$t`dDwm6x85u6?AJ*1H{VlN8|9E#VL>`uq6Wz4 zGetZG`5lYpzmFjk~MB%(;`A~qRwGNx}1j~T#7;J=;my*+*Tv>noo!vTh1nb%L*-kwH_NsfkBOYB)5_9R+e{VUKW?gQP$ zlRn0hIgHSpCiD=pabD>k{^1$Zt9AEE%`JF%u`4G&*?(hx=8jX+>5{xM-t&T(22_wD zM6Ik!T7*_giu7=`L?DkOGS?Yu4dyGTr8_-)vzLLJr5R)0#7&<_1G!U&^wrmU{x&ae z+ZQsUzT8mmp6%6q(Vu?wlI0`G;=|k{cNXR(#FuIqB1?rOko;qbBLva0XyY-A+z?|?g(=K zksGr#<|Mb)kyf}WmL<`<(c3o#MxlC095{27+;Z>==`=fy>#Yk2M^|&3+5StjeSxp! zCELioip=K3I8rDzV1)GmW?o&8aX!;^#VrW|!B4zrLp`ns6;nt{EI$}P^o>UnnhPIZ z1i*@wkqn_owqI%MpH7vnbhQOrDY8GLE1yD9!1n+=2sF4JON%KQ)+Iyl^3>E#bM3K| z6M?zOF-mDb$WSD|WPf*ya}<1V)hIavPOvG+m$==j{xB%RpH*RCCt4RRGd>x7MKJr88ep0hXP@<-th!jq0rWM*%QKrSv{dogoEIR5BX z(4oduuCagNg)ba008Jf-%x348QkS6sEd4wX#=|$|h`Me&W^H>q9~mRn*%}Q1M(k_1 z;fBSf$IeVg?jiw%$HIp4hiA8Dg%d_jw$=*5+RyQ8vf+I(!>m_3_AahkVMIAF=ARwT zv3qiZTFI+~NWK$>La37W7=9ICS&$U%IbTVyPVXaSVYCaJ!74l6G@>e+_zTs6t1Wdt z(=9E=WiB1_W|`LxP1M2ZdVq_(wZP)5%XB{X%$Q1*{e6$Yfx2rbst}qk*U|{ zfj!pY5(9qzLC-yAB50x|qEqAf`B!>-7ge~5If9IJ7iFiDyIeRyMDVz#@GjUz;5>Tn zK`zPaBy*&ZOtkB|V?ua_?zA?E@k}J`PQ*{HtO|B}hnUefH9&@$oweq@1#(6DWsdu?qw0v`nk*_fWv%CwK2Ytf#FjzYiQ1*hscIN`HS_VKBts z8%d}MNPpATdj3cK1tUI@u-PB?t6>8b;OhL6Bn%B<{mX3NhK({i;;eh`cl{W$_4HPL zv;K!X*MRcRo?C;vkiEQHT_V@a^NYc}s==G1Ow|o4a-;Lk*yn)hp0#bK0AbZ)b^j;& zPrW&l+xezf{5D=rPXguO@u??XslI@94N>T{8Q$*gU9C^Zz5vuu?YnPXgP@$hI|677 zLiB2AxPfdS$@4dXNi_!Yp(8&`AGg#5Kxt9*sdy+oju3zDB#JBlTzhKB{+|mgYtgzh zmXc77{rqfs9>SRG1`2;#;5f3CZI?`;ENxD5jb^=g*t$&9Bhn|@DlQS?c0Gd9m2;ij5ua*X_OpL~6 zd&Mpz z2}_Kb=xpI{(|s)Qz3u_;)$4|&;N{XXhs4cmJUl%ZWQw{?wtqsdMjquetiS^de1wk2U0D{-Yx!Iu}mK&&^sFcs0>? z#qhM}4SxW3fMS##ZlW{o<|0wbGwBf(a(Vv*kdF{UDTo zp<5I9sopyd%um^JaTUlXHcZ^~448od+gP*uvjqu`xh5DJFF`g?a?ruHVt)LqS5HFq z=VC_7I-#$Q5-~sUm^%rT$EPBhu=eMCkOn4(yl6a8g6N#!tRoo@9$W5~zYJp5+OJ&K zid_PAS_CRcgzZ+An?M)c@i4+b{42h+;=%s#{8}2wymKT_9sOxKF)<*^lSSZ3NF9p35K@U~x-9UBfkT+JilVDT|uXNlLYe3jtom>us;xYyc(rB3;3sOP?WOP$5 z#=gJ4dr#~40ac1+EX;T;Fj0af3)@v~nmwcQSNZ|gYd+qxLUgc1>N@4s-5qy8x zb5RQVqj`+@h={A51TBe`8(3Gs6th1ro+ZrD;V-}VS2T)2uzP>xld>uBg_)KW4q3{T z`TIV*<9b#`d&*q5LJ2l)H=LYc@3rEBHTtCeq}j{Z@G4xGF7ekN1Iy}9Lqp3E+_7-O zB%KI`lW9+@ZkkIpd>=X~Ozv>&$-ZKDl@B=jG4C+2TVf$E`ti5i3kE9}dO}F*kVYb< z(#B!;ztcB?oc&a|VmE&(8JKqlP;ky&B?fQDrgHiy-XlV$pEiFD;!9u_d6uKv3; zSYarrh^G^wp{0KxymM`Bo+o_D2Zq^&Z!u|gtE=CrG0_U5IzsU#wj=Ul=b8WLhS3zbW7heBwt>itZBmU0JXOW-3_N#&fMef+ox7d+TM35EnR zXb>Yv2$%KtV(y->*QvP+ItB)kEy$&n#3G??o(g9bO;0ak+@Q0#u)z%Uf406um^?y7g~vY7L9A?MKttG z9Tw#8-@*ra(PI8ik&jeXo;ojZ;bq?7hLogzWoGl{=JN=`X~~z|>7Bo`NDHX>;&YgSQX&KcGD-{XIlXnQE1knYeuGzC70p#Jp&hVZ-mvu`;o#=&hl7~S0 zE8WZ~;Kp#OQBHMGgR;(q9GqqL1|hA!|1DOT`c9fQFQCiDq6Vd1X1JaK0jz5p+^>I+ z?Gg>f+L%C%NjwQD7x0aI{aTk=ehVU%f~TfB$PH=(i%KD6$C*kuPy7mOHGm5AB?~b_sVMBTL z7bxmv6k<5H20thzyad#!uN{U@96!b-jmbh@*sR8~wq zoauMo(0k1M3aA(a1Pz>{d*Vp*4EP?{*NXi9Un7XT(DCv7Bs3kq=j{T7B4lu5N`u+6 z&=DWt@y_#C>bV#Z8Q89L>?OyGnoPXDXncg>%+#^ftMuS0+uN%+IOvRh-(|)xB!rqM z4NHM4m>#J2SrbQS5EG7z!|VP{6#05n%`0$)z1Z%+xjey;u{vo1IG1YrpTCA|+huB&upt`RPxG--D z5YO$xd+H~^|3jWt`pmgP-;>_(MFpkzY{ z`@5hrScb0-$SFnd3PZJ@* zI21kyz)%j0RVj1%4(Vm_v@Tf6`7D(pcX?*Az14V-gVh_IE-GGsn(C!Po5hZ)p&5j; zvBqQ0Hn()4>pZA*?QoDvrP3@L@Ts)+fh9B-x#!B|TvN0iXm(6kW?Qj2daD1iI(HG7 zTU)WZRG0(X&+N!nkUK(=0WWZQr1JLHNK~mMX{>$)Jf!H2{li~z3Ff-{Y3{U zI{?>>$|uBF_X2Z*CU;gQb58a<0HU1YFL)+`(yU6Kb=)CeY@;Z&@h>L7W6oC^o?p zLFwMJf44>m4}8$J2tpmeI9J`30s{{YOdNwtQ$k3IG}jKVf#w@J?2!CQ{|e{4Ty%FE zc&}CteWS>WwCOkIi@4+C^Hm5B|EL{sAlH_d_I5|>rc)_vw_iwZ8W;(WtL8G?AYMfq z7F?fVGx4?%e>Vfo3^pnG!(QRLF9zu?SPBY}7JkJKcYtIDK#nA7% zlvF@SD&)Glx(W*@X876)$bV&JX>PyG;D(Sf^G%nE1I8;<@)d?P#@TOnHm}rl+Jbmc z1zm9PS?LXGdx)t69dYp%98SVEP7fhvIwSQM1S8bNYEUfIx*i-jI$34pF&4<96R$9_ zi8WM;U?i)iUFDWz6i1{hPn^1lvU8a}=oEP}urD#?2_ZG)RIx|}F5=SmN=EAi6&n+%k#h6j9)Z{^6OlZHmXAN2VOj&diKI`WgY&l2VCE>%(wog@*2-()|eVzNM^uCGLVSBADN6~)N+Mmu`7*# z3Kk;hWh$8T{9daU0`3IbDn(MHv z3ViDWiRreCSH2pd3-5sZ?W2|9(D%`hC47xa^c=Vd8S(IkkHeCw zB!7E%?G%m+1A#7C8Lg?Vm5gi(99!F+vYqw%r79^+KzfnlmNiPr_nmqU2uco^v4L;2 z0<=U&ZyP1S9t|u3CT~76+v8qiI%b&%?kP3;$58=$g^xUeOjYY*$$iNlx55SSPWxfR zxEP)VI(&Ec#OCz}52lykZ~cx$oO!g64;E^(XOT&Ejp+?)@I&C*Ad%Kef;|bL1awL) z=zsVxhGd^vD$sAmSR7gJHUAg$67Io_KMB`)5|&=Ljx`dJ14ftm!4P zn^d&JjfJm8{@UslHd5^O2-MrZ`#?%su|`h09=b3Ty&8%Nvc6iRenRKUiitUknn~F0 zDa)tg2vL8feU$V$e6j7xvAF*nHk`9=Xd+Nf`sdGgH;2x$cK=C;wbsl0?cTCqRQe)P zS2Vj%ir?6a9efN{!LjG7YtWXPGel~HZMD;2{5D$|-WPZF1)Tw|3YdNzrj@XqV#QDo zeSJ}%e?|q@Q&_E5E&AXVH0<~I|99%Od5y|XZzwWuYUg3I4S6_#mU!8!S8$d}3{UGa z$A}NE(*g8LrF>weX0hgs^nv%-kdqDj4ColzF1`VAB78#ghH@hOV12~SlnkxBCWkLN z`_|eq^NkZu2QM3hv4etn{9qnmdqW1hwz)T`B)^=k<~@cp(X< z@%z;@?!8j*guyK4T^AtEVgN*VfGImfjDPv{<$2i8sp3xg#z=OFOW-Jv5C;&qc7$d0 z_3PIzFK*`0@NP{yZEf?n>h#GeBdMTpMZX=@93PA1v^+Y>fPVy0HKjU2Pa@%L z0y#%LiDYQRo0%vjQnF9EmOp$l3M}5hIM4@TwWTP?macpEz0$yf8BG(F8 zFnF{)BHrBI?kFQ>DUyOqDtJ+!1994>9EL{MQM4FVR%Wu9?T6Hc6+l1F0vT}g;ewXt zgc9m}{U9XR;u7gm{IZZWvOGHURfqIT)$-l7Q*c0KakaF-`Behbw;IJMU6hYe!S$#X z3_p9N$&FC|OmIL=R4@~=YS#`q$=dQ%us-Fz&?@!zSM~ z)JhI^Zl#RR=D{3JDI=ph+tuta53BkLZDdrieb1_iA~k?!_+e_{eo3OCn8-l5n^tV{ zzuq}rSzZUI!)9wW>WCpJbLjum+nYyI+5PY1Cxs^UBr3^}j*w7MDpLxPA~Hp&44KEw zQ=>|fA*qCt$eek0DhZJ};uw-Cgv^=0*L`~4pJ9F1@Av!fch<9>W$C!@eeZkk*M3dc zb&(CfUcLUtB1i{e3`K|q9W-%K!d=;F?pikUe~iq-MKezI?Jp~Lp6Gkm3=q*qrr@hA z9gMnJLQ?I^UaNG!pIF>MV~eCci35ZSVq#-?z%^zC(izgVmF|o&bl=j97`?H=x}RC5 z$a;Gh*9?grIuson`)R21eh5+EIAHr8rKi=!J8Uy8ErLS|Q8rf8&q^-5K;$3-cC09# z-ED`%s3PQq<)+MMk?7Eh8(!<=8R4H;*7J7wPjFZqB}-`MdnxS-0nie~Z(i21v{qK# zSi}>bzg;Ko!yivF7^OATcBa;xi}Q6}_M3}H)v}%{D}ONrv5F)RIQ6H`)#8@QzXqn0M)e%8HWgnGyDjQihq zI_G(9nXJo?fAsj1r*@~~5Cwy)0$o_Z_{6RPo{HXjaog9nhMltY&4$XX4-nA=?U9gh z{bTO(raLi9=v#d=OL)FhLuc#bDBG0$k877&53E$Xq`X&)e0w`t=mI}iHLzYxY#`d$ z$>sbf{o-}22i1my{JTy_lU8dj8j}fb{IeQNoYClTgsd_JH{D~VAlfqw56J*xX;ZL@D}PHIwhQ`&8#3g)Mj{hql}^JZJ6SXl?> zof(uKJgAp$9@~$YZR@UNYZVl%`u0EsBTq<&TT^}|qp$P)hvQGBvaG{LM{A-)gi5v6 zUukYoxI$J$Vk=)kkkHH6B=7&fk!c9jj}uTP=#E+}dzhM<&Y%4Sp^jUMTt8Mp0+zCVtyti;V-^;BARf4AxKqYZ<)@zRBt1!Mk5H(1A%_eT z63I~sVrG^hs)oDn1%>|9aLu(0^DtoGzM0LZ;$P^aw+~!6bx2*@3z<9eZSe`q)A4cR zsg#R*_G}}|KSNJ9!M;#9aoor1prt+fOqxxOb?ow-6R@LNzwGoWE7Q#p{_W%sSjH%K zeE!i%?t3L5<@m6Nz@_ad`$tR7`4g+m0Z~zV@Cj%7a5z@&KRnwL5(M}(`du;&Qy~Zg z1@J9^s%^>48U_zqDJcBM&R+@D>`>R$z(=*0b=cVx&%znDM9Fy2AMOCvt=iLc&rm1 z#9|ScX>xzB!6plU>=R95U(cqn%*|QO)|$6w){b_DBq*+ZtDekDlm<|Wr5&3;w!x?c z15r%0ZmP;GiY=j$?Q&hTJfW%b)@ewKpcxr%rot>hP3JezWm$%YAP*ur)8?4l`&ch9 zDCjzvloJ*$4c|BhDBd(nL50GcvE;8ni5^A zK64?>Ud5>T5G}0!qy7ngYh3o#HSLzIe<1F>WMFD)mD*}bwi{4fESw&7$=NCASL9*R zVlD-d?|&LHc+4x9-Oi2Hjk6t`+(1=PVW7);_N*(Fr=q$nro1tky|!&H--qjaoOw(X z)Kvn0*`Bz0YlVAbQ-1%CkBCFoG!70uxXmBsNq73m#p2_Yoxh?V$~?BFf;#I{W8bryAgbX9c+KwsSR?UfHUd?^Noi)=ck$T#DVxvm?vgv%;%Jgq z`|Jwkbu{jfND5p^aL-tF^lMYZnFPTDRV8FuRX+^ZvPe2o73$de6tw z6fe?IMx~aPX{p`vyK9=)_`E-UE+li0{PDn;$bzq5KLPGA`7=Cji~pP(hjfz6M>DHc z=oaJj{|K7}(r*R)YZIUEn7y^Xjo`1e0$gLi;8n3NIscTIHY%Yfa}&JFs2y@ zc=zXy0kdTydsjJ+Od2eXg{YGK4~hE3ze>@6U3#SUvo%hrugwL^9Zj*L2a(>r+N1bO zvk={FZr~4nsjL`i$I~W?EJDwNfom%wDBvog0O5K|to}yKS*_@y<{JGVzzv_vDRd*P zp`sxs^t;0F%mXiAC1)vZ_GW1biqU_f>H3eY%+U;FQNj_#VAI_6$vH5FL24-N;&qLi z0ozTM#@4Ot&9WNOnlYK^U)lYMEpbErybkmT(1i$gW7{=ZJogB!6-zggsz=~JX^^Kz zM?IXCwH|Bah#m5EiTdy|OOv6^cRP8-xT*A(8|0cMyx!NZt%0AP5!x}Sp>b5`h6h6kqe;$7BnS6NLY3txwq-=9Ktz$$37iupK$NQFTq}u;~Al7pNNOXn#se(sO4+K_0u5q+3cJ|kV9npPtQqPfQgF53|_evbl3yO(U-I0yv zm>!Ts9^jLaDHH`pU6=YJjYr_ombU>sSy^9xQsQCSg$?$yo~b)(#%{+b2ZHQxpkodW z+v@qvnT%lJK&PN+2=ZMB*x%Ewaf7zdgw+wadGQMAMoxtR)>kHk|C&N@ez7?7;ewVLiGw4xBK z5cLfe6bb6pTRvG5g5eJ7|0F1`5B^;FKpTkJ3P;-t=I^=P(GtJQ%L7$^MQ}nNTR~ln zOM9*U+va&(xYCshaAhYuoNw_5Y67DxJMAxmeA&`jOz61rO-?kcbJz3_gc^XLN-uiW~A!XqlEd}Ahnat>VK*mpaK zB?c4uLY7Jje8O+rf#fXC1_C1s5a61#$Dhu}XxZEkg~;2zw`-rEH0s>5y0-pD$5tG+ zh&*9Is}k-vfw$&ctCrrnB%MoM@EwF6?uh}Ier|q7Lx-PBcs7B7&Yw}~T0zknTlrwh z9y5Z-`!OH>>%@Mv??uj7hGWe z4LZ6V9!XO(qw#OgE+t#+T#MvxjE&`PG}jgr>Kbb6l^xs}`9H^-+(cW>)eg7oshjVk zPEu+rAAEGh>fTPLvC?tYb#$5n%XXr_UFzEP|cDR)hX3h2I(wZ~v^NYN?l$+-j=%J@}~r-^Q*4?saH2;pQ#v2=$7VeQ}A}g-D<8 z3=R!A_H6v9?me;k^bIL0I%8&D>Ml`r7gY<1?-s|NPK=3y~`{As6<0 z&~F}1wjOF${^`6h)H-N#&@Sd|6o;=>JF?N%=?=IJeYZhAC{VYR zT>9hfFghUVtJVgwkfG4!8WzEnnVC!7Bl)%QmPlhog?$cO9zWY9VukQnLatJZUwRvB z`o}>fZ!irF>t6h0w`Q%8?4_55HX|+F#UPtd0DI|UYE-d=)vYG-zZR092 z@%g2tX9v?nQK04~0uPo`fsdZ=t+|(f)a8Cm$G)O&j}5dAwjFotT+~^#5_l>yGTSVi zoW}1NN#9mLmeF8Ax9mW4-2EKsYR`6u<9COtwyIZDPo!Mrl6jLMNXo zx%4wsyS5!ZOUIo-pdfO_)(3;?`bPrc-!Q_^yUorXP(4-?{^u-s(5OMXZW@YAMq)l5@A2S@|dXq3(>V=mgMn5UYK>r9jgb;RD{9R(VZWa>TGY+P1=% zjUF<2>ra_GzBsa}F1E0(6`iS{e1=XOwrXEn^<>3!dtcAe0!1q#H{Urv4sGhMjfzx_ z7)|sLjr+*eBv{|sSawoAIJMJfa3z}*<>V`+OD~;d)xC;lx#Lp;tchoxAf6S6TGI%@ zpm+ZKo^4LS{|F9SGckecs)+&HCOg#bO)@=jN9AAW`>N_@L?1YzDw3TYwAmA4?;218 z1Wt3YjcK=mE)|q-mQ9L3nJp)3pj(6jcnmpyPkl*qf;A!g8ojLyI{@W$m@LO-FhikR z?be#Pl++G#8K?|k%yl_{3$+7{^OGZ+i7fxJSG@a+caY>xB3@0*pY7Qk6(D9kbd$^t zmtt>LmH=|M=J5p)?NFnCz9ar>LkVijv{gboGvfEXOSu@GtFy-eFuXLapiBU4Jf<`I zmH^uP_&la_L4hu&*)WXvjS|}4)?6g7r%Jm9CT0hXZqb3sYK2%S1fnZNg8_ZCwnpP@ zJqiSG&zj1FhE{YYLo#{%;4VUW@b95);E-N8kS?gGn|)d(v+O~G%OZT8BuI#@_6eMf zS|+JVu8^E1{|?fKJLBf6cik3~vqhEIb zF__l-{><`Pr3m9lfr#fr_A{Ols}sL6qW7;aPlkn}@jAYL=$3BU8QV(L!|}3X??1W& z4}TWTs1~)(NG&mK;}(5T#emLRbSzIzF~19U^Mv{0|G<4ODJaUe(;uBulgQSjv5}{K znYUNyU5+kP+Mpe%aemDOyA}Jc^8c>!QRq`LBU&-zRrCIHTTiCm)Hiplz(orLaDN`6$=JmPTQ`;cHrZ)rtCTQZlKX-V{Z7+1RsdKe-LWeDDMvo#Qmf^sQl{2 z8fpfC$H@VsVK3*QEs>hkk;2W;(jaPJ4R+%q5(N{EO-TN3pa$4lY*5eCAj1-i2~9%- zV$J@59rnxm;HeYX(!UFfwHmrj#7k7WfHt))VJ{p)dtx}O!wrl!9%$X?4*V6V!W zc^YNwBd?CmyBW$#%2HMu_dWx93!-ec7Hmy!9K(lq2?XCny?uuHcX4BL79KKqbYrO6#M z!FSY_%~jkFZM^sP3Q70e=B4xUEiC{ld_Xe`y;6|mt~IrAAsYN+yylI+IZNJ}yM4P) zkIJ$O)b*)P4<4+i+dkO6nSkNpRWPB+zPA_c&U9rQe9G^jfQVAq`_a^RrN<98zD{OZ zZ5ue*YMHK|Sx(H!`#Bhc`F;~^54d2@%RiP>+}G4_?RMk2pwF>%_71C3%vCFE$Go;w zjvxi8D;qh+eIi_+6F#$5_s&LHIEh=^#$>;_sDjEkct{XFI%}KNG_WPK^-p0u;^XZ;+kF zKX5C_7_zdW^}XFG%yW~Im33uJSo5r6V^-oX?2y8nF$XPsX1ZR;R*#%J?3(>vLM=!` zuad}5{9w3APWPRQ(M2~NLz{bIW>(R>aBV@(1@HIt2aI;Fa#LJt1HX#%fPpdaW3*$fJ%dx5xKu8YdM7(e|bjc=eCnFPsd_mv_ykma->Hv@Jv&H1}E55OL zGW2|l8ofacNSU*+9_Bg>VwyTmbYp)y8Tcq@XKD#6>w#2ZgX4?EUcDt-kL{`PSI^lG zXV8vNnk(hNcp^sHuS3oLz6gVglxAy@C)?qF4$oUH4_nJDaEBf+y$Q+k>Fj2M5b>0_ z$UQz$lIL}2liW|^Wu22`9MCnNv03wG8<{NxbSWx1*_#9+L+^vmh#tO*D)tChIcaP{ z_xY36>pAycQs>VolmtfhvL!q2%|_gN!81O(qH0x#yWns3Ja*ehU5ys}8F40+F-dJ1 zyw&{Z9gp1||0aL@l*^iho(gsPtgG>}^UY%Oxcm47`YOP8WnK8ZPzRK6^m^2CrG)?N z+Q3u6MDQeyoX;|Q*xsPQhD#;S+oSCE>|I~9U9r6y8a1@wZZj;|QE+^*r@2kEMH$aK zcDNn5M%5imN+{0tY&M_y(^AtG#tOJe1sL4N9ln2m|1*i8NErdc-q5H- z69jgWlgAfK~Q`aNms5`1kaHB7VSsIQcA*fSI(+W0!1>6KE87!aJh)OHi( zd91SazD(&^8Ax+;5B%xwbq(Af$@jtHD=*S&=V0{PF1YD4M$0Jr35F=jzvl1wstR88Jkk6zYNlt-d_2Z zy*-lCHtL36DJ_yNR{{qr^R9npS(PBa`mD2caDIO5&@O@fYNq5MZ{(44+Ac3R?}4Ls zh3KNOzb?ysWazo$wWY=O6FFf*qMc0)<)QgtJ~~uRw=b+t%F<^G&|lE`T23eHNW64@VoKc z_0*ZZ@M9Kxem6wTQL1|?4)lmMT12HEXt4#>pm9xfs3zf*{4>6)IP=ymsEC<24#~x6 zBSE(;SlmZ8*IfY8irMt-0WnW(^Tt6UrrYm-#iOS%YIz3=6<{0aTX!CAv8Kl4x4d($ z^&=Fn+)?S7s@MQ{5s*ZGr^yq-S9u7!B3(L|Z9YZa3}4YBRcN2iyRdLI{COzRLmc&R z)v67_+u2|D3UT7V`tBzPU*E-Dqr3{ zeNV8i-s|my9Uj-?q8PkqB)<2B^EQrW+Hkz}Eleq4$DgHAiqAC&_qn<=M;4&S~r zXQ$Ng*UMEzy6i5gymlDu=wDke2wRSvvnBZ2pXbhuZaf{nn{OA4KkH5r%HE9Kq5m*1`fRfW#u>ayn;xxO=_1dAM|Y{(K`ldALnLPjki9cBglgWucOsIMuN0jmgRH z*liUPudRJd_f+yc-Ak`||52C2ID*NUil@?r3DXwvpk2ri70*%fXHSwgLV z_IJ{c9e{qKYLNl039EDkaN%ss%-2xMN=!)~kv_qMw-QVs#2BULPbWM=OYP%#8=9J8 z-yBi8zm>h!1W9*xabG@v0I6EvB@lh^Eu)WEfK}&oi(CzSH;=5CCfR9#Csn5S1^_S6 z%pj-VUmL7LfC}x=DqWhAwAZU7+hyd@iG}{hdGJ z4#H?JT4K;sea?UBY>v&lGvaLFdMuxx$Jix*i$%Ncm$$>QyIWjvxs&ut2M?kOWTBu|D^!9&g~=CfD1Ub+1e`b3V!$J(O`f@H>aH03TU1dJS}g)KQoR3ZT^B+X-0}Z@7tcY*3q?pEUhIx zbxiqr$0kJ6*W9t*yh&{b;;=M9(Q|!FS(&SfJE|j&$zGp&*a^umX#GxS6HlvyK{*RS zwP_)h7r)%7qHcaTA0{%In^m zpCc&bkVvtz%<#QJOgxVfXx7)~SE$e`H;^R7Tu;RcqXuYVlAa|n3sI9jYc=Y?*P6OFGZ=W1)?#46|aDaqc*1q4XTeHy~)B&3mlsq-;ihlJ6$ zoqEuKJZ7k=o~v!OvQ60|{Cu8u<(BW|o7hgN^qH9C-8_lZI@@M5^h%P`H`BmLw7F1L z4a{_mTW^?U8JY>zo)3z)fsL>3wC%nyHx-hcP_?$SH7ngjKKNjIM!U>LxZG?dY)U;W ze-m@WTq@)HX!>|$vb6{i$+C3eS4Gmki`T3)O|!}Bl727`$AbIX;L)SQ*ev#Q^RtiS!!fsi5dy0Pi! zrQ3ctvy($3h9ifd^1@@J@u$`P5I|Nz?Xjm<0J@6=`he6<@LX}A?LZ9Odh56CBkAW~ zeUHAQG_-u2$gu^v9NKld9P$Vn8EQsI(6O4g!*;8b8e+kr!EnJt=&k*$QFJAJPI?Kx zuDlEAc>4*w5ku1Xl37Kq@ucI}Eyn(H%MJG555e`a$;*3JvT@}aPfJwf(jHQVqpQ&S zIELQ6??gw&!t4a0S${g>C@aJbn!UJ03V%bFuf?d_K;dqWF*gK>(e@1$?sag{V;=m9 zzVBHx`lpNWMb^OubYSKYDzN7-foiE~Z2dzN@QhEZ1XxKFSEHX4D1q_==(SV?^>U=^ zk2P)H^TWB2*w_z|xk4xpY&&en9qBo)q(7quKsV8UL0KHOweFy<4thEuL$K~Tk7mmN z`nt+=9l3qmUnGB-pYQY$q2kF`^q~v7H9KR0}kDolbTPa+y9(?k5V$EeoyTL+R)>9Ee3`oM|pu?)Q!BBlU z7WaZ%a$a3$b>Z}ODC~h!-Ob0ggaKmHGjRi5Vi&=+c{EZcB1vEHME3l2mS@h7$456+ zS3kBm#8Zzh_K$49g3HN01~Oe9@_#ExC+=*R%J&4R`Ud51}c%5AUBk zz5j=?$?V49#?Lxk=+RBcXOApCw8sDz( z$TnG5yYm|vd_kLT3-i+dtS^<=rp>h*djnqb=|R`()^Fc_P70SwGAM02)gYFMk55*1 zHXZ4L^d1=v~~KuTSjA|Js5R!Zon5iQ&;>x)Z)S6xFe+46H2w%8IcR0*IbX_)t84knNnB zIBFYl(DHFZ0g*kHe+Wk58KMA=KCnu+QZP4U7Ex)r{!$Q7oRC>q#K1*fRtc z-IFH{pu@?)ufrh7^Oz?m%p~ zUBel1H?L6$ihWW$$#EOd2$dXV5Y}5*B(N^;=r=aar=0aiieHTgT z85#SpR6n0a<@>)nx%zUllCRW6Ld89(*`&9(7g{DiurOM)NveoL@c$;x(nCZ(rKeXenql%#L@>ku#Rfj5$^iA_9bj^DVY zw4iaDU#z92qa>SFOTJsKjv@Siw@EcM^?i?2SnZRTt*Kqvi<-*{_3J-kpkdGhwwh{!UhjnmIw!jwnPUy+9KsiCtA(A9rEmIq-MTDTXxJNV#P{^%+wy!zN7nh zik_2x&?GkcR^iJ^zen5^=JfQ(G+FI>RtH9=jAGmX&`*zO7rJ*d#m!vfcb~PeR#plU z!OG$}02yC3tIDcJCWAwRg(BKfS`tUBB!b^`f1efb=6<$YmRZ|sSF|EDzb=D#eAA{N z6CB7^=lC{q^uqyNL*`p4<>BiS#Zz}}Mfz6>yMlyM0RpH+^4B8^zJZH|e3dgMIEJW= z?=Bk>Xx9VZo?PoFy8cnlrP6R8NS_IL@u1xj>~B&TpFQgt{r72EYN7?Q<3sKQ0A_De znync0i2VLp9U&1iKMe*LbUvhOYNCaQ#-sA@9Tu{(hrhjtM$`>Dlqi8`p&#$@zsqw5 zXNqRvMsL%c<>$`g|6NB(OEVR_}!z$p$a&99acznB1q z_4~Gv?Mn~w_{k9wCO9zB&U{%9ijZfipasfcG_y1cY*drIzJvYXi$j6+>Tny>jMa%Q z!Sr4N6Y}pp`y49x>LO|WEj%7udirazD@;tTr@H+7=d30!i9zTl(?Bn{1){pkhNiJD z^1&9JmjipweDPn~)S;!jM-Tpu8u^q8Lwy4wOeH60ycttFk1~gZ3tfAD#XW|+sUbeS zB3Wj*aQx~fS3DRAOe4mS4z4agFXDLH<<#_gU2UTFZ?qe&Dz9AW;8>DjC?}@?-T8Q> zVwyWHXN)us%RIE+w9k^PIywIwef=86uITc|5keDf8@L`*p#Jh-7d?L-5-of9YqbgI zk;l(Cfk&xm`hR{QGBMly*(%=k!JY)RE{&V7jHV1vT%nS1GW%OvV|6(i*P8Sszaa4o zuGH1@@}Z<`_?B@8tywHB8h@+Xckk|ngLpaW3gOFHoe%2KlU-3;Pxq^RXntNV=ksH@ zi!2rJA&6OF|9AM}aMBi;+JpsyV7Bh^2ghHc{yQD${P6Rly^dAEA=h68fSQIU)M9o* zoG(6;MNzeX3{U$BJ{Ue!VvJ!=*N72JGLUi@VkG~LA3dkQ5E#}E3a*>k`YFFhwr&an zJHhxfLrq>T&Y0MvMJF*1P@7+AV3}{^tINt_4z#OvzIa{BMo{1qc9MMYFRPZ3jxQFp zTihpBC6wNW+bU)i-@}tQrF;Z5`8-h|QXl@|`QP6`liZUPD?J$@N9Aw?sV={(NJ9u3 z5e-F~f=UTrx$<-H@b#7NBBj6=tp?Ao&ueyi_?R9$u_?#hw8%u+9?OA+c>X~$>WClU zn_g#Sf$a9XVy}*{u*ZL=a=^jCxWGk(NT>1e(8d?1!h;|>dq+0MlcUlcR4NQ|@2g_3 z5{fET&zt7OCM^Uq#I@1QMqK-H2}8^DP=8yo+|6Tzq{uNk!$WuZKk)9(FZqv?<#**4 z%RJceGR5K*Ki4qkJ3F~Q_RaxjMNoO0@V*&R78Nla3{5H_jlrm;HuoHqvpRG>;1mh$~n6RHO8EpQtj1h*KU;54F zH5*!+;*0`OenW7ESMA67@7kcp_%{y9-C$LkgB2rSj@X{j_$~2tdRPsA-qBveH6+k# zcg}h!-AlUh(NCDiZwn0=EoZ!r7%I}IujqE6p$T*aynGTl#4ORm$GhcQhyg$VegA)+ z(wAtB{=4sgoG{%DK*r?HuE|KbJFf!h)YS|BoeLyqh((o+iu%~4w!_~+>99;!Z1{d z-Yw5rv64)*yHg1ND7qr8j6hos45*B*!JMj{hEbF&#u}mKDC@D4)^G3mIr0xN+na*m zEo};dHN)D04$ZO7;$ync8a#+$;x4Y~KLEg{!O`{mtt_+bHPsJY1?8Ddzp5AjPnwPEYO#{|Yty#8Y-+5O$C zJUn!PJ(C;GPC^5TboSi2?%$DTzp=28fX$DZb%aC+a%QVvRx5vX+jARa!S%`6*;`12 zyv?1NTOWia5MFT;-{BdI!9TP}8m1H8r) z=g+gHr>8&fhHw|&cm#3XHAcG<6EdZ=gS+IGC;?>}ujV)HEaC9YfK zF{5ytE^hi30&T^dTent0P_22O{VqR5)fP(sq5FTo4DO&XYH!3FEV_;Q@04x|68?KO z5mTgy-ZZ(Lif7K0+Imm01EX@ZNKoze|HGL6=Wy?DphfS0UtRv+rAnLgtHkU7_mhw1 zpL)WwPG10hfc zn$Zgm6}d<(u8>oLjoDpxob>jvd4b+G3ljV7z<&7MC+Q0FyK*WiFnQ(gf6@erx#3lg;DbQV zF;FG!bN;2X*Pt-B!6|z2$%-b<9M_}TZwGTrsd>j6>+k-kDs-F3EcTkW0DPqx9r&4> zb}3xVuO3D2Q%FESi4fw;1QtP~000BAu*%m(%S-bH&Ftbf&k-3M9!^M1sVA91@y5qf zx-YLkGI46~tFi|I%FK;23uoD-+%pd$q)$#ssebJ9iX-W$4HqcA51G|+gR(lZR(g5C zaTsE8O>nN0d;9t#5_IM2=lahMCu0IY8e9?|f!GDHy z*gm=`v6NDXCzEg) znanrH$`Bw*`ll5|{hurJDtHi#vc5hMb zU#rED<0_$jR|^NkNM?3+_VO&kwUJ4Mn6xyt%BreHLgXj45=-fVb>kZ~=}}9FZ+ZJ1 z@q||NcG58DF)b)oh;eiQpPX)Uta-XLbL+ll*)SB=bR8Wtnw@~Y$$IpNWjN?0r#Z#H zEO}xg`PJjcTapba-mOgHZMEKP55t$c4?)hGwoafv(kP(#!+yxK4MK>SdNz5y-@C_h zouq95NUYs``pI{vD2}mQdLJ*3f!ASTq*lm4Q$@P%C!}q#E-)g@`(!?Hi086^$V^)t z#pW%uF7;#GPD-^KfOM_yLr$s8<;BirYs%g&0ebkOxJrTyF&PZ3x=>AxHGk#|ihfrq zL$NarHo094$=ay~!)1jpU&>ZiS2qEu5jj1CjDZ=@whh7sztX2GFLymLf96znrs2oX zh4%^FTWtH+*7+q{v$!xb)5#1gD=Vw%>D>obmpU}=Sbcm(= zo*Z?um^WFRpDDYwhG2=}Vr%p0*Kult+EU=7M?Fp`tQNeo?JWL=Wi@XiX%i7;%fY*Q zHxC~l=Lp>*_gB=i#Bs8j4jE0m6 z{h1hkAl`B8(xsiB0ywC|?9bSrOI>VfX)$Q7tf&}6gA6g5mdz=*D|hO*l_f!Z=`2%*@Q*LHWM++qb*7*KPfQRW_4f4R!f@_{Qbo)-2LB zELCP8+z@vw%sVk}-fZ7u5_V~6VRpDCN9VBVHK#JByQz-!i2-_r5rT41(~{$OGs`6) zAYc#{7aL1QLQ0CM-O%2Cf54HxZ8Yt(I~%!d{uJ{(ecG+1rRB>Z;d~L)&;QN>PP^By z>E~xg%aAi9G>aF0Y^9hjOFz`>?{0!Vm7Z-8>qhB6-0o9@y%S5WZJsin20n{UI6ZiF z?ONFod!_{Dv#7X+*rv$&-0bT!C6mb)@$6ulRDeD$p+os9D=TZ!S<>dZkBA1~)!2T7 zm}8-EKQV7B88>`2JuB?+Sut+As;%s_wM6w^c4^PUa17Es zW=5umN`LHMd0fKA9h*VMX3R6?pLjS%PWxaE6yIFV39m zvSF~bmLdyfk5lqb5f&MKwJ=drW(a9EPY8X1t9sclIGFkLlYJ}Dt!+)&Ks>CVhzKKz zcn=WB=pKWq`|()M<~_tc3lq9r+vt`Jh?@Y;B|)!*wf&F^%w1igwg?L=DulJyfu?6vX>!@Q#4lPSVX<`z*E#JfNh@6S-a!U&LKu8o&wQ!VUN` z<3`R^*Ho6P@tmX)Bj z!|Za=Z{XvNkf1w^*Cyfhu3K9h1emGGD&c%Mm3=m@W*I~fMbLm}+LmxlBjD6!eHxNgg!TDUT z$-0w&@U@$e+}-d>;Q<}#^fR$FcQ$Z*z&ZTDGPP(6NPs16I0Ez>91e824+&jxZZR1y zOMdfabxxgTh-KnxpJhKZEjpkEh=q18ywuO~eW@95Vn>n!4f@#C+>2ab$&u{}VUrUSgD}LwBx)*amnEhL(`Q}5m z9sY#6UZ>X#$KcS=e$%Qa(;*!?M{IzZCK%Gy)YPs+P-0bkzLVv02i%FE;NW`m$SSa6 zuE;fPXl!bpA{_NN)hL+45Bpap9U)91PC`lTyLQte>rb-ec?xk2;q<+K$Zg!%8v9V1 zl6pm{=Aq0|?qYg`gp0yd{Ug(dRDdnBg!S{yn)O{GI;q%8KE-ZFyZuEoFO~RaEk_XR zCRx8gOOMDPH>Z_xCPsj;sx?~vekP7=%V(duWBX^bJ!Wi#1O=4|2e=j92)q{s|7|bN zL`k;~rNf7jhf3*JxZ+}Gy6x9Zdg5SsatKA7s7~5hOQ@9`B(*}%pfJT_C^-j>9BiJ& z=5vOfg!5$6o^MW=Wh(9DZEH(f?BSECL#{gui?X`kiO>@ojtNm2&;Rv_-F8wppy*#$ z#d|8zFv~FR^=o;Y8@X7RVOfS-)|H)ZNRqIL{Lyvc{pap2sc=<3zBcFJk&ND^t*w0r zL>ntfOP~;kenw`|EpPBx5P5ClkZf#jHYyL?2`^-$fom7N6*Zg5&}T8LIYG74bAl<` zwle@q6J2mvKLNOVJ4EH^n)U*hvKXb1)>PWfJz+ zfDa0B_%1xe@zIW=F=)56TYi0%l$10b%1+`Guj!tKSJY*=Jhv5eu9A~O%^Mjv?0ARm z=S!TcU&K2nbx(zs4O>zYGgiEV&H@qOdBQ>%giYk;J*2{I2iMWuYA}n)i@!ujY$XmHL>mF-4Y8M6(oI}{+~fOQCbP7Ovd}9+O)cwg znq3}f&9dC@K4}J*e);#a(e8VMu_dSrqN91pyLbN(r_3_B&x|}EJR%p$JjY1!yiMu~ z!V|5@A^ZXmZoen$A;j`0N{IUr`gy`Y8i)o*FXpzR)b@wkx#UY&`zHt|nrD1Bas^=8 zl62DO;v(Jmq~!H&zu?vzs#ENFgs3Y(+-6il>U5$s&}+?XR*}!W!)iotp@G`yp$Q( ziRE0@ty?X@tbe@;(X&aj*BIcF35_F-A#-13Wib&Xi+V^OD#Ih$XL(5% z1ubWg3YZ~O`*pR*Ovdz?>9GC%xc~NRiOT-pPmu8c%z4m?ZM6Ikk3OxuMf=kKm0bBZ z6%w=VZ$^;#gZkSqH+7;`)9XOL=#B}nhrGl^H6jyU3nf90MxCwcf_5ppA`J-3&)Qg>Qbmt?kR`mywVBwOq_kwGiGu4`Ka$7~% zQOmNhu(0bX--Hyyp0Sb;$%u^@HU7h^aBSm_c z1Qyx9>2~{9MwZk{5z6|o%G;rW+CL^6R8K;?#2&}B(fb?pgy3DO&!Xkzua?d5rQa>i zbUq~F`*}fO;ZG2Ixc2Lp6(M2^heF=m+`M+twKr6IpLu=E^gr<+=>WIQ#}gZ+8VB;V zFW3WN$JxF#17C%(ZA1(Z2=;~jMxVYz)^;lY3j`l;iPlIG_c^AewS^#U#x2fa-{ds zG=-0kkD!HlkBR!ls>K2C2u>*B_c@Y|p@3WskCU*HuHzMaV<;dSq2_-W>5(+=#kpzf z;_~!}P|SdMK=+?A-6g-(y6V?6?!Vu`TeYG7%+Nl2vGki4oi|-#Gj4pt&;QTM%L=Ou zhONn7o5&x1{z1TYHrAaSPv$6h7cVYb`*^SQ**)v6Z{1neKPxSh z+EX&&ae8-%F2%bD^3SiZeaeTX60)T3TI_N0EV(iwlA$ zvqp+4iS&CWoz;CMiA1`M!$%s&i*j&}-Qk`wSC;#wRRv*m*RQ2q$cs)%5uWIHeU_hG z(4Im?x5{8+d@!Ilk1GbQCcJue;{C;aS5k|MrHg!)4LBe2y($I~dR9&jVst%Ss6&$n z5RL>JPfL1CCvi%9zQKx&kB<-kaiL1LuPFUjt*V zJDbB=W8@I>dPhr{9Y^f=^+dCgtr zIN~6O7^&i!IBN9I=TAqa&y97T?6|T5x85krOd`G8hErJ?FZ;{vS}JJGmK-A|`xnBshonG~%gotKEh`k!Cf=9`d^(Adxr;#ioM zcfQL`2Usgx|JtW99B{`c%J1I^-KCbx^C1Z$Om4&%E z8bbmZ7#W==f2tkX5uhezxOYnfJoT!AU~C1lFqMHs+O`==Vaik7rtR#B6D4o@NTh2N!#c$j z%0-OiH)VEmAf)e0c$}G70&KLYf05^`IwCNGWR|tJCl_WNr>3Vrn+jG5t-@y3R#OY&2~w1u zTS>Y$c3+!{+pc|)m4yxzP2_RH$OCZi#@&lf?G73I^=q%rT`@_?r15bF + document.addEventListener('DOMContentLoaded', function() { + showFilteredRecipes('all'); + + function showFilteredRecipes(filter) { + const recipes = document.querySelectorAll('.sphx-glr-thumbcontainer'); + recipes.forEach(recipe => { + if (filter === 'all') { + recipe.style.display = 'inline-block'; + } else if (recipe.classList.contains(filter)) { + recipe.style.display = 'inline-block'; + } else { + recipe.style.display = 'none'; + } + }); + } + + const filterButtons = document.querySelectorAll('.filter-menu button'); + filterButtons.forEach(button => { + button.addEventListener('click', () => { + const filter = button.getAttribute('data-filter'); + showFilteredRecipes(filter); + }); + }); + }); + + +.. raw:: html + + + +
+ + + + + + + + + + + + + +
+ +.. raw:: html + + + + +.. raw:: html + +
+ +.. thumbnail-parent-div-open + +.. raw:: html + +
+ +.. only:: html + + .. image:: /recipes/images/thumb/sphx_glr_plot_22_recipe_thumb.png + :alt: + + :doc:`/recipes/plot_22_recipe` + +.. raw:: html + +
Plotting a wind rose as a scatter plot
+
+ + +.. thumbnail-parent-div-close + +.. raw:: html + +
+ + +.. toctree:: + :hidden: + + /recipes/plot_22_recipe + + +.. only:: html + + .. container:: sphx-glr-footer sphx-glr-footer-gallery + + .. container:: sphx-glr-download sphx-glr-download-python + + :download:`Download all examples in Python source code: recipes_python.zip ` + + .. container:: sphx-glr-download sphx-glr-download-jupyter + + :download:`Download all examples in Jupyter notebooks: recipes_jupyter.zip ` + + +.. only:: html + + .. rst-class:: sphx-glr-signature + + `Gallery generated by Sphinx-Gallery `_ diff --git a/recipes-docs/source/recipes-source/index.rst b/recipes-docs/source/recipes-source/index.rst new file mode 100644 index 0000000000..646c716435 --- /dev/null +++ b/recipes-docs/source/recipes-source/index.rst @@ -0,0 +1,7 @@ +cf-python recipes +================= + +.. toctree:: + :maxdepth: 2 + + recipes/gallery diff --git a/recipes-docs/source/recipes-source/plot_01_recipe.py b/recipes-docs/source/recipes-source/plot_01_recipe.py new file mode 100644 index 0000000000..93b9dde417 --- /dev/null +++ b/recipes-docs/source/recipes-source/plot_01_recipe.py @@ -0,0 +1,77 @@ +""" +Calculating global mean temperature timeseries +============================================== + +In this recipe we will calculate and plot monthly and annual global mean temperature timeseries. +""" + +# %% +# 1. Import cf-python and cf-plot: + +import cfplot as cfp + +import cf + +# %% +# 2. Read the field constructs: + +f = cf.read("~/recipes/cru_ts4.06.1901.2021.tmp.dat.nc") +print(f) + +# %% +# 3. Select near surface temperature by index and look at its contents: + +temp = f[1] +print(temp) + +# %% +# 4. Select latitude and longitude dimensions by identities, with two different techniques: + +lon = temp.coordinate("long_name=longitude") +lat = temp.coordinate("Y") + +# %% +# 5. Print the description of near surface temperature using the dump method to show properties of all constructs: + +temp.dump() + +# %% +# 6. Latitude and longitude dimension coordinate cell bounds are absent, which are created and set: + +a = lat.create_bounds() +lat.set_bounds(a) +lat.dump() + +# %% + +b = lon.create_bounds() +lon.set_bounds(b) +lon.dump() + +# %% + +print(b.array) + +# %% +# 7. Time dimension coordinate cell bounds are similarly created and set for cell sizes of one calendar month: + +time = temp.coordinate("long_name=time") +c = time.create_bounds(cellsize=cf.M()) +time.set_bounds(c) +time.dump() + +# %% +# 8. Calculate and plot the area weighted mean surface temperature for each time: + +global_avg = temp.collapse("area: mean", weights=True) +cfp.lineplot(global_avg, color="red", title="Global mean surface temperature") + +# %% +# 9. Calculate and plot the annual global mean surface temperature: + +annual_global_avg = global_avg.collapse("T: mean", group=cf.Y()) +cfp.lineplot( + annual_global_avg, + color="red", + title="Annual global mean surface temperature", +) diff --git a/recipes-docs/source/recipes-source/plot_02_recipe.py b/recipes-docs/source/recipes-source/plot_02_recipe.py new file mode 100644 index 0000000000..902c796788 --- /dev/null +++ b/recipes-docs/source/recipes-source/plot_02_recipe.py @@ -0,0 +1,96 @@ +""" +Calculating and plotting the global average temperature anomalies +================================================================= + +In this recipe we will calculate and plot the global average temperature anomalies. +""" + +# %% +# 1. Import cf-python and cf-plot: + +import cfplot as cfp + +import cf + +# %% +# 2. Read the field constructs: + +f = cf.read("~/recipes/cru_ts4.06.1901.2021.tmp.dat.nc") +print(f) + +# %% +# 3. Select near surface temperature by index and look at its contents: + +temp = f[1] +print(temp) + +# %% +# 4. Select latitude and longitude dimensions by identities, with two different techniques: + +lon = temp.coordinate("long_name=longitude") +lat = temp.coordinate("Y") + +# %% +# 5. Print the description of near surface temperature to show properties of all constructs: + +temp.dump() + +# %% +# 6. Latitude and longitude dimension coordinate cell bounds are absent, which are created and set: + +a = lat.create_bounds() +lat.set_bounds(a) +lat.dump() + +# %% + +b = lon.create_bounds() +lon.set_bounds(b) +lon.dump() + +# %% + +print(b.array) + +# %% +# 7. Time dimension coordinate cell bounds are similarly created and set for cell sizes of one calendar month: + +time = temp.coordinate("long_name=time") +c = time.create_bounds(cellsize=cf.M()) +time.set_bounds(c) +time.dump() + +# %% +# 8. Calculate the area weighted mean surface temperature for each time using the collapse method: + +global_avg = temp.collapse("area: mean", weights=True) + +# %% +# 9. Calculate the annual global mean surface temperature: + +annual_global_avg = global_avg.collapse("T: mean", group=cf.Y()) + +# %% +# 10. The temperature values are averaged for the climatological period of 1961-1990 by defining a subspace within these years using `cf.wi` query instance over subspace and doing a statistical collapse with the collapse method: + +annual_global_avg_61_90 = annual_global_avg.subspace( + T=cf.year(cf.wi(1961, 1990)) +) +print(annual_global_avg_61_90) + +# %% + +temp_clim = annual_global_avg_61_90.collapse("T: mean") +print(temp_clim) + +# %% +# 11. The temperature anomaly is then calculated by subtracting these climatological temperature values from the annual global average temperatures and plotted: + +temp_anomaly = annual_global_avg - temp_clim +cfp.lineplot( + temp_anomaly, + color="red", + title="Global Average Temperature Anomaly (1901-2021)", + ylabel="1961-1990 climatology difference ", + yunits="degree Celcius", +) diff --git a/recipes-docs/source/recipes-source/plot_03_recipe.py b/recipes-docs/source/recipes-source/plot_03_recipe.py new file mode 100644 index 0000000000..25633b1c54 --- /dev/null +++ b/recipes-docs/source/recipes-source/plot_03_recipe.py @@ -0,0 +1,35 @@ +""" +Plotting global mean temperatures spatially +=========================================== + +In this recipe, we will plot the global mean temperature spatially. +""" + +# %% +# 1. Import cf-python and cf-plot: + +import cfplot as cfp + +import cf + +# %% +# 2. Read the field constructs: + +f = cf.read("~/recipes/cru_ts4.06.1901.2021.tmp.dat.nc") +print(f) + +# %% +# 3. Select near surface temperature by index and look at its contents: + +temp = f[1] +print(temp) + +# %% +# 4. Average the monthly mean surface temperature values by the time axis using the collapse method: + +global_avg = temp.collapse("mean", axes="long_name=time") + +# %% +# 5. Plot the global mean surface temperatures: + +cfp.con(global_avg, lines=False, title="Global mean surface temperature") diff --git a/recipes-docs/source/recipes-source/plot_04_recipe.py b/recipes-docs/source/recipes-source/plot_04_recipe.py new file mode 100644 index 0000000000..a50296bb9e --- /dev/null +++ b/recipes-docs/source/recipes-source/plot_04_recipe.py @@ -0,0 +1,40 @@ +""" +Comparing two datasets with different resolutions using regridding +================================================================== + +In this recipe, we will regrid two different datasets with different resolutions. An example use case could be one where the observational dataset with a higher resolution needs to be regridded to that of the model dataset so that they can be compared with each other. +""" + +# %% +# 1. Import cf-python: + +import cf + +# %% +# 2. Read the field constructs: + +obs = cf.read("~/recipes/cru_ts4.06.1901.2021.tmp.dat.nc", dask_chunks=None) +print(obs) + +# %% + +model = cf.read( + "~/recipes/tas_Amon_HadGEM3-GC3-1_hist-1p0_r3i1p1f2_gn_185001-201412.nc" +) +print(model) + +# %% +# 3. Select observation and model temperature fields by identity and index respectively, and look at their contents: + +obs_temp = obs.select_field("long_name=near-surface temperature") +print(obs_temp) + +# %% + +model_temp = model[0] +print(model_temp) + +# %% +# 4. Regrid observational data to that of the model data and create a new low resolution observational data using bilinear interpolation: +obs_temp_regrid = obs_temp.regrids(model_temp, method="linear") +print(obs_temp_regrid) diff --git a/recipes-docs/source/recipes-source/plot_05_recipe.py b/recipes-docs/source/recipes-source/plot_05_recipe.py new file mode 100644 index 0000000000..e40f1c23ad --- /dev/null +++ b/recipes-docs/source/recipes-source/plot_05_recipe.py @@ -0,0 +1,64 @@ +""" +Plotting wind vectors overlaid on precipitation data +==================================================== + +In this recipe we will plot wind vectors, derived from northward and eastward wind components, over precipitation data. +""" + +# %% +# 1. Import cf-python and cf-plot: + +import cfplot as cfp + +import cf + +# %% +# 2. Read the field constructs: + +f1 = cf.read("~/recipes/northward.nc") +print(f1) + +# %% + +f2 = cf.read("~/recipes/eastward.nc") +print(f2) + +# %% + +f3 = cf.read("~/recipes/monthly_precipitation.nc") +print(f3) + +# %% +# 3. Select wind vectors and precipitation data by index and look at their contents: +v = f1[0] +print(v) + +# %% + +u = f2[0] +print(u) + +# %% + +pre = f3[0] +print(pre) + +# %% +# 4. Plot the wind vectors on top of precipitation data for June 1995 by creating a subspace with a date-time object and using `cfplot.con `_. Here `cfplot.gopen `_ is used to define the parts of the plot area, which is closed by `cfplot.gclose `_; `cfplot.cscale `_ is used to choose one of the colour maps amongst many available; `cfplot.levs `_ is used to set the contour levels for precipitation data; and `cfplot.vect `_ is used to plot the wind vectors for June 1995: +june_95 = cf.year(1995) & cf.month(6) +cfp.gopen() +cfp.cscale("precip4_11lev") +cfp.levs(step=100) +cfp.con( + pre.subspace(T=june_95), + lines=False, + title="June 1995 monthly global precipitation", +) +cfp.vect( + u=u.subspace(T=june_95), + v=v.subspace(T=june_95), + key_length=10, + scale=35, + stride=5, +) +cfp.gclose() diff --git a/recipes-docs/source/recipes-source/plot_06_recipe.py b/recipes-docs/source/recipes-source/plot_06_recipe.py new file mode 100644 index 0000000000..91ca4ab100 --- /dev/null +++ b/recipes-docs/source/recipes-source/plot_06_recipe.py @@ -0,0 +1,65 @@ +""" +Converting from rotated latitude-longitude to regular latitude-longitude +======================================================================== + +In this recipe, we will be regridding from a rotated latitude-longitude source domain to a regular latitude-longitude destination domain. +""" + +# %% +# 1. Import cf-python, cf-plot and numpy: + +import cfplot as cfp + +import cf + +# %% +# 2. Read the field constructs using read function: + +f = cf.read("~/recipes/au952a.pd20510414.pp") +print(f) + +# %% +# 3. Select the field by index and print its description to show properties of all constructs: + +gust = f[0] +gust.dump() + +# %% +# 4. Access the time coordinate of the gust field and retrieve the datetime values of the time coordinate: + +print(gust.coordinate("time").datetime_array) + +# %% +# 5. Create a new instance of the `cf.dt` class with a specified year, month, day, hour, minute, second and microsecond. Then store the result in the variable ``test``: +test = cf.dt(2051, 4, 14, 1, 30, 0, 0) +print(test) + +# %% +# 6. Plot the wind gust by creating a subspace for the specified variable ``test`` using `cfplot.con `_. Here `cfplot.mapset `_ is used to set the mapping parameters like setting the map resolution to 50m: +cfp.mapset(resolution="50m") +cfp.con(gust.subspace(T=test), lines=False) + +# %% +# 7. To see the rotated pole data on the native grid, the above steps are repeated and projection is set to rotated in `cfplot.mapset `_: +cfp.mapset(resolution="50m", proj="rotated") +cfp.con(gust.subspace(T=test), lines=False) + +# %% +# 8. Create dimension coordinates for the destination grid with the latitude and +# longitude values for Europe. `cf.Domain.create_regular +# `_ +# method is used to +# create a regular grid with longitudes and latitudes. Spherical regridding is +# then performed on the gust variable by passing the target domain as argument. +# The method also takes an argument ``'linear'`` which specifies the type of +# regridding method to use. The description of the ``regridded_data`` is finally +# printed to show properties of all its constructs: + +target_domain = cf.Domain.create_regular((-25, 45, 10), (32, 72, 10)) +regridded_data = gust.regrids(target_domain, "linear") +regridded_data.dump() + +# %% +# 9. Step 6 is similarly repeated for the ``regridded_data`` to plot the wind gust on a regular latitude-longitude domain: +cfp.mapset(resolution="50m") +cfp.con(regridded_data.subspace(T=test), lines=False) diff --git a/recipes-docs/source/recipes-source/plot_07_recipe.py b/recipes-docs/source/recipes-source/plot_07_recipe.py new file mode 100644 index 0000000000..02908e83ee --- /dev/null +++ b/recipes-docs/source/recipes-source/plot_07_recipe.py @@ -0,0 +1,65 @@ +""" +Plotting members of a model ensemble +==================================== + +In this recipe, we will plot the members of a model ensemble. +""" + +# %% +# 1. Import cf-python and cf-plot: + +import cfplot as cfp + +import cf + +# %% +# 2. Read the field constructs using read function and store it in the variable ``f``. The * in the filename is a wildcard character which means the function reads all files in the directory that match the specified pattern. [0:5] selects the first five elements of the resulting list: + +f = cf.read("~/recipes/realization/PRMSL.1941_mem*.nc")[0:5] +print(f) + +# %% +# 3. The description of one of the fields from the list shows ``'realization'`` as a property by which the members of the model ensemble are labelled: + +f[1].dump() + +# %% +# 4. An ensemble of the members is then created by aggregating the data in ``f`` along a new ``'realization'`` axis using the cf.aggregate function, and storing the result in the variable ``ensemble``. ``'relaxed_identities=True'`` allows for missing coordinate identities to be inferred. [0] selects the first element of the resulting list. ``id%realization`` now shows as an auxiliary coordinate for the ensemble: + +ensemble = cf.aggregate( + f, dimension=("realization",), relaxed_identities=True +)[0] +print(ensemble) + + +# %% +# 5. To see the constructs for the ensemble, print the *constructs* attribute: + +print(ensemble.constructs) + +# %% +# 6. Loop over the realizations in the ensemble using the *range* function and the *domain_axis* to determine the size of the realization dimension. For each realization, extract a subspace of the ensemble using the *subspace* method and the ``'id%realization'`` keyword argument along a specific latitude and longitude and plot the realizations from the 4D field using `cfplot.lineplot `_. +# A moving average of the ensemble along the time axis, with a window size of 90 (i.e. an approximately 3-month moving average) is calculated using the *moving_window* method. The ``mode='nearest'`` parameter is used to specify how to pad the data outside of the time range. The *squeeze* method removes any dimensions of size 1 from the field to produce a 2D field: + +cfp.gopen() + +for realization in range(1, ensemble.domain_axis("id%realization").size + 1): + cfp.lineplot( + ensemble.subspace( + **{"id%realization": realization}, latitude=[0], longitude=[0] + ).squeeze(), + label=f"Member {realization}", + linewidth=1.0, + ) + +cfp.lineplot( + ensemble.moving_window( + method="mean", window_size=90, axis="T", mode="nearest" + )[0, :, 0, 0].squeeze(), + label="Ensemble mean", + linewidth=2.0, + color="black", + title="Model Ensemble Pressure", +) + +cfp.gclose() diff --git a/recipes-docs/source/recipes-source/plot_08_recipe.py b/recipes-docs/source/recipes-source/plot_08_recipe.py new file mode 100644 index 0000000000..63427f62a7 --- /dev/null +++ b/recipes-docs/source/recipes-source/plot_08_recipe.py @@ -0,0 +1,143 @@ +""" +Plotting statistically significant temperature trends with stippling +==================================================================== + +In this recipe, we will analyse and plot temperature trends from the HadCRUT.5.0.1.0 dataset for two different time periods. The plotted maps also include stippling, which is used to highlight areas where the temperature trends are statistically significant. +""" + +# %% +# 1. Import cf-python, cf-plot, numpy and scipy.stats: + +import cfplot as cfp +import cf + +import numpy as np +import scipy.stats as stats + + +# %% +# 2. Three functions are defined: + +# %% +# * ``linear_trend(data, time_axis)``: This function calculates the linear regression slope and p-value for the input data along the time axis. It takes two arguments: ``'data'``, which represents the temperature anomalies or any other data you want to analyse, and ``'time_axis'``, which represents the corresponding time points for the data. The function uses the `stats.linregress `_ method from the `scipy.stats `_ library to calculate the slope and p-value of the linear regression. It returns these two values as a tuple: + + +def linear_trend(data, time_axis): + slope, _, _, p_value, _ = stats.linregress(time_axis, data) + return slope, p_value + + +# %% +# * ``create_trend_stipple_obj(temp_data, input_data)``: This function creates a new object with the input data provided and *collapses* the time dimension by taking the mean. It takes two arguments: ``'temp_data'``, which represents the temperature data object, and ``'input_data'``, which is the data to be set in the new object. The function creates a copy of the ``'temp_data'`` object by selecting the first element using index and *squeezes* it to remove the size 1 axis. It then sets the input data with the ``'Y'`` (latitude) and ``'X'`` (longitude) axes, and then *collapses* the time dimension using the ``"T: mean"`` operation: + + +def create_trend_stipple_obj(temp_data, input_data): + trend_stipple_obj = temp_data[0].squeeze() + trend_stipple_obj.set_data(input_data, axes=["Y", "X"]) + return trend_stipple_obj + + +# %% +# * ``process_subsets(subset_mask)``: This function processes the subsets of data by applying the ``linear_trend`` function along a specified axis. It takes one argument, ``'subset_mask'``, which is a boolean mask representing the time points to be considered in the analysis. The function first extracts the masked subset of data and then applies the ``linear_trend`` function along the time axis (axis 0) using the `numpy.ma.apply_along_axis `_ function. The result is an array containing the slope and p-value for each grid point in the dataset: + + +def process_subsets(subset_mask): + subset_data = masked_data[subset_mask, :, :] + return np.ma.apply_along_axis( + linear_trend, 0, subset_data, time_axis[subset_mask] + ) + + +# %% +# 3. Read the field constructs: + +temperature_data = cf.read( + "~/recipes/HadCRUT.5.0.1.0.analysis.anomalies.ensemble_mean.nc" +)[0] +print(temperature_data) + +# %% +# 4. Calculate the annual mean temperature anomalies. The ``'weights=True'`` argument is used take the varying lengths of months into account which ensures that the calculated mean is more accurate. A masked array is created for the annual mean temperature anomalies, masking any invalid values: + +annual_temperature = temperature_data.collapse( + "T: mean", weights=True, group=cf.Y() +) +time_axis = annual_temperature.coordinate("T").year.array +masked_data = np.ma.masked_invalid(annual_temperature.array) + +# %% +# 5. Define two time periods for analysis: 1850-2020 and 1980-2020, along with a significance level (alpha) of 0.05: + +time_periods = [(1850, 2020, "sub_1850_2020"), (1980, 2020, "sub_1980_2020")] +alpha = 0.05 +results = {} + +# %% +# 6. Loop through the time periods, processing the subsets, calculating trend p-values, and creating stipple objects. For each time period, the script calculates the trends and p-values using the ``process_subsets`` function. If the p-value is less than the significance level (alpha = 0.05), a stippling mask is created. The script then creates a new object for the trend and stippling mask using the ``create_trend_stipple_obj`` function: + +for start, end, prefix in time_periods: + subset_mask = (time_axis >= start) & (time_axis <= end) + subset_trend_pvalue = process_subsets(subset_mask) + results[prefix + "_trend_pvalue"] = subset_trend_pvalue + results[prefix + "_stipple"] = subset_trend_pvalue[1] < alpha + results[prefix + "_trend"] = create_trend_stipple_obj( + temperature_data, subset_trend_pvalue[0] + ) + results[prefix + "_stipple_obj"] = create_trend_stipple_obj( + temperature_data, results[prefix + "_stipple"] + ) + +# %% +# 7. Create two plots - one for the 1850-2020 time period and another for the 1980-2020 time period using `cfplot.con `_. +# The results are multiplied by 10 so that each plot displays the temperature trend in K/decade with stippling to indicate areas where the trend is statistically significant (p-value < 0.05). +# Here `cfplot.gopen `_ is used to define the parts of the plot area with two rows and one column, and setting the bottom margin to 0.2. +# It is closed by `cfplot.gclose `_; +# `cfplot.gpos `_ is used to set the plotting position of both the plots; +# `cfplot.mapset `_ is used to set the map projection to Robinson; +# `cfplot.cscale `_ is used to choose one of the colour maps amongst many available; +# `cfplot.levs `_ is used to set the contour levels; +# and `cfplot.stipple `_ is used to add stippling to show statistically significant areas: + +cfp.gopen(rows=2, columns=1, bottom=0.2) + +cfp.gpos(1) +cfp.mapset(proj="robin") +cfp.cscale("temp_19lev") +cfp.levs(min=-1, max=1, step=0.1) +cfp.con( + results["sub_1850_2020_trend"] * 10, + lines=False, + colorbar=None, + title="Temperature Trend 1850-2020", +) +cfp.stipple( + results["sub_1850_2020_stipple_obj"], + min=1, + max=1, + size=5, + color="k", + marker=".", +) + +cfp.gpos(2) +cfp.mapset(proj="robin") +cfp.cscale("temp_19lev") +cfp.levs(min=-1, max=1, step=0.1) +cfp.con( + results["sub_1980_2020_trend"] * 10, + lines=False, + title="Temperature Trend 1980-2020", + colorbar_position=[0.1, 0.1, 0.8, 0.02], + colorbar_orientation="horizontal", + colorbar_title="K/decade", +) +cfp.stipple( + results["sub_1980_2020_stipple_obj"], + min=1, + max=1, + size=5, + color="k", + marker=".", +) + +cfp.gclose() diff --git a/recipes-docs/source/recipes-source/plot_09_recipe.py b/recipes-docs/source/recipes-source/plot_09_recipe.py new file mode 100644 index 0000000000..42fe49cbd9 --- /dev/null +++ b/recipes-docs/source/recipes-source/plot_09_recipe.py @@ -0,0 +1,82 @@ +""" +Plotting a joint histogram +========================== + +In this recipe, we will be creating a joint histogram of PM2.5 mass +concentration and 2-metre temperature. +""" + +# %% +# 1. Import cf-python, numpy and matplotlib.pyplot: + +import matplotlib.pyplot as plt +import numpy as np + +import cf + +# %% +# 2. Read the field constructs using read function: +f = cf.read("~/recipes/levtype_sfc.nc") +print(f) + +# %% +# 3. Select the PM2.5 mass concentration and 2-metre temperature fields by +# index and print the description to show properties of all constructs: +pm25_field = f[2] +pm25_field.dump() + +# %% + +temp_field = f[3] +temp_field.dump() + +# %% +# 4. Convert the units to degree celsius for the temperature field: +temp_field.units = "degC" +temp_field.get_property("units") + +# %% +# 5. Digitize the PM2.5 mass concentration and 2-metre temperature fields. +# This step counts the number of values in each of the 10 equally sized bins +# spanning the range of the values. The ``'return_bins=True'`` argument makes +# sure that the calculated bins are also returned: +pm25_indices, pm25_bins = pm25_field.digitize(10, return_bins=True) +temp_indices, temp_bins = temp_field.digitize(10, return_bins=True) + +# %% +# 6. Create a joint histogram of the digitized fields: +joint_histogram = cf.histogram(pm25_indices, temp_indices) + +# %% +# 7. Get histogram data from the ``joint_histogram``: +histogram_data = joint_histogram.array + +# %% +# 8. Calculate bin centres for PM2.5 and temperature bins: +pm25_bin_centers = [(a + b) / 2 for a, b in pm25_bins] +temp_bin_centers = [(a + b) / 2 for a, b in temp_bins] + + +# %% +# 9. Create grids for PM2.5 and temperature bins using `numpy.meshgrid +# `_: +temp_grid, pm25_grid = np.meshgrid(temp_bin_centers, pm25_bin_centers) + +# %% +# 10. Plot the joint histogram using `matplotlib.pyplot.pcolormesh +# `_. Use `cf.Field.unique +# `_ to get +# the unique data array values and show the bin boundaries as ticks on the +# plot. ``'style='plain'`` in `matplotlib.axes.Axes.ticklabel_format +# `_ +# disables the scientific notation on the y-axis: +plt.pcolormesh(temp_grid, pm25_grid, histogram_data, cmap="viridis") +plt.xlabel("2-metre Temperature (degC)") +plt.ylabel("PM2.5 Mass Concentration (kg m**-3)") +plt.xticks(temp_bins.unique().array, rotation=45) +plt.yticks(pm25_bins.unique().array) +plt.colorbar(label="Frequency") +plt.title("Joint Histogram of PM2.5 and 2-metre Temperature", y=1.05) +plt.ticklabel_format(style="plain") +plt.show() diff --git a/recipes-docs/source/recipes-source/plot_10_recipe.py b/recipes-docs/source/recipes-source/plot_10_recipe.py new file mode 100644 index 0000000000..8c1a3f65dd --- /dev/null +++ b/recipes-docs/source/recipes-source/plot_10_recipe.py @@ -0,0 +1,96 @@ +""" +Calculating and plotting the relative vorticity +=============================================== + +Vorticity, the microscopic measure of rotation in a fluid, is a vector field +defined as the curl of velocity +`(James R. Holton, Gregory J. Hakim, An Introduction to Dynamic Meteorology, +2013, Elsevier : Academic Press p95-125) +`_. +In this recipe, we will be calculating and plotting the relative vorticity +from the wind components. +""" + +# %% +# 1. Import cf-python and cf-plot: + +import cfplot as cfp + +import cf + +# %% +# 2. Read the field constructs: +f = cf.read("~/recipes/ERA5_monthly_averaged_pressure_levels.nc") +print(f) + +# %% +# 3. Select wind components and look at their contents: +u = f.select_field("eastward_wind") +print(u) + +# %% + +v = f.select_field("northward_wind") +print(v) + +# %% +# 4. Create a date-time object for the required time period: +jan_2023 = cf.year(2023) & cf.month(1) + +# %% +# 5. The relative vorticity is calculated using `cf.curl_xy +# `_ and +# plotted using `cfplot.con `_. +# The ``with cf.relaxed_identities(True)`` context manager statement prevents +# the curl operation broadcasting across the two ``expver`` dimensions because +# it can't be certain that they are the same as they lack the standardised +# metadata. Setting +# ``cf.relaxed_identities(True)`` allows the ``long_name`` to be treated +# as standardised metadata. Since the horizontal coordinates are latitude and +# longitude, the +# `cf.curl_xy `_ +# function automatically accounts for the Earth's spherical geometry when +# calculating the spatial derivatives in the horizontal directions, and for this +# it requires the Earth's radius. In this case the radius is not stored in the +# wind fields, so must be provided by setting ``radius="earth"`` keyword +# parameter. While plotting, the relative vorticity is subspaced for January +# 2023 and one of the `experiment versions` using the dictionary unpacking +# operator (``**``) as there is an equal to sign in the identifier +# (``"long_name=expver"``): + +with cf.relaxed_identities(True): + rv = cf.curl_xy(u, v, radius="earth") + +cfp.con( + rv.subspace(T=jan_2023, **{"long_name=expver": 1}), + lines=False, + title="Relative Vorticity", +) + +# %% +# 6. Although the X axis is cyclic, it is not recognised as such, owing to the +# fact that the longitude coordinate bounds are missing. This results in +# discontinuities in the calculated vorticity field on the plot at the +# wrap-around location of 0 degrees east. The cyclicity could either be set on +# the field itself or just in the curl command by setting ``'x_wrap=True'`` +# while calculating the relative vorticity. Setting ``rv.units = "s-1"``, +# ensures that the units of the relative vorticity field are consistent with +# the calculation and the physical interpretation of the quantity: + +print(v.coordinate("X").has_bounds()) + +# %% + +with cf.relaxed_identities(True): + rv = cf.curl_xy(u, v, x_wrap=True, radius="earth") + +rv.units = "s-1" +print(rv) + +# %% + +cfp.con( + rv.subspace(T=jan_2023, **{"long_name=expver": 1}), + lines=False, + title="Relative Vorticity", +) diff --git a/recipes-docs/source/recipes-source/plot_11_recipe.py b/recipes-docs/source/recipes-source/plot_11_recipe.py new file mode 100644 index 0000000000..b147e51809 --- /dev/null +++ b/recipes-docs/source/recipes-source/plot_11_recipe.py @@ -0,0 +1,92 @@ +""" +Plotting the Warming Stripes +============================ + +In this recipe, we will plot the `Warming Stripes (Climate Stripes) +`_ created by +Professor Ed Hawkins at NCAS, University of Reading. Here we will use the +ensemble mean of the +`HadCRUT.5.0.1.0 analysis gridded data +`_ for +the same. + +""" + +# %% +# 1. Import cf-python and matplotlib.pyplot: + +import matplotlib.pyplot as plt + +import cf + +# %% +# 2. Read the field constructs: +temperature_data = cf.read( + "~/recipes/HadCRUT.5.0.1.0.analysis.anomalies.ensemble_mean.nc" +)[0] +print(temperature_data) + +# %% +# 3. Calculate the annual mean temperature anomalies. The ``'weights=True'`` +# argument is used to take the varying lengths of months into account which +# ensures that the calculated mean is more accurate: +annual_temperature = temperature_data.collapse( + "T: mean", weights=True, group=cf.Y() +) + +# %% +# 4. Select the data from 1850 to 2022: +period = annual_temperature.subspace(T=cf.year(cf.wi(1850, 2022))) + +# %% +# 5. Calculate the global average temperature for each year: +global_temperature = period.collapse("X: Y: mean") + +# %% +# 6. Get the global average temperature and squeeze it to remove the size 1 axis: +global_avg_temp = global_temperature.array.squeeze() + +# %% +# 7. Create a normalisation function that maps the interval from the minimum to +# the maximum temperature to the interval [0, 1] for colouring: +norm_global = plt.Normalize(global_avg_temp.min(), global_avg_temp.max()) + +# %% +# 8. Set the colormap instance: +cmap = plt.get_cmap("RdBu_r") + +# %% +# 9. Create the figure and the axes for the global plot. Loop over the selected +# years, plot a colored vertical stripe for each and remove the axes: +fig_global, ax_global = plt.subplots(figsize=(10, 2)) + +for i in range(global_avg_temp.shape[0]): + ax_global.axvspan( + xmin=i - 0.5, xmax=i + 0.5, color=cmap(norm_global(global_avg_temp[i])) + ) + +ax_global.axis("off") + +plt.show() + +# %% +# 10. For the regional warming stripes, steps 5 to 9 are repeated for the +# specific region. Here, we define the bounding box for UK by subspacing over +# a domain spanning 49.9 to 59.4 degrees north and -10.5 to 1.8 degrees east: +uk_temperature = period.subspace(X=cf.wi(-10.5, 1.8), Y=cf.wi(49.9, 59.4)) +uk_avg_temperature = uk_temperature.collapse("X: Y: mean") +uk_avg_temp = uk_avg_temperature.array.squeeze() +norm_uk = plt.Normalize(uk_avg_temp.min(), uk_avg_temp.max()) + +# %% + +fig_uk, ax_uk = plt.subplots(figsize=(10, 2)) + +for i in range(uk_avg_temp.shape[0]): + ax_uk.axvspan( + xmin=i - 0.5, xmax=i + 0.5, color=cmap(norm_uk(uk_avg_temp[i])) + ) + +ax_uk.axis("off") + +plt.show() diff --git a/recipes-docs/source/recipes-source/plot_12_recipe.py b/recipes-docs/source/recipes-source/plot_12_recipe.py new file mode 100644 index 0000000000..b09db0b29f --- /dev/null +++ b/recipes-docs/source/recipes-source/plot_12_recipe.py @@ -0,0 +1,117 @@ +""" +Using mask to plot Aerosol Optical Depth +======================================== + +In this recipe, we will make use of a +`masked array +`_ +to plot the `high-quality` retrieval of Aerosol Optical Depth (AOD) from all other +retrievals. + +""" + +# %% +# 1. Import cf-python, cf-plot and matplotlib.pyplot: + +import matplotlib.pyplot as plt +import cfplot as cfp + +import cf + +# %% +# 2. Read the field constructs: +fl = cf.read( + "~/recipes/JRR-AOD_v3r0_npp_s202012310752331_e202012310753573_c202100000000000.nc" +) +print(fl) + +# %% +# 3. Select AOD from the field list by identity and look at the contents: +aod = fl.select_field("long_name=AOT at 0.55 micron for both ocean and land") +print(aod) + +# %% +# 4. Select AOD retrieval quality by index and look at the quality flags: +quality = fl[13] +print(quality) + +# %% +# 5. Select latitude and longitude dimensions by identities, with two different +# techniques: +lon = aod.coordinate("long_name=Longitude") +lat = aod.coordinate("Y") + +# %% +# 6. Plot the AOD for all the retrievals using +# `cfplot.con `_. Here the argument +# ``'ptype'`` specifies the type of plot to use (latitude-longitude here) and +# the argument ``'lines=False'`` does not draw contour lines: +cfp.con(f=aod.array, x=lon.array, y=lat.array, ptype=1, lines=False) + +# %% +# 7. Create a mask for AOD based on the quality of the retrieval. The +# ``'__ne__'`` method is an implementation of the ``!=`` operator. It is used to +# create a mask where all the `high-quality` AOD points (with the flag 0) are +# marked as ``False``, and all the other data points (medium quality, low +# quality, or no retrieval) are marked as ``True``: +mask = quality.array.__ne__(0) + +# %% +# 8. Apply the mask to the AOD dataset. The ``'where'`` function takes the +# mask as an input and replaces all the values in the AOD dataset that +# correspond to ``True`` in the mask with a masked value using `cf.masked +# `_. +# In this case, all AOD values that are not of `high-quality` (since they were +# marked as ``True`` in the mask) are masked. This means that the ``high`` +# variable contains only the AOD data that was retrieved with `high-quality`: +high = aod.where(mask, cf.masked) + +# %% +# 9. Now plot both the AOD from `high-quality` retrieval and all other retrievals +# using `cfplot.con `_. Here: +# +# - `cfplot.gopen `_ is used to +# define the parts of the plot area, specifying that the figure should have +# 1 row and 2 columns, which is closed by +# `cfplot.gclose `_; +# - `plt.suptitle `_ +# is used to add a title for the whole figure; +# - the subplots for plotting are selected using +# `cfplot.gpos `_ after which +# `cfplot.mapset `_ is used to +# set the map limits and resolution for the subplots; +# - and as cf-plot stores the plot in a plot object with the name +# ``cfp.plotvars.plot``, country borders are added using normal +# `Cartopy operations `_ +# on the ``cfp.plotvars.mymap`` object: +import cartopy.feature as cfeature + +cfp.gopen(rows=1, columns=2, bottom=0.2) +plt.suptitle("AOD for both ocean and land", fontsize=20) +cfp.gpos(1) +cfp.mapset(resolution="50m", lonmin=68, lonmax=98, latmin=7, latmax=36) +cfp.con( + f=aod.array, + x=lon.array, + y=lat.array, + ptype=1, + lines=False, + title="All retrievals", + colorbar=None, +) +cfp.plotvars.mymap.add_feature(cfeature.BORDERS) +cfp.gpos(2) +cfp.mapset(resolution="50m", lonmin=68, lonmax=98, latmin=7, latmax=36) +cfp.con( + f=high.array, + x=lon.array, + y=lat.array, + ptype=1, + lines=False, + title="High quality retrieval", + colorbar_position=[0.1, 0.20, 0.8, 0.02], + colorbar_orientation="horizontal", + colorbar_title="AOD at 0.55 $\mu$", +) +cfp.plotvars.mymap.add_feature(cfeature.BORDERS) +cfp.gclose() diff --git a/recipes-docs/source/recipes-source/plot_13_recipe.py b/recipes-docs/source/recipes-source/plot_13_recipe.py new file mode 100644 index 0000000000..bf0398713e --- /dev/null +++ b/recipes-docs/source/recipes-source/plot_13_recipe.py @@ -0,0 +1,268 @@ +""" +Calculate and plot the Niño 3.4 Index +===================================== + +In this recipe, we will calculate and plot the sea surface temperature (SST) +anomaly in the Niño 3.4 region. According to `NCAR Climate Data Guide +`_, +the Niño 3.4 anomalies may be thought of as representing the average equatorial +SSTs across the Pacific from about the dateline to the South American coast. +The Niño 3.4 index typically uses a 5-month running mean, and El Niño or La +Niña events are defined when the Niño 3.4 SSTs exceed +/- 0.4 degrees Celsius for a +period of six months or more. + +""" + +# %% +# 1. Import cf-python and cf-plot, as well as some other libraries for use +# in next steps. + +import cartopy.crs as ccrs +import matplotlib.patches as mpatches + +import cfplot as cfp + +import cf + + +# %% +# 2. Read and select the SST by index and look at its contents: +sst = cf.read("~/recipes/ERA5_monthly_averaged_SST.nc")[0] +print(sst) + +# %% +# 3. Set the units from Kelvin to degrees Celsius: +sst.Units = cf.Units("degreesC") + +# %% +# 4. SST is subspaced for the Niño 3.4 region (5N-5S, 170W-120W) and as the +# dataset is using longitudes in 0-360 degrees East format, they are subtracted +# from 360 to convert them: +region = sst.subspace(X=cf.wi(360 - 170, 360 - 120), Y=cf.wi(-5, 5)) + +# %% +# 5. Plot the various Niño regions using cf-plot. Here: +# +# - `cfplot.gopen `_ is used to +# define the parts of the plot area, which is closed by +# `cfplot.gclose `_; +# - `cfplot.mapset `_ is used to +# set the map limits and projection; +# - `cfplot.setvars `_ is used to +# set various attributes of the plot, like setting the land colour to grey; +# - `cfplot.cscale `_ is used to +# choose one of the colour maps amongst many available; +# - `cfplot.con `_ plots contour data +# from the ``region`` subspace at a specific time with no contour lines and a +# title; +# - next, four Niño regions and labels are defined using +# `Matplotlib's Rectangle `_ +# and +# `Text `_ +# function with cf-plot plot object (``cfp.plotvars.plot``): + +cfp.gopen() +cfp.mapset(proj="cyl", lonmin=0, lonmax=360, latmin=-90, latmax=90) +cfp.setvars(land_color="grey") +cfp.cscale(scale="scale1") +cfp.con( + region.subspace(T=cf.dt(2022, 12, 1, 0, 0, 0, 0)), + lines=False, + title="Niño Index Regions", +) + +# Niño 3.4 region(5N-5S, 170W-120W): +rectangle = mpatches.Rectangle( + (-170, -5), + 50, + 10, + fill=False, + linewidth=1, + edgecolor="black", + transform=ccrs.PlateCarree(), +) +cfp.plotvars.mymap.add_patch(rectangle) +cfp.plotvars.mymap.text( + -145, + 7, + "3.4", + horizontalalignment="center", + fontsize=14, + weight="bold", + transform=ccrs.PlateCarree(), +) + +# Niño 1+2 region (0-10S, 90W-80W): +rectangle = mpatches.Rectangle( + (-90, 0), + 10, + 10, + hatch="**", + fill=False, + linewidth=1, + edgecolor="black", + alpha=0.3, + transform=ccrs.PlateCarree(), +) +cfp.plotvars.mymap.add_patch(rectangle) +cfp.plotvars.mymap.text( + -85, + 3, + "1+2", + horizontalalignment="center", + fontsize=8, + weight="bold", + transform=ccrs.PlateCarree(), +) + +# Niño 3 region (5N-5S, 150W-90W): +rectangle = mpatches.Rectangle( + (-150, -5), + 60, + 10, + hatch="xxx", + fill=False, + linewidth=1, + edgecolor="black", + alpha=0.3, + transform=ccrs.PlateCarree(), +) +cfp.plotvars.mymap.add_patch(rectangle) +cfp.plotvars.mymap.text( + -120, + -3, + "3", + horizontalalignment="center", + fontsize=14, + weight="bold", + transform=ccrs.PlateCarree(), +) + +# Niño 4 region (5N-5S, 160E-150W): +rectangle = mpatches.Rectangle( + (-200, -5), + 50, + 10, + hatch="oo", + fill=False, + linewidth=1, + edgecolor="black", + alpha=0.3, + transform=ccrs.PlateCarree(), +) +cfp.plotvars.mymap.add_patch(rectangle) +cfp.plotvars.mymap.text( + -175, + -3, + "4", + horizontalalignment="center", + fontsize=14, + weight="bold", + transform=ccrs.PlateCarree(), +) +cfp.gclose() + +# %% +# 6. Calculate the Niño 3.4 index and standardise it to create an anomaly index. +# The `collapse `_ +# method is used to calculate the mean over the longitude (X) and latitude (Y) +# dimensions: +nino34_index = region.collapse("X: Y: mean") + +# %% +# 7. The result, ``nino34_index``, represents the average SST in the defined +# Niño 3.4 region for each time step. In the variable ``base_period``, +# ``nino34_index`` is subset to only include data from the years 1961 to 1990. +# This period is often used as a reference period for calculating anomalies. +# The variables ``climatology`` and ``std_dev`` include the mean and the +# standard deviation over the time (T) dimension of the ``base_period`` data +# respectively: +base_period = nino34_index.subspace(T=cf.year(cf.wi(1961, 1990))) +climatology = base_period.collapse("T: mean") +std_dev = base_period.collapse("T: sd") + +# %% +# 8. The line for variable ``nino34_anomaly`` calculates the standardised +# anomaly for each time step in the ``nino34_index`` data. It subtracts the +# ``climatology`` from the ``nino34_index`` and then divides by the ``std_dev``. +# The resulting ``nino34_anomaly`` data represents how much the SST in the Niño +# 3.4 region deviates from the 1961-1990 average, in units of standard +# deviations. This is a common way to quantify climate anomalies like El Niño +# and La Niña events: +nino34_anomaly = (nino34_index - climatology) / std_dev + +# %% +# 9. A moving average of the ``nino34_anomaly`` along the time axis, with a +# window size of 5 (i.e. an approximately 5-month moving average) is calculated +# using the +# `moving_window `_ +# method. The ``mode='nearest'`` parameter is used to specify how to pad the +# data outside of the time range. The resulting ``nino34_rolling`` variable +# represents a smoothed version of the ``nino34_anomaly`` data. It removes +# short-term fluctuations and highlights longer-term trends or cycles: +nino34_rolling = nino34_anomaly.moving_window( + method="mean", window_size=5, axis="T", mode="nearest" +) + +# %% +# 10. Define El Niño and La Niña events by creating Boolean masks to identify +# El Niño and La Niña events. Now plot SST anomalies in the Niño 3.4 region over +# time using cf-plot. Here: +# +# - `cfplot.gset `_ sets the limits +# of the x-axis (years from 1940 to 2022) and y-axis (anomalies from -3 +# degrees C to 3 degrees C) for the plot; +# - `cfplot.gopen `_ is used to +# define the parts of the plot area, which is closed by +# `cfplot.gclose `_; +# - `cfplot.lineplot `_ plots +# the rolling Niño 3.4 index over time; +# - a zero line and also horizontal dashed lines are drawn for El Niño and +# La Niña thresholds using +# `Matplotlib's axhline `_ +# with cf-plot plot object (``cfp.plotvars.plot``); +# - `fill_between `_ +# from Matplotlib is used with cf-plot plot object (``cfp.plotvars.plot``) +# to fill the area between the Niño 3.4 index and the El Niño/La Niña +# thresholds; +# - similarly, +# `cfplot.plotvars.plot.legend `_ +# is used to add a legend in the end: +elnino = nino34_rolling >= 0.4 +lanina = nino34_rolling <= -0.4 + +cfp.gset(xmin="1940-1-1", xmax="2022-12-31", ymin=-3, ymax=3) + +cfp.gopen(figsize=(10, 6)) +cfp.lineplot( + nino34_rolling, + color="black", + title="SST Anomaly in Niño 3.4 Region (5N-5S, 120-170W)", + ylabel="Temperature anomaly ($\degree C$)", + xlabel="Year", +) +cfp.plotvars.plot.axhline( + 0.4, color="red", linestyle="--", label="El Niño Threshold" +) +cfp.plotvars.plot.axhline( + -0.4, color="blue", linestyle="--", label="La Niña Threshold" +) +cfp.plotvars.plot.axhline(0, color="black", linestyle="-", linewidth=1) +cfp.plotvars.plot.fill_between( + nino34_rolling.coordinate("T").array, + 0.4, + nino34_rolling.array.squeeze(), + where=elnino.squeeze(), + color="red", + alpha=0.3, +) +cfp.plotvars.plot.fill_between( + nino34_rolling.coordinate("T").array, + -0.4, + nino34_rolling.array.squeeze(), + where=lanina.squeeze(), + color="blue", + alpha=0.3, +) +cfp.plotvars.plot.legend(frameon=False, loc="lower center", ncol=2) +cfp.gclose() diff --git a/recipes-docs/source/recipes-source/plot_14_recipe.py b/recipes-docs/source/recipes-source/plot_14_recipe.py new file mode 100644 index 0000000000..957a3f38d7 --- /dev/null +++ b/recipes-docs/source/recipes-source/plot_14_recipe.py @@ -0,0 +1,181 @@ +""" +Overlay Geopotential height contours over Temperature anomalies +=============================================================== + +In this recipe, we will overlay Geopotential height contours over Temperature +anomalies to help analyse meteorological conditions during July 2018, +specifically focusing on the significant concurrent extreme events that occurred +during the 2018 boreal spring/summer season in the Northern Hemisphere. + +""" + +# %% +# 1. Import cf-python and cf-plot: +import cfplot as cfp + +import cf + +# %% +# 2. Read and select the 200 hPa geopotential by index and look at its contents: +gp = cf.read("~/recipes/ERA5_monthly_averaged_z200.nc")[0] +print(gp) + +# %% +# 3. Convert the geopotential data to geopotential height by dividing it by the +# acceleration due to gravity (approximated as 9.81 :math:`m \cdot {s}^{-2}`): +gph = gp / 9.81 + +# %% +# 4. Subset the geopotential height to extract data specifically for July 2018, +# a significant month due to heat extremes and heavy rainfall: +gph_july = gph.subspace(T=cf.month(7) & cf.year(2018)).squeeze() + +# %% +# 5. Plot contour lines of this geopotential height for July 2018. Here: +# +# - `cfplot.gopen `_ is used to +# define the parts of the plot area, which is closed by +# `cfplot.gclose `_; +# - `cfplot.mapset `_ is used to +# set the map projection to North Polar Stereographic; +# - `cfplot.setvars `_ is used to +# set various attributes of the plot, like setting the thickness of the lines +# that represent continents; +# - `cfplot.con `_ plots the contour +# lines representing the 200 hPa geopotential height values without filling +# between the contour lines (``fill=False``) and no colour bar +# (``colorbar=False``); +# - `cfplot.levs `_ is used to +# specify two contour levels, 12000 and 12300 m, corresponding to the +# approximate polar-front jet and subtropical jet respectively; +# - `cfplot.con `_ is again used to +# plot the contour lines for polar-front jet and subtropical jet with a +# thicker line width; +# - `cfp.plotvars.mymap.stock_img() `_ +# then finally visualises the Earth's surface in cf-plot's +# ``cfp.plotvars.mymap`` plot object: +cfp.gopen() +cfp.mapset(proj="npstere") +cfp.setvars(continent_thickness=0.5) + +cfp.con( + f=gph_july, + fill=False, + lines=True, + line_labels=False, + colors="black", + linewidths=1, + colorbar=False, +) + +cfp.levs(manual=[12000, 12300]) +cfp.con( + f=gph_july, + fill=False, + lines=True, + colors="black", + linewidths=3.0, + colorbar=False, +) + +cfp.plotvars.mymap.stock_img() +cfp.gclose() + +# %% +# 6. Read and select the 2-metre temperature by index and look at its contents: +t2m = cf.read("~/recipes/ERA5_monthly_averaged_t2m.nc")[0] +print(t2m) + +# %% +# 7. Set the units from Kelvin to degrees Celsius: +t2m.Units = cf.Units("degreesC") + +# %% +# 8. Extract a subset for July across the years for ``t2m``: +t2m_july = t2m.subspace(T=cf.month(7)) + +# %% +# 9. The 2-meter temperature climatology is then calculated for the month of +# July over the period from 1981 to 2010, which provides a baseline against +# which anomalies in later years are compared: +t2m_july_climatology = t2m_july.subspace( + T=cf.year(cf.wi(1981, 2010)) +).collapse("T: mean") + +# %% +# 10. Calculate the temperature anomaly for the month of July in the year 2018 +# relative to the climatological baseline (``t2m_july_climatology``). This +# indicates how much the temperatures for that month in that year deviated from +# the long-term average for July across the 1981-2010 period: +t2m_july_anomaly_2018 = ( + t2m_july.subspace(T=cf.year(2018)).squeeze() - t2m_july_climatology +) + +# %% +# 11. +# The July 2018 season experienced extreme heat in many parts of the Northern +# Hemisphere. This period's extreme events were related to unusual +# meteorological conditions, particularly abnormalities in the jet stream. To +# provide an insight into the atmospheric conditions, the temperature anomalies +# and the geopotential height contours are plotted using cf-plot. Here: +# +# - `cfplot.gopen `_ is used to +# define the parts of the plot area, which is closed by +# `cfplot.gclose `_; +# - `cfplot.mapset `_ is used to +# set the map projection to Robinson; +# - `cfplot.setvars `_ is used to +# set various attributes of the plot, like setting the thickness of the lines +# that represent continents and master title properties; +# - `cfplot.levs `_ is used to +# specify the contour levels for temperature anomalies, starting from -2 to 2 +# with an interval of 0.5; +# - `cfplot.cscale `_ is used to +# choose one of the colour maps amongst many available; +# - `cfplot.con `_ plots contour fill +# of temperature anomalies without contour lines (``lines=False``); +# - `cfplot.levs() `_ is used to +# reset contour levels to default after which the steps to plot the contour +# lines representing the 200 hPa geopotential height values, the approximate +# polar-front jet and subtropical jet from Step 5 are repeated: +cfp.gopen() +cfp.mapset(proj="robin") +cfp.setvars( + continent_thickness=0.5, + master_title="July 2018", + master_title_fontsize=22, + master_title_location=[0.53, 0.83], +) + +cfp.levs(min=-2, max=2, step=0.5) +cfp.cscale("temp_19lev") +cfp.con( + f=t2m_july_anomaly_2018, + lines=False, + colorbar_title="Temperature anomaly relative to 1981-2010 ($\degree C$)", + colorbar_fontsize=13, + colorbar_thick=0.04, +) + +cfp.levs() +cfp.con( + f=gph_july, + fill=False, + lines=True, + line_labels=False, + colors="black", + linewidths=1, + colorbar=False, +) + +cfp.levs(manual=[12000, 12300]) +cfp.con( + f=gph_july, + fill=False, + lines=True, + colors="black", + linewidths=3.0, + colorbar=False, +) + +cfp.gclose() diff --git a/recipes-docs/source/recipes-source/plot_15_recipe.py b/recipes-docs/source/recipes-source/plot_15_recipe.py new file mode 100644 index 0000000000..d20f622826 --- /dev/null +++ b/recipes-docs/source/recipes-source/plot_15_recipe.py @@ -0,0 +1,163 @@ +""" +Resampling Land Use Flags to a Coarser Grid +=========================================== + +In this recipe, we will compare the land use distribution in different countries +using a land use data file and visualize the data as a histogram. This will help +to understand the proportion of different land use categories in each country. + +The land use data is initially available at a high spatial resolution of +approximately 100 m, with several flags defined with numbers representing the +type of land use. Regridding the data to a coarser resolution of approximately +25 km would incorrectly represent the flags on the new grids. + +To avoid this, we will resample the data to the coarser resolution by +aggregating the data within predefined spatial regions or bins. This approach +will give a dataset where each 25 km grid cell contains a histogram of land use +flags, as determined by the original 100 m resolution data. It retains the +original spatial extent of the data while reducing its spatial complexity. +Regridding, on the other hand, involves interpolating the data onto a new grid, +which can introduce artefacts and distortions in the data. + +""" + +import cartopy.io.shapereader as shpreader +import matplotlib.pyplot as plt +import numpy as np + +# %% +# 1. Import the required libraries. We will use Cartopy's ``shapereader`` to +# work with shapefiles that define country boundaries: +import cf + +# %% +# 2. Read and select land use data by index and see properties of +# all constructs: +f = cf.read("~/recipes/output.tif.nc")[0] +f.dump() + + +# %% +# 3. Define a function to extract data for a specific country: +# +# - The ``extract_data`` function is defined to extract land use data for a +# specific country, specified by the ``country_name`` parameter. +# - It uses the `Natural Earth `_ +# shapefile to get the bounding coordinates of the selected country. +# - The `shpreader.natural_earth `_ +# function is called to access the Natural +# Earth shapefile of country boundaries with a resolution of 10 m. +# - The `shpreader.Reader `_ +# function reads the shapefile, and the selected country's record is retrieved +# by filtering the records based on the ``NAME_LONG`` attribute. +# - The bounding coordinates are extracted using the ``bounds`` attribute of the +# selected country record. +# - The land use data file is then read and subset using these bounding +# coordinates with the help of the ``subspace`` function. The subset data is +# stored in the ``f`` variable. + + +def extract_data(country_name): + shpfilename = shpreader.natural_earth( + resolution="10m", category="cultural", name="admin_0_countries" + ) + reader = shpreader.Reader(shpfilename) + country = [ + country + for country in reader.records() + if country.attributes["NAME_LONG"] == country_name + ][0] + lon_min, lat_min, lon_max, lat_max = country.bounds + + f = cf.read("~/recipes/output.tif.nc")[0] + f = f.subspace(X=cf.wi(lon_min, lon_max), Y=cf.wi(lat_min, lat_max)) + + return f + + +# %% +# 4. Define a function to plot a histogram of land use distribution for a +# specific country: +# +# - The `digitize `_ +# function of the ``cf.Field`` object is called to convert the land use data +# into indices of bins. It takes an array of bins (defined by +# the `np.linspace `_ function) +# and the ``return_bins=True`` parameter, which returns the actual bin values +# along with the digitized data. +# - The `np.linspace `_ +# function is used to create an array of evenly spaced bin edges from 0 to 50, +# with 51 total values. This creates bins of width 1. +# - The ``digitized`` variable contains the bin indices for each data point, +# while the bins variable contains the actual bin values. +# - The `cf.histogram `_ +# function is called on the digitized data to create a histogram. This +# function returns a field object with the histogram data. +# - The `squeeze `_ +# function applied to the histogram ``array`` extracts the histogram data as a NumPy +# array and removes any single dimensions. +# - The ``total_valid_sub_cells`` variable calculates the total number of valid +# subcells (non-missing data points) by summing the histogram data. +# - The last element of the bin_counts array is removed with slicing +# (``bin_counts[:-1]``) to match the length of the ``bin_indices`` array. +# - The ``percentages`` variable calculates the percentage of each bin by +# dividing the ``bin_counts`` by the ``total_valid_sub_cells`` and multiplying +# by 100. +# - The ``bin_indices`` variable calculates the centre of each bin by averaging +# the bin edges. This is done by adding the ``bins.array[:-1, 0]`` and +# ``bins.array[1:, 0]`` arrays and dividing by 2. +# - The ``ax.bar`` function is called to plot the histogram as a bar chart on +# the provided axis. The x-axis values are given by the ``bin_indices`` array, +# and the y-axis values are given by the ``percentages`` array. +# - The title, x-axis label, y-axis label, and axis limits are set based on the +# input parameters. + + +def plot_histogram(field, ax, label, ylim, xlim): + digitized, bins = field.digitize(np.linspace(0, 50, 51), return_bins=True) + + h = cf.histogram(digitized) + bin_counts = h.array.squeeze() + + total_valid_sub_cells = bin_counts.sum() + + bin_counts = bin_counts[:-1] + + percentages = bin_counts / total_valid_sub_cells * 100 + + bin_indices = (bins.array[:-1, 0] + bins.array[1:, 0]) / 2 + + ax.bar(bin_indices, percentages, label=label) + ax.set_title(label) + ax.set_xlabel("Land Use Flag") + ax.set_ylabel("Percentage") + ax.set_ylim(ylim) + ax.set_xlim(xlim) + + +# %% +# 5. Define the countries of interest: +countries = ["Ireland", "Belgium", "Switzerland"] + +# %% +# 6. Set up the figure and axes for plotting the histograms: +# +# - The ``plt.subplots`` function is called to set up a figure with three +# subplots, with a figure size of 8 inches by 10 inches. +# - A loop iterates over each country in the countries list and for each +# country, the ``extract_data`` function is called to extract its land use +# data. +# - The ``plot_histogram`` function is then called to plot the histogram of land +# use distribution on the corresponding subplot. +# - The ``plt.tight_layout`` function is called to ensure that the subplots are +# properly spaced within the figure and finally, the ``plt.show`` function +# displays the figure with the histograms. +fig, axs = plt.subplots(3, 1, figsize=(8, 10)) + +for i, country in enumerate(countries): + ax = axs[i] + data = extract_data(country) + plot_histogram(data, ax, label=country, ylim=(0, 50), xlim=(0, 50)) + +plt.tight_layout() +plt.show() diff --git a/recipes-docs/source/recipes-source/plot_16_recipe.py b/recipes-docs/source/recipes-source/plot_16_recipe.py new file mode 100644 index 0000000000..00a7f2e1d0 --- /dev/null +++ b/recipes-docs/source/recipes-source/plot_16_recipe.py @@ -0,0 +1,57 @@ +""" +Plotting contour subplots with different projections +==================================================== + +In this recipe, we will plot the same data using different projections +as subplots to illustrate visually some available possibilities. + +""" + +# %% +# 1. Import cf-python and cf-plot: + +import cfplot as cfp + +import cf + +# %% +# 2. Read the field in: +f = cf.read("~/recipes/ggap.nc")[0] + +# %% +# 3. List the projection types to use. Here we are using +# Cylindrical/Default, North Pole Stereographic, South Pole Stereographic, +# Mollweide, Mercator and Robinson. However there are several other choices +# possible, see: +# https://ncas-cms.github.io/cf-plot/build/user_guide.html#appendixc. Our +# chosen list is: +projtypes = ["cyl", "npstere", "spstere", "moll", "merc", "robin"] + +# %% +# 4. Create the file with subplots. If changing the number of subplots, +# ensure the number of rows * number of columns = the number of projections. +# Here we are doing 6 projections so 2 x 3 is fine. Then loop through the +# list of projection types and plot each as a sub-plot: +cfp.gopen(rows=2, columns=3, bottom=0.2) +for i, proj in enumerate(projtypes): + # gpos has 1 added to the index because it takes 1 as its first value + cfp.gpos(i + 1) + cfp.mapset(proj=proj) + + # For the final plot only, add a colour bar to cover all the sub-plots + if i == len(projtypes) - 1: + cfp.con( + f.subspace(pressure=850), + lines=False, + title=proj, + colorbar_position=[0.1, 0.1, 0.8, 0.02], + colorbar_orientation="horizontal", + ) + else: + cfp.con( + f.subspace(pressure=850), + lines=False, + title=proj, + colorbar=False, + ) +cfp.gclose() diff --git a/recipes-docs/source/recipes-source/plot_17_recipe.py b/recipes-docs/source/recipes-source/plot_17_recipe.py new file mode 100644 index 0000000000..c94769e2ba --- /dev/null +++ b/recipes-docs/source/recipes-source/plot_17_recipe.py @@ -0,0 +1,109 @@ +""" +Plotting contour subplots with different colour maps/scales +=========================================================== + +In this recipe, we will plot the same data with different colour maps from +three categories in separate subplots to illustrate the importance of +choosing a suitable one for given data. To avoid unintended bias and +misrepresentation, or lack of accessibility, a careful choice must be made. +""" + +# %% +# 1. Import cf-python and cf-plot: + +import matplotlib.pyplot as plt +import cfplot as cfp + +import cf + +# %% +# 2. Read the field in: +f = cf.read("~/recipes/ggap.nc")[0] + +# %% +# 3. Choose a set of predefined colour scales to view. These can be chosen +# from the selection at: +# https://ncas-cms.github.io/cf-plot/build/colour_scales.html or you +# can define your own, see: +# https://ncas-cms.github.io/cf-plot/build/colour_scales.html#user-defined-colour-scales. +# Here we take three colour scales each from three different general +# categories, to showcase some differences in representation. +# Note colour scale levels can be adjusted using 'cscale' and keywords such as: +# cfp.cscale(, ncols=16, below=2, above=14) + +# %% +# a. Perceptually uniform colour scales, with no zero value, see: +# https://ncas-cms.github.io/cf-plot/build/colour_scales.html#perceptually-uniform-colour-scales. +colour_scales_pu = ["viridis", "magma", "plasma"] + +# %% +# b. NCAR Command Language colour scales enhanced to help with colour +# blindness, see: +# https://ncas-cms.github.io/cf-plot/build/colour_scales.html#ncar-command-language-enhanced-to-help-with-colour-blindness. +# These colour maps are better for accessibility. +colour_scales_ncl = ["posneg_1", "GreenMagenta16", "StepSeq25"] + +# %% +# c. Orography/bathymetry colour scales, see: +# https://ncas-cms.github.io/cf-plot/build/colour_scales.html#orography-bathymetry-colour-scales. +# These are used to show the shape/contour of land masses, but bear in mind the +# data we show here is air temperature so doesn't represent this and +# therefore it is not a good choice in this case: +colour_scales_ob = ["wiki_1_0_2", "wiki_2_0", "wiki_2_0_reduced"] + + +# %% +# 4. We plot each category of colourmap in a given columns of the subplot, +# but given the 'gpos' function positions subplots from left to right, row by +# row from the top, we need to interleave the values in a list. We can use +# zip to do this: +colour_scales_columns = [ + cscale + for category in zip(colour_scales_pu, colour_scales_ncl, colour_scales_ob) + for cscale in category +] + + +# %% +# 5. Create the figure and give it an overall title. Ensure the +# number of rows * number of columns = number of colour scales. +# Then we loop through all the different colour maps defined and plot +# as subplots, with each category in the same column, labelling each column +# with the colour scale category: +cfp.gopen(rows=3, columns=3, bottom=0.1, top=0.85) +plt.suptitle( + ( + "Air temperature (K) at 850 mbar pressure shown in different " + "categories of colour scale" + ), + fontsize=18, +) +for i, colour_scale in enumerate(colour_scales_columns): + cfp.gpos(i + 1) + cfp.cscale(colour_scale, ncols=15) + + # For the topmost plots, label the column with the colour scale category + # using the 'title' argument, otherwise don't add a title. + # Ensure the order the titles are written in corresponds to the + # order unzipped in step 4, so the columns match up correctly. + if i == 0: + set_title = "Perceptually uniform\ncolour maps" + elif i == 1: + set_title = "NCL colour maps enhanced to \nhelp with colour blindness" + elif i == 2: + set_title = "Orography/bathymetry\ncolour maps" + else: + set_title = "" + + cfp.con( + f.subspace(pressure=850), + title=set_title, + lines=False, + axes=False, + colorbar_drawedges=False, + colorbar_title=f"Shown in '{colour_scale}'", + colorbar_fraction=0.04, + colorbar_thick=0.02, + colorbar_fontsize=11, + ) +cfp.gclose() diff --git a/recipes-docs/source/recipes-source/plot_18_recipe.py b/recipes-docs/source/recipes-source/plot_18_recipe.py new file mode 100644 index 0000000000..f0eae36e35 --- /dev/null +++ b/recipes-docs/source/recipes-source/plot_18_recipe.py @@ -0,0 +1,143 @@ +""" +Calculating the Pearson correlation coefficient between datasets +================================================================ + +In this recipe, we will take two datasets, one for an independent variable +(in this example elevation) and one for a dependent variable (snow +cover over a particular day), regrid them to the same resolution then +calculate the correlation coefficient, to get a measure of the relationship +between them. + +""" + +# %% +# 1. Import cf-python, cf-plot and other required packages: +import matplotlib.pyplot as plt +import scipy.stats.mstats as mstats +import cfplot as cfp + +import cf + + +# %% +# 2. Read the data in and unpack the Fields from FieldLists using indexing. +# In our example We are investigating the influence of the land height on +# the snow cover extent, so snow cover is the dependent variable. The snow +# cover data is the +# 'Snow Cover Extent 2017-present (raster 500 m), Europe, daily – version 1' +# sourced from the Copernicus Land Monitoring Service which is described at: +# https://land.copernicus.eu/en/products/snow/snow-cover-extent-europe-v1-0-500m +# and the elevation data is the 'NOAA NGDC GLOBE topo: elevation data' dataset +# which can be sourced from the IRI Data Library, or details found, at: +# http://iridl.ldeo.columbia.edu/SOURCES/.NOAA/.NGDC/.GLOBE/.topo/index.html. +orog = cf.read("~/recipes/1km_elevation.nc")[0] +snow = cf.read("~/recipes/snowcover")[0] + +# %% +# 3. Choose the day of pre-aggregated snow cover to investigate. We will +# take the first datetime element corresponding to the first day from the +# datasets, 1st January 2024, but by changing the indexing you can explore +# other days by changing the index. We also get the string corresponding to +# the date, to reference later: +snow_day = snow[0] +snow_day_dt = snow_day.coordinate("time")[0].data +snow_day_daystring = f"{snow_day_dt.datetime_as_string[0].split(' ')[0]}" + +# %% +# 4. Choose the region to consider to compare the relationship across, +# which must be defined across both datasets, though not necessarily on the +# same grid since we regrid to the same grid next and subspace to the same +# area for both datasets ready for comparison in the next steps. By changing +# the latitude and longitude points in the tuple below, you can change the +# area that is used: +region_in_mid_uk = ((-3.0, -1.0), (52.0, 55.0)) +sub_orog = orog.subspace( + longitude=cf.wi(*region_in_mid_uk[0]), latitude=cf.wi(*region_in_mid_uk[1]) +) +sub_snow = snow_day.subspace( + longitude=cf.wi(*region_in_mid_uk[0]), latitude=cf.wi(*region_in_mid_uk[1]) +) + +# %% +# 5. Ensure data quality, since the standard name here corresponds to a +# unitless fraction, but the values are in the tens, so we need to +# normalise these to all lie between 0 and 1 and change the units +# appropriately: +sub_snow = (sub_snow - sub_snow.minimum()) / (sub_snow.range()) +sub_snow.override_units("1", inplace=True) + +# %% +# 6. Regrid the data so that they lie on the same grid and therefore each +# array structure has values with corresponding geospatial points that +# can be statistically compared. Here the elevation field is regridded to the +# snow field since the snow is higher-resolution, but the other way round is +# possible by switching the field order: +regridded_orog = sub_orog.regrids(sub_snow, method="linear") + +# %% +# 7. Squeeze the snow data to remove the size 1 axes so we have arrays of +# the same dimensions for each of the two fields to compare: +sub_snow = sub_snow.squeeze() + +# %% +# 8. Finally, perform the statistical calculation by using the SciPy method +# to find the Pearson correlation coefficient for the two arrays now they are +# in comparable form. Note we need to use 'scipy.stats.mstats' and not +# 'scipy.stats' for the 'pearsonr' method, to account for masked +# data in the array(s) properly: +coefficient = mstats.pearsonr(regridded_orog.array, sub_snow.array) +print(f"The Pearson correlation coefficient is: {coefficient}") + +# %% +# 9. Make a final plot showing the two arrays side-by-side and quoting the +# determined Pearson correlation coefficient to illustrate the relationship +# and its strength visually. We use 'gpos' to position the plots in two +# columns and apply some specific axes ticks and labels for clarity. +cfp.gopen( + rows=1, + columns=2, + top=0.85, + user_position=True, +) + +# Joint configuration of the plots, including adding an overall title +plt.suptitle( + ( + "Snow cover compared to elevation for the same area of the UK " + f"aggregated across\n day {snow_day_daystring} with correlation " + "coefficient (on the same grid) of " + f"{coefficient.statistic:.4g} (4 s.f.)" + ), + fontsize=17, +) +cfp.mapset(resolution="10m") +cfp.setvars(ocean_color="white", lake_color="white") +label_info = { + "xticklabels": ("3W", "2W", "1W"), + "yticklabels": ("52N", "53N", "54N", "55N"), + "xticks": (-3, -2, -1), + "yticks": (52, 53, 54, 55), +} + +# Plot the two contour plots as columns +cfp.gpos(1) +cfp.cscale("wiki_2_0_reduced", ncols=11) +cfp.con( + regridded_orog, + lines=False, + title="Elevation (from 1km-resolution orography)", + colorbar_drawedges=False, + **label_info, +) +cfp.gpos(2) +# Don't add extentions on the colourbar since it can only be 0 to 1 inclusive +cfp.levs(min=0, max=1, step=0.1, extend="neither") +cfp.cscale("precip_11lev", ncols=11, reverse=1) +cfp.con( + sub_snow, + lines=False, + title="Snow cover extent (from satellite imagery)", + colorbar_drawedges=False, + **label_info, +) +cfp.gclose() diff --git a/recipes-docs/source/recipes-source/plot_19_recipe.py b/recipes-docs/source/recipes-source/plot_19_recipe.py new file mode 100644 index 0000000000..02d493dc21 --- /dev/null +++ b/recipes-docs/source/recipes-source/plot_19_recipe.py @@ -0,0 +1,108 @@ +""" +Plotting per-season trends in global sea surface temperature extrema +==================================================================== + +In this recipe we find the area-based extrema of global sea surface +temperature per month and, because it is very difficult to +interpret for trends when in a monthly form, we calculate and plot +on top of this the mean across each season for both the minima and the +maxima. +""" + +# %% +# 1. Import cf-python, cf-plot and other required packages: +import matplotlib.pyplot as plt +import cfplot as cfp + +import cf + +# %% +# 2. Read the dataset in extract the SST Field from the FieldList: +f = cf.read("~/recipes/ERA5_monthly_averaged_SST.nc") +sst = f[0] # this gives the sea surface temperature (SST) + +# %% +# 3. Collapse the SST data by area extrema (extrema over spatial dimensions): +am_max = sst.collapse("area: maximum") # equivalent to "X Y: maximum" +am_min = sst.collapse("area: minimum") # equivalent to "X Y: minimum" + +# %% +# 4. Reduce all timeseries down to just 1980+ since there are some data +# quality issues before 1970 and also this window is about perfect size +# for viewing the trends without the line plot becoming too cluttered: +am_max = am_max.subspace(T=cf.ge(cf.dt("1980-01-01"))) +am_min = am_min.subspace(T=cf.ge(cf.dt("1980-01-01"))) + +# %% +# 5. Create a mapping which provides the queries we need to collapse on +# the four seasons, along with our description of them, as a value, with +# the key of the string encoding the colour we want to plot these +# trend lines in. This structure will be iterated over to make our plot: +colours_seasons_mapping = { + "red": (cf.mam(), "Mean across MAM: March, April and May"), + "blue": (cf.jja(), "Mean across JJA: June, July and August"), + "green": (cf.son(), "Mean across SON: September, October and November"), + "purple": (cf.djf(), "Mean across DJF: December, January and February"), +} + +# %% +# 6. Create and open the plot file. Put maxima subplot at top since these +# values are higher, given increasing x axis. +# Note we set limits manually with 'gset' only to +# allow space so the legend doesn't overlap the data, which isn't +# possible purely from positioning it anywhere within the default plot. +# Otherwise cf-plot handles this for us. To plot the per-season means +# of the maxima, we loop through the season query mapping and do a +# "T: mean" collapse setting the season as the grouping: +cfp.gopen( + rows=2, + columns=1, + bottom=0.1, + top=0.85, +) +cfp.gpos(1) +cfp.gset(xmin="1980-01-01", xmax="2022-12-01", ymin=304, ymax=312) +for colour, season_query in colours_seasons_mapping.items(): + query_on_season, season_description = season_query + am_max_collapse = am_max.collapse("T: mean", group=query_on_season) + cfp.lineplot( + am_max_collapse, + color=colour, + markeredgecolor=colour, + marker="o", + label=season_description, + title="Maxima per month or season", + ) +cfp.lineplot( + am_max, + color="grey", + xlabel="", + label="All months", +) +# Create and add minima subplot below the maxima one. Just like for the +# maxima case, we plot per-season means by looping through the season query +# mapping and doing a "T: mean" collapse setting the season as the grouping +cfp.gpos(2) +cfp.gset(xmin="1980-01-01", xmax="2022-12-01", ymin=269, ymax=272) +for colour, season_query in colours_seasons_mapping.items(): + query_on_season, season_description = season_query + am_min_collapse = am_min.collapse("T: mean", group=query_on_season) + cfp.lineplot( + am_min_collapse, + color=colour, + markeredgecolor=colour, + marker="o", + xlabel="", + title="Minima per month or season", + ) +cfp.lineplot( + am_min, + color="grey", +) +# Add an overall title to the plot and close the file to save it +plt.suptitle( + "Global mean sea surface temperature (SST) monthly\nminima and maxima " + "showing seasonal means of these extrema", + fontsize=18, +) +cfp.gclose() diff --git a/recipes-docs/source/recipes-source/plot_20_recipe.py b/recipes-docs/source/recipes-source/plot_20_recipe.py new file mode 100644 index 0000000000..95e02b01d7 --- /dev/null +++ b/recipes-docs/source/recipes-source/plot_20_recipe.py @@ -0,0 +1,97 @@ +""" +Calculating and plotting the divergence of sea currents +======================================================= + +In this recipe, we will calculate the divergence of depth-averaged +currents in the Irish Sea, then plot the divergence as a contour +fill plot underneath the vectors themselves in the form of a vector plot. +""" + +# %% +# 1. Import cf-python and cf-plot: +import cfplot as cfp + +import cf + +# %% +# 2. Read the fields in. This dataset consists of depth-averaged eastward and +# northward current components plus the sea surface height above sea level and +# is a gridded dataset, with grid resolution of 1.85 km, covering the entire +# Irish Sea area. It was found via the CEDA Archive at the location of: +# https://catalogue.ceda.ac.uk/uuid/1b89e025eedd49e8976ee0721ec6e9b5, with +# DOI of https://dx.doi.org/10.5285/031e7ca1-9710-280d-e063-6c86abc014a0: +f = cf.read("~/recipes/POLCOMS_WAM_ZUV_01_16012006.nc") + +# %% +# 3. Get the separate vector components, which are stored as separate fields. +# The first, 'u', corresponds to the eastward component and the second, 'v', +# the northward component: +u = f[0] +v = f[1] + +# %% +# 4. Squeeze the fields to remove the size 1 axes in each case: +u = u.squeeze() +v = v.squeeze() + +# %% +# 5. Consider the currents at a set point in time. To do this we +# select one of the 720 datetime sample points in the fields to +# investigate, in this case by subspacing to pick out a particular +# datetime value we saw within the time coordinate data of the field (but +# you could also use indexing or filtering to select a specific value). +# Once we subspace to one datetime, we squeeze out the size 1 time axis +# in each case: +chosen_time = "2006-01-15 23:30:00" # 720 choices to pick from, try this one! +u_1 = u.subspace(T=cf.dt(chosen_time)) +v_1 = v.subspace(T=cf.dt(chosen_time)) +u_1 = u_1.squeeze() +v_1 = v_1.squeeze() + +# %% +# 6. +# When inspecting the u and v fields using cf inspection methods such as +# from print(u_1.data) and u_1.data.dump(), for example, we can see that there are +# lots of -9999 values in their data array, apparently used as a +# fill/placeholder value, including to indicate undefined data over the land. +# In order for these to not skew the data and dominate the plot, we need +# to mask values matching this, so that only meaningful values remain. +u_2 = u_1.where(cf.lt(-9000), cf.masked) +v_2 = v_1.where(cf.lt(-9000), cf.masked) + +# %% +# 7. Calculate the divergence using the 'div_xy' function operating on the +# vector eastward and northward components as the first and second argument +# respectively. We need to calculate this for the latitude-longitude plane +# of the Earth, defined in spherical polar coordinates, so we must specify +# the Earth's radius for the appropriate calculation: +div = cf.div_xy(u_2, v_2, radius="earth") + +# %% +# 8. First we configure the overall plot by +# making the map higher resolution, to show the coastlines of the UK and +# Ireland in greater detail, and changing the colourmap to better reflect +# the data which can be positive or negative, i.e. has 0 as the 'middle' +# value of significance, so should use a diverging colour map. +cfp.mapset(resolution="10m") +cfp.cscale("ncl_default", ncols=21) + +# %% +# 9. Now generate the final plot. Plot the current vectors, noting we had +# to play around with the 'stride' and 'scale' parameter values to adjust +# the vector spacing and size so that the vector field is best represented +# and visible without over-cluttering the plot. Finally we plot the +# divergence as a contour plot without any lines showing. This compound +# plot is saved on one canvas using 'gopen' and 'gclose' to wrap the two +# plotting calls: +cfp.gopen() +cfp.vect(u=u_2, v=v_2, stride=6, scale=3, key_length=1) +cfp.con( + div, + lines=False, + title=( + f"Depth-averaged Irish Sea currents at {chosen_time} with " + "their divergence" + ), +) +cfp.gclose() diff --git a/recipes-docs/source/recipes-source/plot_22_recipe.py b/recipes-docs/source/recipes-source/plot_22_recipe.py new file mode 100644 index 0000000000..20a77eeb86 --- /dev/null +++ b/recipes-docs/source/recipes-source/plot_22_recipe.py @@ -0,0 +1,121 @@ +""" +Plotting a wind rose as a scatter plot +====================================== + +Given a file containing northerly and easterly wind components, we can +calculate the magnitude and bearing of the resultant wind at each point +in the region and plot them using a scatter plot on a polar grid to +create a wind rose representing wind vectors in the given area. +""" + +# %% +# 1. Import cf-python, Dask.array, NumPy, and Matplotlib: +import dask.array as da +import matplotlib.pyplot as plt +import numpy as np + +import cf + +# %% +# 2. Read the field constructs and load the wind speed component fields: + +f = cf.read("~/recipes/data1.nc") +print(f) + +U = f[2].squeeze() # Easterly wind speed component +V = f[3].squeeze() # Northerly wind speed component + + +# %% +# 3. Set a bounding region for the data and discard readings outside of it: + +tl = (41, 72) # (long, lat) of top left of bounding box. +br = (65, 46) # (long, lat) of bottom right of bounding box. + +U_region = U.subspace(X=cf.wi(tl[0], br[0]), Y=cf.wi(br[1], tl[1])) +V_region = V.subspace(X=cf.wi(tl[0], br[0]), Y=cf.wi(br[1], tl[1])) + +# %% +# 4. Select measurements for a specific pressure using the subspace method, +# then use squeeze to remove the size 1 axis: + +U_sub = U_region.subspace(pressure=500.0) +V_sub = V_region.subspace(pressure=500.0) + +U_sub.squeeze(inplace=True) +V_sub.squeeze(inplace=True) + +# %% +# 5. Calculate the magnitude of each resultant vector using Dask's hypot +# function: + +magnitudes = da.hypot(U_sub.data, V_sub.data) + +# %% +# 6. Calculate the angle of the resultant vector (relative to an Easterly ray) +# using Dask's arctan2 function, then convert to a clockwise bearing: + +azimuths = da.arctan2(V_sub.data, U_sub.data) + +bearings = ((np.pi / 2) - azimuths) % (np.pi * 2) + +# %% +# 7. Flatten the two dimensions of each array for plotting with Matplotlib: + +bearings_flattened = da.ravel(bearings) + +magnitudes_flattened = da.ravel(magnitudes) + +# %% +# 8. Draw the scatter plot using Matplotlib: + +plt.figure(figsize=(5, 6)) +# sphinx_gallery_start_ignore +fig = plt.gcf() +# sphinx_gallery_end_ignore +ax = plt.subplot(polar=True) +ax.set_theta_zero_location("N") # Place 0 degrees at the top. +ax.set_theta_direction(-1) # Arrange bearings clockwise around the plot. +# sphinx_gallery_start_ignore +plt.close() +# sphinx_gallery_end_ignore + +# %% +# 9. Draw a scatter plot on the polar plot, using the wind direction bearing as +# the angle and the magnitude of the resultant wind speed as the distance +# from the pole. +ax.scatter(bearings_flattened.compute(), magnitudes_flattened.compute(), s=1.2) +# sphinx_gallery_start_ignore +plt.close() +# sphinx_gallery_end_ignore + +# %% +# 10. Label the axes and add a title. + +# sphinx_gallery_start_ignore +plt.figure(fig) +# sphinx_gallery_end_ignore +plt.title( + f"Wind Rose Scatter Plot\nLat: {br[1]}°-{tl[1]}°, Long: {tl[0]}°-{br[0]}°" +) + +ax.set_xlabel("Bearing [°]") + +ax.set_ylabel("Speed [m/s]", rotation=45, labelpad=30, size=8) + +ax.yaxis.set_label_coords(0.45, 0.45) + +ax.yaxis.set_tick_params(which="both", labelrotation=45, labelsize=8) + +ax.set_rlabel_position(45) +# sphinx_gallery_start_ignore +plt.close() +# sphinx_gallery_end_ignore + +# %% +# 11. Display the plot. + +# sphinx_gallery_start_ignore +plt.figure(fig) +# sphinx_gallery_end_ignore +plt.show() diff --git a/recipes-docs/source/recipes-source/plot_23_recipe.py b/recipes-docs/source/recipes-source/plot_23_recipe.py new file mode 100644 index 0000000000..4ae11c3863 --- /dev/null +++ b/recipes-docs/source/recipes-source/plot_23_recipe.py @@ -0,0 +1,174 @@ +""" +Combining cf and Matplotlib plots in one figure +=============================================== + +Presently, cf-plot has very few exposed interfaces to its Matplotlib and +Cartopy backend. This makes it difficult to combine plots from the three +in one figure, but not impossible. + +A combined cf and Matplotlib plot can be achieved by amending the figure +stored at ``cfp.plotvars.master_plot``, and then redrawing it with the +new subplots. +""" + +# %% +# 1. Import cf-python, cf-plot, Matplotlib, NumPy, and Dask.array: + +# sphinx_gallery_start_ignore +# sphinx_gallery_thumbnail_number = 2 +# sphinx_gallery_end_ignore + +import matplotlib.pyplot as plt +import numpy as np +import dask.array as da + +import cfplot as cfp +import cf + + +# %% +# 2. Read example data field constructs, and set region for our plots: + +f = cf.read(f"~/recipes/data1.nc") + +u = f.select_by_identity("eastward_wind")[0] +v = f.select_by_identity("northward_wind")[0] +t = f.select_by_identity("air_temperature")[0] + +# Subspace to get values for a specified pressure, here 500 mbar +u = u.subspace(pressure=500) +v = v.subspace(pressure=500) +t = t.subspace(pressure=500) + +lonmin, lonmax, latmin, latmax = 10, 120, -30, 30 + +# %% [markdown] +# +# Outlining the figure with cf-plot +# --------------------------------- +# + +# %% +# 1. Set desired dimensions for our final figure: + +rows, cols = 2, 2 + +# %% +# 2. Create a figure of set dimensions with ``cfp.gopen()``, then set the +# position of the cf plot: + + +cfp.gopen(rows, cols) + +pos = 2 # Second position in the figure + +cfp.gpos(pos) +# sphinx_gallery_start_ignore +plt.close() +# sphinx_gallery_end_ignore + +# %% +# 3. Create a simple vector plot: + +cfp.mapset(lonmin=lonmin, lonmax=lonmax, latmin=latmin, latmax=latmax) +cfp.vect(u=u, v=v, key_length=10, scale=120, stride=4) + +# %% [markdown] +# +# Creating our Matplotlib plots +# ----------------------------- +# + +# %% +# 1. Access the newly-created figure: + +fig = cfp.plotvars.master_plot + +# %% +# 2. Reduce fields down to our test data for a wind rose scatter plot: + +# Limit to specific geographic region +u_region = u.subspace(X=cf.wi(lonmin, lonmax), Y=cf.wi(latmin, latmax)) +v_region = v.subspace(X=cf.wi(lonmin, lonmax), Y=cf.wi(latmin, latmax)) +t_region = t.subspace(X=cf.wi(lonmin, lonmax), Y=cf.wi(latmin, latmax)) + +# Remove size 1 axes +u_squeeze = u_region.squeeze() +v_squeeze = v_region.squeeze() +t_squeeze = t_region.squeeze() + +# Flatten to one dimension for plot +u_f = da.ravel(u_squeeze.data) +v_f = da.ravel(v_squeeze.data) +t_f = da.ravel(t_squeeze.data) + +# %% +# 3. Perform calculations to create appropriate plot data: + +mag_f = da.hypot(u_f, v_f) # Wind speed magnitude + +azimuths_f = da.arctan2(v_f, u_f) +rad_f = ((np.pi / 2) - azimuths_f) % (np.pi * 2) # Wind speed bearing + +# Normalise temperature data into a range appropriate for setting point sizes (1-10pt). +temp_scaled = 1 + (t_f - t_f.min()) / (t_f.max() - t_f.min()) * (10 - 1) + +# %% +# 4. Add Matplotlib subplot to our existing cf figure: + +pos = 1 # First position in the figure + +ax = fig.add_subplot(rows, cols, pos, polar=True) +ax.set_theta_zero_location("N") +ax.set_theta_direction(-1) + +ax.scatter( + rad_f.compute(), + mag_f.compute(), + s=temp_scaled.compute(), + c=temp_scaled.compute(), + alpha=0.5, +) + +ax.set_xlabel("Bearing [°]") +ax.set_ylabel("Speed [m/s]", rotation=45, labelpad=30, size=8) +ax.yaxis.set_label_coords(0.45, 0.45) +ax.yaxis.set_tick_params(which="both", labelrotation=45, labelsize=8) +ax.set_rlabel_position(45) + +# %% +# 5. Create and add a third plot, for example: + +x = np.linspace(0, 10, 100) +y = np.sin(x) + +pos = 3 # Third position in the figure + +ax1 = fig.add_subplot(rows, cols, pos) + +ax1.plot(x, y, label="sin(x)") +ax1.legend() + +# %% [markdown] +# +# Drawing the new figure +# ---------------------- +# + +# %% +# 1. Draw final figure: + +fig = plt.figure(fig) +fig.tight_layout() +fig.show() + +# %% [markdown] +# +# Summary +# ------- +# +# In summary, to use other plotting libraries with cf-plot, you must first +# create your figure with cf-plot with placeholders for your other plots, +# then add subplots by accessing the ``cfp.plotvars.master_plot`` object, +# and finally redraw the figure containing the new plots. +