From 99173ad53432b5f38569795902e1f8012bbe316b Mon Sep 17 00:00:00 2001 From: Willi Rath Date: Sun, 20 Dec 2020 12:48:40 +0100 Subject: [PATCH 01/11] Don't mention ball tree in docstring --- src/xoak/accessor.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/xoak/accessor.py b/src/xoak/accessor.py index ba127ed..e8162c0 100644 --- a/src/xoak/accessor.py +++ b/src/xoak/accessor.py @@ -227,7 +227,7 @@ def _get_pos_indexers(self, indices, indexers): def sel( self, indexers: Mapping[Hashable, Any] = None, **indexers_kwargs: Any ) -> Union[xr.Dataset, xr.DataArray]: - """Selection based on a ball tree index. + """Selection based on an tree index. The index must have been already built using `xoak.set_index()`. From ea3eed970de3717402c1af26f853bbca341b9b6e Mon Sep 17 00:00:00 2001 From: Willi Rath Date: Sun, 20 Dec 2020 12:51:31 +0100 Subject: [PATCH 02/11] Expose indexers --- src/xoak/accessor.py | 30 +++++++++++++++++++++++++++--- 1 file changed, 27 insertions(+), 3 deletions(-) diff --git a/src/xoak/accessor.py b/src/xoak/accessor.py index e8162c0..4e87580 100644 --- a/src/xoak/accessor.py +++ b/src/xoak/accessor.py @@ -224,10 +224,10 @@ def _get_pos_indexers(self, indices, indexers): return pos_indexers - def sel( + def get_indexers( self, indexers: Mapping[Hashable, Any] = None, **indexers_kwargs: Any ) -> Union[xr.Dataset, xr.DataArray]: - """Selection based on an tree index. + """Indexers based on an tree index. The index must have been already built using `xoak.set_index()`. @@ -253,11 +253,35 @@ def sel( indices = self._query(indexers) if not isinstance(indices, np.ndarray): - # TODO: remove (see todo below) + # TODO: remove (see TODO in self.sel below) indices = indices.compute() pos_indexers = self._get_pos_indexers(indices, indexers) + return pos_indexers + + def sel( + self, indexers: Mapping[Hashable, Any] = None, **indexers_kwargs: Any + ) -> Union[xr.Dataset, xr.DataArray]: + """Selection based on an tree index. + + The index must have been already built using `xoak.set_index()`. + + It behaves mostly like :meth:`xarray.Dataset.sel` and + :meth:`xarray.DataArray.sel` methods, with some limitations: + + - Orthogonal indexing is not supported + - For vectorized (point-wise) indexing, you need to supply xarray + objects + - Use it for nearest neighbor lookup only (it implicitly + assumes method="nearest") + + This triggers :func:`dask.compute` if the given indexers and/or the index + coordinates are chunked. + + """ + pos_indexers = self.get_indexers(indexers, **indexers_kwargs) + # TODO: issue in xarray. 1-dimensional xarray.Variables are always considered # as OuterIndexer, while we want here VectorizedIndexer # This would also allow lazy selection From 2cc512daa52a3df4f95535a3ead0b4461b1ab5aa Mon Sep 17 00:00:00 2001 From: Willi Rath Date: Sun, 20 Dec 2020 16:11:59 +0100 Subject: [PATCH 03/11] Add two-step test via indexers --- src/xoak/tests/test_scipy_adapters.py | 9 +++++++++ 1 file changed, 9 insertions(+) diff --git a/src/xoak/tests/test_scipy_adapters.py b/src/xoak/tests/test_scipy_adapters.py index a15f119..f12b224 100644 --- a/src/xoak/tests/test_scipy_adapters.py +++ b/src/xoak/tests/test_scipy_adapters.py @@ -13,6 +13,15 @@ def test_scipy_kdtree(xyz_dataset, xyz_indexer, xyz_expected): xr.testing.assert_equal(ds_sel.load(), xyz_expected.load()) +def test_scipy_kdtree_via_indexer(xyz_dataset, xyz_indexer, xyz_expected): + xyz_dataset.xoak.set_index(['x', 'y', 'z'], 'scipy_kdtree') + + indexers = xyz_dataset.xoak.get_indexers(x=xyz_indexer.xx, y=xyz_indexer.yy, z=xyz_indexer.zz) + ds_sel = xyz_dataset.isel(indexers) + + xr.testing.assert_equal(ds_sel.load(), xyz_expected.load()) + + def test_scipy_kdtree_options(): ds = xr.Dataset(coords={'x': ('points', [1, 2]), 'y': ('points', [1, 2])}) From eedb5b9cd9cf8df574cb1e9ab0f89881eb2a7880 Mon Sep 17 00:00:00 2001 From: Willi Rath Date: Sun, 20 Dec 2020 16:31:43 +0100 Subject: [PATCH 04/11] Document new method --- doc/api.rst | 2 ++ 1 file changed, 2 insertions(+) diff --git a/doc/api.rst b/doc/api.rst index 97fa1d7..25a1799 100644 --- a/doc/api.rst +++ b/doc/api.rst @@ -36,6 +36,7 @@ properties listed below. Proper use of this accessor should be like: Dataset.xoak.set_index Dataset.xoak.sel + Dataset.xoak.get_indexers DataArray.xoak -------------- @@ -58,6 +59,7 @@ The accessor above is also registered for :py:class:`xarray.DataArray`. DataArray.xoak.set_index DataArray.xoak.sel + DataArray.xoak.get_indexers Indexes ------- From 407608806619733dad1cf5ba3a3cdeef60595536 Mon Sep 17 00:00:00 2001 From: Willi Rath Date: Tue, 22 Dec 2020 11:08:32 +0100 Subject: [PATCH 05/11] Rename to public query method --- src/xoak/accessor.py | 11 ++--- src/xoak/tests/conftest.py | 57 +++++++++++++++++++++++-- src/xoak/tests/test_s2_adapters.py | 8 ++++ src/xoak/tests/test_scipy_adapters.py | 2 +- src/xoak/tests/test_sklearn_adapters.py | 12 ++++-- 5 files changed, 76 insertions(+), 14 deletions(-) diff --git a/src/xoak/accessor.py b/src/xoak/accessor.py index 4e87580..8395007 100644 --- a/src/xoak/accessor.py +++ b/src/xoak/accessor.py @@ -224,16 +224,13 @@ def _get_pos_indexers(self, indices, indexers): return pos_indexers - def get_indexers( + def query( self, indexers: Mapping[Hashable, Any] = None, **indexers_kwargs: Any ) -> Union[xr.Dataset, xr.DataArray]: - """Indexers based on an tree index. + """Directly query the underlying tree index. The index must have been already built using `xoak.set_index()`. - It behaves mostly like :meth:`xarray.Dataset.sel` and - :meth:`xarray.DataArray.sel` methods, with some limitations: - - Orthogonal indexing is not supported - For vectorized (point-wise) indexing, you need to supply xarray objects @@ -258,7 +255,7 @@ def get_indexers( pos_indexers = self._get_pos_indexers(indices, indexers) - return pos_indexers + return xr.Dataset(pos_indexers) def sel( self, indexers: Mapping[Hashable, Any] = None, **indexers_kwargs: Any @@ -280,7 +277,7 @@ def sel( coordinates are chunked. """ - pos_indexers = self.get_indexers(indexers, **indexers_kwargs) + pos_indexers = self.query(indexers, **indexers_kwargs) # TODO: issue in xarray. 1-dimensional xarray.Variables are always considered # as OuterIndexer, while we want here VectorizedIndexer diff --git a/src/xoak/tests/conftest.py b/src/xoak/tests/conftest.py index c131b68..e985904 100644 --- a/src/xoak/tests/conftest.py +++ b/src/xoak/tests/conftest.py @@ -23,7 +23,11 @@ def indexer_array_lib(request): @pytest.fixture( - params=[(('d1',), (200,)), (('d1', 'd2'), (20, 10)), (('d1', 'd2', 'd3'), (4, 10, 5))], + params=[ + (('d1',), (200,)), + (('d1', 'd2'), (20, 10)), + (('d1', 'd2', 'd3'), (4, 10, 5)), + ], scope='session', ) def dataset_dims_shape(request): @@ -31,7 +35,11 @@ def dataset_dims_shape(request): @pytest.fixture( - params=[(('i1',), (100,)), (('i1', 'i2'), (10, 10)), (('i1', 'i2', 'i3'), (2, 10, 5))], + params=[ + (('i1',), (100,)), + (('i1', 'i2'), (10, 10)), + (('i1', 'i2', 'i3'), (2, 10, 5)), + ], scope='session', ) def indexer_dims_shape(request): @@ -64,6 +72,34 @@ def query_brute_force(dataset, dataset_dims_shape, indexer, indexer_dims_shape, return dataset.isel(indexers=pos_indexers) +def query_brute_force_indexers( + dataset, dataset_dims_shape, indexer, indexer_dims_shape, metric='euclidean' +): + """Find indexers for nearest neighbors using brute-force approach.""" + + # for lat/lon coordinate, assume they are ordered lat, lon!! + X = np.stack([np.ravel(c) for c in indexer.coords.values()]).T + Y = np.stack([np.ravel(c) for c in dataset.coords.values()]).T + + if metric == 'haversine': + X = np.deg2rad(X) + Y = np.deg2rad(Y) + + positions, _ = pairwise_distances_argmin_min(X, Y, metric=metric) + + dataset_dims, dataset_shape = dataset_dims_shape + indexer_dims, indexer_shape = indexer_dims_shape + + u_positions = list(np.unravel_index(positions.ravel(), dataset_shape)) + + pos_indexers = { + dim: xr.Variable(indexer_dims, ind.reshape(indexer_shape)) + for dim, ind in zip(dataset_dims, u_positions) + } + + return xr.Dataset(pos_indexers) + + @pytest.fixture(scope='session') def geo_dataset(dataset_dims_shape, dataset_array_lib): """Dataset with coords lon and lat on a grid of different shapes.""" @@ -93,7 +129,22 @@ def geo_indexer(indexer_dims_shape, indexer_array_lib): @pytest.fixture(scope='session') def geo_expected(geo_dataset, dataset_dims_shape, geo_indexer, indexer_dims_shape): return query_brute_force( - geo_dataset, dataset_dims_shape, geo_indexer, indexer_dims_shape, metric='haversine' + geo_dataset, + dataset_dims_shape, + geo_indexer, + indexer_dims_shape, + metric='haversine', + ) + + +@pytest.fixture(scope='session') +def geo_expected_indexers(geo_dataset, dataset_dims_shape, geo_indexer, indexer_dims_shape): + return query_brute_force_indexers( + geo_dataset, + dataset_dims_shape, + geo_indexer, + indexer_dims_shape, + metric='haversine', ) diff --git a/src/xoak/tests/test_s2_adapters.py b/src/xoak/tests/test_s2_adapters.py index fd47a68..b541261 100644 --- a/src/xoak/tests/test_s2_adapters.py +++ b/src/xoak/tests/test_s2_adapters.py @@ -16,6 +16,14 @@ def test_s2point(geo_dataset, geo_indexer, geo_expected): xr.testing.assert_equal(ds_sel.load(), geo_expected.load()) +def test_s2point_via_query(geo_dataset, geo_indexer, geo_expected): + geo_dataset.xoak.set_index(['lat', 'lon'], 's2point') + ds_indexer = geo_dataset.xoak.query(lat=geo_indexer.latitude, lon=geo_indexer.longitude) + ds_sel = geo_dataset.xoak.isel(ds_indexer) + + xr.testing.assert_equal(ds_sel.load(), geo_expected.load()) + + def test_s2point_sizeof(): ds = xr.Dataset(coords={'lat': ('points', [0.0, 10.0]), 'lon': ('points', [-5.0, 5.0])}) points = np.array([[0.0, -5.0], [10.0, 5.0]]) diff --git a/src/xoak/tests/test_scipy_adapters.py b/src/xoak/tests/test_scipy_adapters.py index f12b224..70776c1 100644 --- a/src/xoak/tests/test_scipy_adapters.py +++ b/src/xoak/tests/test_scipy_adapters.py @@ -16,7 +16,7 @@ def test_scipy_kdtree(xyz_dataset, xyz_indexer, xyz_expected): def test_scipy_kdtree_via_indexer(xyz_dataset, xyz_indexer, xyz_expected): xyz_dataset.xoak.set_index(['x', 'y', 'z'], 'scipy_kdtree') - indexers = xyz_dataset.xoak.get_indexers(x=xyz_indexer.xx, y=xyz_indexer.yy, z=xyz_indexer.zz) + indexers = xyz_dataset.xoak.query(x=xyz_indexer.xx, y=xyz_indexer.yy, z=xyz_indexer.zz) ds_sel = xyz_dataset.isel(indexers) xr.testing.assert_equal(ds_sel.load(), xyz_expected.load()) diff --git a/src/xoak/tests/test_sklearn_adapters.py b/src/xoak/tests/test_sklearn_adapters.py index 0a019c5..f4c5ec6 100644 --- a/src/xoak/tests/test_sklearn_adapters.py +++ b/src/xoak/tests/test_sklearn_adapters.py @@ -38,10 +38,13 @@ def test_sklearn_balltree_options(): assert ds.xoak._index._index_adapter._index_options == {'leaf_size': 10} -def test_sklearn_geo_balltree(geo_dataset, geo_indexer, geo_expected): +def test_sklearn_geo_balltree(geo_dataset, geo_indexer, geo_expected, geo_expected_indexers): geo_dataset.xoak.set_index(['lat', 'lon'], 'sklearn_geo_balltree') - ds_sel = geo_dataset.xoak.sel(lat=geo_indexer.latitude, lon=geo_indexer.longitude) + ds_indexers = geo_dataset.xoak.query(lat=geo_indexer.latitude, lon=geo_indexer.longitude) + xr.testing.assert_equal(ds_indexers.load(), geo_expected_indexers.load()) + + ds_sel = geo_dataset.xoak.sel(lat=geo_indexer.latitude, lon=geo_indexer.longitude) xr.testing.assert_equal(ds_sel.load(), geo_expected.load()) @@ -52,4 +55,7 @@ def test_sklearn_geo_balltree_options(): # sklearn tree classes init options are not exposed as class properties # user-defined metric should be ignored - assert ds.xoak._index._index_adapter._index_options == {'leaf_size': 10, 'metric': 'haversine'} + assert ds.xoak._index._index_adapter._index_options == { + 'leaf_size': 10, + 'metric': 'haversine', + } From 10edb013b1c9a3917b93551311510d17659b1ac8 Mon Sep 17 00:00:00 2001 From: Willi Rath Date: Tue, 22 Dec 2020 11:16:10 +0100 Subject: [PATCH 06/11] Fix s2 test --- src/xoak/tests/test_s2_adapters.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/xoak/tests/test_s2_adapters.py b/src/xoak/tests/test_s2_adapters.py index b541261..2c1ecb5 100644 --- a/src/xoak/tests/test_s2_adapters.py +++ b/src/xoak/tests/test_s2_adapters.py @@ -19,7 +19,7 @@ def test_s2point(geo_dataset, geo_indexer, geo_expected): def test_s2point_via_query(geo_dataset, geo_indexer, geo_expected): geo_dataset.xoak.set_index(['lat', 'lon'], 's2point') ds_indexer = geo_dataset.xoak.query(lat=geo_indexer.latitude, lon=geo_indexer.longitude) - ds_sel = geo_dataset.xoak.isel(ds_indexer) + ds_sel = geo_dataset.isel(ds_indexer) xr.testing.assert_equal(ds_sel.load(), geo_expected.load()) From 1537e29fa3c870dafd96319a8f261190e671ef89 Mon Sep 17 00:00:00 2001 From: Willi Rath Date: Tue, 22 Dec 2020 11:22:04 +0100 Subject: [PATCH 07/11] Add kwargs to adapter base clas query --- src/xoak/index/base.py | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/src/xoak/index/base.py b/src/xoak/index/base.py index 540f202..b0c7f31 100644 --- a/src/xoak/index/base.py +++ b/src/xoak/index/base.py @@ -42,7 +42,7 @@ def build(self, points: np.ndarray) -> Index: raise NotImplementedError() @abc.abstractmethod - def query(self, index: Index, points: np.ndarray) -> Tuple[np.ndarray, np.ndarray]: + def query(self, index: Index, points: np.ndarray, **kwargs) -> Tuple[np.ndarray, np.ndarray]: """Query points/samples, Parameters @@ -53,6 +53,8 @@ def query(self, index: Index, points: np.ndarray) -> Tuple[np.ndarray, np.ndarra Two-dimensional array of points/samples (rows) and their corresponding coordinate labels (columns) to query. + Additional keywords are passed to the query method of the underlying index. + Returns ------- distances : ndarray of shape (n_points) From bc55f899494e8238cf6fc64bcb7e16ffe7e9dfc3 Mon Sep 17 00:00:00 2001 From: Willi Rath Date: Tue, 22 Dec 2020 12:50:50 +0100 Subject: [PATCH 08/11] Add query_kwargs everywhere --- src/xoak/accessor.py | 17 ++++++++++------- src/xoak/index/base.py | 10 ++++++---- src/xoak/index/s2_adapters.py | 2 +- src/xoak/index/scipy_adapters.py | 6 ++++-- src/xoak/index/sklearn_adapters.py | 18 ++++++++++++------ src/xoak/tests/test_index_base.py | 14 +++++++------- src/xoak/tests/test_scipy_adapters.py | 11 +++++++++++ 7 files changed, 51 insertions(+), 27 deletions(-) diff --git a/src/xoak/accessor.py b/src/xoak/accessor.py index 8395007..dd49bb1 100644 --- a/src/xoak/accessor.py +++ b/src/xoak/accessor.py @@ -1,4 +1,4 @@ -from typing import Any, Hashable, Iterable, List, Mapping, Tuple, Type, Union +from typing import Any, Dict, Hashable, Iterable, List, Mapping, Optional, Tuple, Type, Union import numpy as np import xarray as xr @@ -133,12 +133,12 @@ def index(self) -> Union[None, Index, Iterable[Index]]: index_wrappers = dask.compute(*self._index) return [wrp.index for wrp in index_wrappers] - def _query(self, indexers): + def _query(self, indexers, query_kwargs=None): X = coords_to_point_array([indexers[c] for c in self._index_coords]) if isinstance(X, np.ndarray) and isinstance(self._index, XoakIndexWrapper): # directly call index wrapper's query method - res = self._index.query(X) + res = self._index.query(X, query_kwargs) results = res['indices'][:, 0] else: @@ -168,7 +168,7 @@ def _query(self, indexers): shape = (chunk_npoints, 1) for idx in indexes: - dlyd = dask.delayed(idx.query)(chunk) + dlyd = dask.delayed(idx.query)(chunk, query_kwargs=query_kwargs) res_chunk_idx.append( da.from_delayed(dlyd, shape, dtype=XoakIndexWrapper._query_result_dtype) ) @@ -225,7 +225,10 @@ def _get_pos_indexers(self, indices, indexers): return pos_indexers def query( - self, indexers: Mapping[Hashable, Any] = None, **indexers_kwargs: Any + self, + indexers: Mapping[Hashable, Any] = None, + query_kwargs: Optional[Dict] = None, + **indexers_kwargs: Any ) -> Union[xr.Dataset, xr.DataArray]: """Directly query the underlying tree index. @@ -246,8 +249,8 @@ def query( 'The index(es) has/have not been built yet. Call `.xoak.set_index()` first' ) - indexers = either_dict_or_kwargs(indexers, indexers_kwargs, 'xoak.sel') - indices = self._query(indexers) + indexers = either_dict_or_kwargs(indexers, indexers_kwargs, 'xoak.query') + indices = self._query(indexers, query_kwargs) if not isinstance(indices, np.ndarray): # TODO: remove (see TODO in self.sel below) diff --git a/src/xoak/index/base.py b/src/xoak/index/base.py index b0c7f31..e8e5970 100644 --- a/src/xoak/index/base.py +++ b/src/xoak/index/base.py @@ -1,7 +1,7 @@ import abc import warnings from contextlib import suppress -from typing import Any, Dict, List, Mapping, Tuple, Type, TypeVar, Union +from typing import Any, Dict, List, Mapping, Optional, Tuple, Type, TypeVar, Union import numpy as np @@ -42,7 +42,9 @@ def build(self, points: np.ndarray) -> Index: raise NotImplementedError() @abc.abstractmethod - def query(self, index: Index, points: np.ndarray, **kwargs) -> Tuple[np.ndarray, np.ndarray]: + def query( + self, index: Index, points: np.ndarray, query_kwargs: Any + ) -> Tuple[np.ndarray, np.ndarray]: """Query points/samples, Parameters @@ -233,8 +235,8 @@ def __init__( def index(self): return self._index - def query(self, points: np.ndarray) -> np.ndarray: - distances, positions = self._index_adapter.query(self._index, points) + def query(self, points: np.ndarray, query_kwargs: Optional[Dict]) -> np.ndarray: + distances, positions = self._index_adapter.query(self._index, points, query_kwargs) result = np.empty(shape=points.shape[0], dtype=self._query_result_dtype) result['distances'] = distances.ravel().astype(np.double) diff --git a/src/xoak/index/s2_adapters.py b/src/xoak/index/s2_adapters.py index f971eca..101b6b2 100644 --- a/src/xoak/index/s2_adapters.py +++ b/src/xoak/index/s2_adapters.py @@ -23,7 +23,7 @@ def build(self, points): self._points_nbytes = points.nbytes return S2PointIndex(points) - def query(self, s2index, points): + def query(self, s2index, points, query_kwargs=None): return s2index.query(points) def __sizeof__(self): diff --git a/src/xoak/index/scipy_adapters.py b/src/xoak/index/scipy_adapters.py index eb35d92..4b2bf03 100644 --- a/src/xoak/index/scipy_adapters.py +++ b/src/xoak/index/scipy_adapters.py @@ -13,5 +13,7 @@ def __init__(self, **kwargs): def build(self, points): return cKDTree(points, **self.index_options) - def query(self, kdtree, points): - return kdtree.query(points) + def query(self, kdtree, points, query_kwargs=None): + if query_kwargs is None: + query_kwargs = {} + return kdtree.query(points, **query_kwargs) diff --git a/src/xoak/index/sklearn_adapters.py b/src/xoak/index/sklearn_adapters.py index e913132..9ecbf22 100644 --- a/src/xoak/index/sklearn_adapters.py +++ b/src/xoak/index/sklearn_adapters.py @@ -14,8 +14,10 @@ def __init__(self, **kwargs): def build(self, points): return KDTree(points, **self._index_options) - def query(self, kdtree, points): - return kdtree.query(points) + def query(self, kdtree, points, query_kwargs=None): + if query_kwargs is None: + query_kwargs = {} + return kdtree.query(points, **query_kwargs) @register_default('sklearn_balltree') @@ -28,8 +30,10 @@ def __init__(self, **kwargs): def build(self, points): return BallTree(points, **self._index_options) - def query(self, btree, points): - return btree.query(points) + def query(self, btree, points, query_kwargs=None): + if query_kwargs is None: + query_kwargs = {} + return btree.query(points, **query_kwargs) @register_default('sklearn_geo_balltree') @@ -54,5 +58,7 @@ def __init__(self, **kwargs): def build(self, points): return BallTree(np.deg2rad(points), **self._index_options) - def query(self, btree, points): - return btree.query(np.deg2rad(points)) + def query(self, btree, points, query_kwargs=None): + if query_kwargs is None: + query_kwargs = {} + return btree.query(np.deg2rad(points), **query_kwargs) diff --git a/src/xoak/tests/test_index_base.py b/src/xoak/tests/test_index_base.py index 2974921..4a29710 100644 --- a/src/xoak/tests/test_index_base.py +++ b/src/xoak/tests/test_index_base.py @@ -16,7 +16,7 @@ def __init__(self, points, option=1): self.points = points self.option = option - def query(self, points): + def query(self, points, query_kwargs): distances = np.zeros((points.shape[0])) indices = np.ones((points.shape[0])) @@ -30,8 +30,8 @@ def __init__(self, **kwargs): def build(self, points): return DummyIndex(points, **self.index_kwargs) - def query(self, index, points): - return index.query(points) + def query(self, index, points, query_kwargs): + return index.query(points, query_kwargs) def test_index_adapter_base(): @@ -39,8 +39,8 @@ class IndexAdapterSubclass(IndexAdapter): def build(self, points): return super().build(points) - def query(self, index, points): - return super().query(index, points) + def query(self, index, points, query_kwargs): + return super().query(index, points, query_kwargs) adapter = IndexAdapterSubclass() @@ -48,7 +48,7 @@ def query(self, index, points): adapter.build(np.zeros((10, 2))) with pytest.raises(NotImplementedError): - adapter.query(None, np.zeros((10, 2))) + adapter.query(None, np.zeros((10, 2)), None) def test_index_registery_constructor(): @@ -130,7 +130,7 @@ def test_xoak_index_wrapper(): assert isinstance(wrapper2.index, DummyIndex) assert wrapper2.index.option == 1 - results = wrapper.query(np.zeros((5, 2))).ravel() + results = wrapper.query(np.zeros((5, 2)), query_kwargs=None).ravel() assert results['distances'].dtype == np.double assert results['indices'].dtype == np.intp diff --git a/src/xoak/tests/test_scipy_adapters.py b/src/xoak/tests/test_scipy_adapters.py index 70776c1..81f3080 100644 --- a/src/xoak/tests/test_scipy_adapters.py +++ b/src/xoak/tests/test_scipy_adapters.py @@ -28,3 +28,14 @@ def test_scipy_kdtree_options(): ds.xoak.set_index(['x', 'y'], 'scipy_kdtree', leafsize=10) assert ds.xoak.index.leafsize == 10 + + +def test_scipy_kdtree_query_k_points_explicit(xyz_dataset, xyz_indexer, xyz_expected): + xyz_dataset.xoak.set_index(['x', 'y', 'z'], 'scipy_kdtree') + + indexers = xyz_dataset.xoak.query( + x=xyz_indexer.xx, y=xyz_indexer.yy, z=xyz_indexer.zz, query_kwargs={'k': 1} + ) + ds_sel = xyz_dataset.isel(indexers) + + xr.testing.assert_equal(ds_sel.load(), xyz_expected.load()) From 6458f50b5d240e6d6c7cf0f56c47498ee0961b1d Mon Sep 17 00:00:00 2001 From: Willi Rath Date: Mon, 1 Mar 2021 11:33:58 +0100 Subject: [PATCH 09/11] Revert kwarg exposure --- src/xoak/accessor.py | 17 +++++++---------- src/xoak/index/base.py | 12 ++++-------- src/xoak/index/s2_adapters.py | 2 +- src/xoak/index/scipy_adapters.py | 6 ++---- src/xoak/index/sklearn_adapters.py | 18 ++++++------------ src/xoak/tests/test_index_base.py | 14 +++++++------- src/xoak/tests/test_scipy_adapters.py | 11 ----------- 7 files changed, 27 insertions(+), 53 deletions(-) diff --git a/src/xoak/accessor.py b/src/xoak/accessor.py index dd49bb1..8395007 100644 --- a/src/xoak/accessor.py +++ b/src/xoak/accessor.py @@ -1,4 +1,4 @@ -from typing import Any, Dict, Hashable, Iterable, List, Mapping, Optional, Tuple, Type, Union +from typing import Any, Hashable, Iterable, List, Mapping, Tuple, Type, Union import numpy as np import xarray as xr @@ -133,12 +133,12 @@ def index(self) -> Union[None, Index, Iterable[Index]]: index_wrappers = dask.compute(*self._index) return [wrp.index for wrp in index_wrappers] - def _query(self, indexers, query_kwargs=None): + def _query(self, indexers): X = coords_to_point_array([indexers[c] for c in self._index_coords]) if isinstance(X, np.ndarray) and isinstance(self._index, XoakIndexWrapper): # directly call index wrapper's query method - res = self._index.query(X, query_kwargs) + res = self._index.query(X) results = res['indices'][:, 0] else: @@ -168,7 +168,7 @@ def _query(self, indexers, query_kwargs=None): shape = (chunk_npoints, 1) for idx in indexes: - dlyd = dask.delayed(idx.query)(chunk, query_kwargs=query_kwargs) + dlyd = dask.delayed(idx.query)(chunk) res_chunk_idx.append( da.from_delayed(dlyd, shape, dtype=XoakIndexWrapper._query_result_dtype) ) @@ -225,10 +225,7 @@ def _get_pos_indexers(self, indices, indexers): return pos_indexers def query( - self, - indexers: Mapping[Hashable, Any] = None, - query_kwargs: Optional[Dict] = None, - **indexers_kwargs: Any + self, indexers: Mapping[Hashable, Any] = None, **indexers_kwargs: Any ) -> Union[xr.Dataset, xr.DataArray]: """Directly query the underlying tree index. @@ -249,8 +246,8 @@ def query( 'The index(es) has/have not been built yet. Call `.xoak.set_index()` first' ) - indexers = either_dict_or_kwargs(indexers, indexers_kwargs, 'xoak.query') - indices = self._query(indexers, query_kwargs) + indexers = either_dict_or_kwargs(indexers, indexers_kwargs, 'xoak.sel') + indices = self._query(indexers) if not isinstance(indices, np.ndarray): # TODO: remove (see TODO in self.sel below) diff --git a/src/xoak/index/base.py b/src/xoak/index/base.py index e8e5970..540f202 100644 --- a/src/xoak/index/base.py +++ b/src/xoak/index/base.py @@ -1,7 +1,7 @@ import abc import warnings from contextlib import suppress -from typing import Any, Dict, List, Mapping, Optional, Tuple, Type, TypeVar, Union +from typing import Any, Dict, List, Mapping, Tuple, Type, TypeVar, Union import numpy as np @@ -42,9 +42,7 @@ def build(self, points: np.ndarray) -> Index: raise NotImplementedError() @abc.abstractmethod - def query( - self, index: Index, points: np.ndarray, query_kwargs: Any - ) -> Tuple[np.ndarray, np.ndarray]: + def query(self, index: Index, points: np.ndarray) -> Tuple[np.ndarray, np.ndarray]: """Query points/samples, Parameters @@ -55,8 +53,6 @@ def query( Two-dimensional array of points/samples (rows) and their corresponding coordinate labels (columns) to query. - Additional keywords are passed to the query method of the underlying index. - Returns ------- distances : ndarray of shape (n_points) @@ -235,8 +231,8 @@ def __init__( def index(self): return self._index - def query(self, points: np.ndarray, query_kwargs: Optional[Dict]) -> np.ndarray: - distances, positions = self._index_adapter.query(self._index, points, query_kwargs) + def query(self, points: np.ndarray) -> np.ndarray: + distances, positions = self._index_adapter.query(self._index, points) result = np.empty(shape=points.shape[0], dtype=self._query_result_dtype) result['distances'] = distances.ravel().astype(np.double) diff --git a/src/xoak/index/s2_adapters.py b/src/xoak/index/s2_adapters.py index 101b6b2..f971eca 100644 --- a/src/xoak/index/s2_adapters.py +++ b/src/xoak/index/s2_adapters.py @@ -23,7 +23,7 @@ def build(self, points): self._points_nbytes = points.nbytes return S2PointIndex(points) - def query(self, s2index, points, query_kwargs=None): + def query(self, s2index, points): return s2index.query(points) def __sizeof__(self): diff --git a/src/xoak/index/scipy_adapters.py b/src/xoak/index/scipy_adapters.py index 4b2bf03..eb35d92 100644 --- a/src/xoak/index/scipy_adapters.py +++ b/src/xoak/index/scipy_adapters.py @@ -13,7 +13,5 @@ def __init__(self, **kwargs): def build(self, points): return cKDTree(points, **self.index_options) - def query(self, kdtree, points, query_kwargs=None): - if query_kwargs is None: - query_kwargs = {} - return kdtree.query(points, **query_kwargs) + def query(self, kdtree, points): + return kdtree.query(points) diff --git a/src/xoak/index/sklearn_adapters.py b/src/xoak/index/sklearn_adapters.py index 9ecbf22..e913132 100644 --- a/src/xoak/index/sklearn_adapters.py +++ b/src/xoak/index/sklearn_adapters.py @@ -14,10 +14,8 @@ def __init__(self, **kwargs): def build(self, points): return KDTree(points, **self._index_options) - def query(self, kdtree, points, query_kwargs=None): - if query_kwargs is None: - query_kwargs = {} - return kdtree.query(points, **query_kwargs) + def query(self, kdtree, points): + return kdtree.query(points) @register_default('sklearn_balltree') @@ -30,10 +28,8 @@ def __init__(self, **kwargs): def build(self, points): return BallTree(points, **self._index_options) - def query(self, btree, points, query_kwargs=None): - if query_kwargs is None: - query_kwargs = {} - return btree.query(points, **query_kwargs) + def query(self, btree, points): + return btree.query(points) @register_default('sklearn_geo_balltree') @@ -58,7 +54,5 @@ def __init__(self, **kwargs): def build(self, points): return BallTree(np.deg2rad(points), **self._index_options) - def query(self, btree, points, query_kwargs=None): - if query_kwargs is None: - query_kwargs = {} - return btree.query(np.deg2rad(points), **query_kwargs) + def query(self, btree, points): + return btree.query(np.deg2rad(points)) diff --git a/src/xoak/tests/test_index_base.py b/src/xoak/tests/test_index_base.py index 4a29710..2974921 100644 --- a/src/xoak/tests/test_index_base.py +++ b/src/xoak/tests/test_index_base.py @@ -16,7 +16,7 @@ def __init__(self, points, option=1): self.points = points self.option = option - def query(self, points, query_kwargs): + def query(self, points): distances = np.zeros((points.shape[0])) indices = np.ones((points.shape[0])) @@ -30,8 +30,8 @@ def __init__(self, **kwargs): def build(self, points): return DummyIndex(points, **self.index_kwargs) - def query(self, index, points, query_kwargs): - return index.query(points, query_kwargs) + def query(self, index, points): + return index.query(points) def test_index_adapter_base(): @@ -39,8 +39,8 @@ class IndexAdapterSubclass(IndexAdapter): def build(self, points): return super().build(points) - def query(self, index, points, query_kwargs): - return super().query(index, points, query_kwargs) + def query(self, index, points): + return super().query(index, points) adapter = IndexAdapterSubclass() @@ -48,7 +48,7 @@ def query(self, index, points, query_kwargs): adapter.build(np.zeros((10, 2))) with pytest.raises(NotImplementedError): - adapter.query(None, np.zeros((10, 2)), None) + adapter.query(None, np.zeros((10, 2))) def test_index_registery_constructor(): @@ -130,7 +130,7 @@ def test_xoak_index_wrapper(): assert isinstance(wrapper2.index, DummyIndex) assert wrapper2.index.option == 1 - results = wrapper.query(np.zeros((5, 2)), query_kwargs=None).ravel() + results = wrapper.query(np.zeros((5, 2))).ravel() assert results['distances'].dtype == np.double assert results['indices'].dtype == np.intp diff --git a/src/xoak/tests/test_scipy_adapters.py b/src/xoak/tests/test_scipy_adapters.py index 81f3080..70776c1 100644 --- a/src/xoak/tests/test_scipy_adapters.py +++ b/src/xoak/tests/test_scipy_adapters.py @@ -28,14 +28,3 @@ def test_scipy_kdtree_options(): ds.xoak.set_index(['x', 'y'], 'scipy_kdtree', leafsize=10) assert ds.xoak.index.leafsize == 10 - - -def test_scipy_kdtree_query_k_points_explicit(xyz_dataset, xyz_indexer, xyz_expected): - xyz_dataset.xoak.set_index(['x', 'y', 'z'], 'scipy_kdtree') - - indexers = xyz_dataset.xoak.query( - x=xyz_indexer.xx, y=xyz_indexer.yy, z=xyz_indexer.zz, query_kwargs={'k': 1} - ) - ds_sel = xyz_dataset.isel(indexers) - - xr.testing.assert_equal(ds_sel.load(), xyz_expected.load()) From a7e4c818a33d27c37dbfadc5977e901ebc0a9204 Mon Sep 17 00:00:00 2001 From: Willi Rath Date: Mon, 10 May 2021 21:15:04 +0200 Subject: [PATCH 10/11] Add query example --- doc/examples/index.rst | 1 + doc/examples/query.ipynb | 382 +++++++++++++++++++++++++++++++++++++++ 2 files changed, 383 insertions(+) create mode 100644 doc/examples/query.ipynb diff --git a/doc/examples/index.rst b/doc/examples/index.rst index bc48c62..a4035f7 100644 --- a/doc/examples/index.rst +++ b/doc/examples/index.rst @@ -7,3 +7,4 @@ Examples introduction dask_support custom_indexes + query diff --git a/doc/examples/query.ipynb b/doc/examples/query.ipynb new file mode 100644 index 0000000..fe728c1 --- /dev/null +++ b/doc/examples/query.ipynb @@ -0,0 +1,382 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Querying for indices\n", + "\n", + "This notebook briefly shows how to use Xoak's `query()` method to build more complicated queries." + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "metadata": {}, + "outputs": [], + "source": [ + "import numpy as np\n", + "import xarray as xr\n", + "import xoak\n", + "\n", + "xr.set_options(display_style='text');" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Load example dataset" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [ + { + "data": { + "text/html": [ + "
<xarray.Dataset>\n",
+       "Dimensions:  (lat: 25, lon: 53, time: 2920)\n",
+       "Coordinates:\n",
+       "  * lat      (lat) float32 75.0 72.5 70.0 67.5 65.0 ... 25.0 22.5 20.0 17.5 15.0\n",
+       "  * lon      (lon) float32 200.0 202.5 205.0 207.5 ... 322.5 325.0 327.5 330.0\n",
+       "  * time     (time) datetime64[ns] 2013-01-01 ... 2014-12-31T18:00:00\n",
+       "    lat2d    (lat, lon) float32 75.0 75.0 75.0 75.0 75.0 ... 15.0 15.0 15.0 15.0\n",
+       "    lon2d    (lat, lon) float32 200.0 202.5 205.0 207.5 ... 325.0 327.5 330.0\n",
+       "Data variables:\n",
+       "    air      (time, lat, lon) float32 241.2 242.5 243.5 ... 296.5 296.2 295.7\n",
+       "Attributes:\n",
+       "    Conventions:  COARDS\n",
+       "    title:        4x daily NMC reanalysis (1948)\n",
+       "    description:  Data is from NMC initialized reanalysis\\n(4x/day).  These a...\n",
+       "    platform:     Model\n",
+       "    references:   http://www.esrl.noaa.gov/psd/data/gridded/data.ncep.reanaly...
" + ], + "text/plain": [ + "\n", + "Dimensions: (lat: 25, lon: 53, time: 2920)\n", + "Coordinates:\n", + " * lat (lat) float32 75.0 72.5 70.0 67.5 65.0 ... 25.0 22.5 20.0 17.5 15.0\n", + " * lon (lon) float32 200.0 202.5 205.0 207.5 ... 322.5 325.0 327.5 330.0\n", + " * time (time) datetime64[ns] 2013-01-01 ... 2014-12-31T18:00:00\n", + " lat2d (lat, lon) float32 75.0 75.0 75.0 75.0 75.0 ... 15.0 15.0 15.0 15.0\n", + " lon2d (lat, lon) float32 200.0 202.5 205.0 207.5 ... 325.0 327.5 330.0\n", + "Data variables:\n", + " air (time, lat, lon) float32 241.2 242.5 243.5 ... 296.5 296.2 295.7\n", + "Attributes:\n", + " Conventions: COARDS\n", + " title: 4x daily NMC reanalysis (1948)\n", + " description: Data is from NMC initialized reanalysis\\n(4x/day). These a...\n", + " platform: Model\n", + " references: http://www.esrl.noaa.gov/psd/data/gridded/data.ncep.reanaly..." + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "ds = xr.tutorial.load_dataset(\"air_temperature\")\n", + "\n", + "# create artificial 2d coords\n", + "ds.coords[\"lat2d\"], ds.coords[\"lon2d\"] = xr.broadcast(\n", + " ds.coords[\"lat\"], ds.coords[\"lon\"]\n", + ")\n", + "\n", + "display(ds)" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": {}, + "outputs": [], + "source": [ + "ds.xoak.set_index(['lat2d', 'lon2d'], 'sklearn_geo_balltree')" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Create trajectory dataset" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "metadata": {}, + "outputs": [ + { + "data": { + "text/html": [ + "
<xarray.Dataset>\n",
+       "Dimensions:      (along_track: 23)\n",
+       "Coordinates:\n",
+       "  * along_track  (along_track) int64 0 1 2 3 4 5 6 7 ... 15 16 17 18 19 20 21 22\n",
+       "Data variables:\n",
+       "    latitude     (along_track) float64 75.0 72.27 69.55 ... 20.45 17.73 15.0\n",
+       "    longitude    (along_track) float64 200.0 205.9 211.8 ... 318.2 324.1 330.0
" + ], + "text/plain": [ + "\n", + "Dimensions: (along_track: 23)\n", + "Coordinates:\n", + " * along_track (along_track) int64 0 1 2 3 4 5 6 7 ... 15 16 17 18 19 20 21 22\n", + "Data variables:\n", + " latitude (along_track) float64 75.0 72.27 69.55 ... 20.45 17.73 15.0\n", + " longitude (along_track) float64 200.0 205.9 211.8 ... 318.2 324.1 330.0" + ] + }, + "execution_count": 4, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "ds_trajectory = xr.Dataset(\n", + " {\n", + " 'latitude': ('along_track', np.linspace(75.0, 15.0, 23)),\n", + " 'longitude': ('along_track', np.linspace(200, 330, 23)),\n", + " },\n", + " coords={\"along_track\": (\"along_track\", np.arange(23))},\n", + ")\n", + "\n", + "ds_trajectory" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Select and check for distances to grid points" + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "metadata": {}, + "outputs": [ + { + "data": { + "text/html": [ + "
<xarray.Dataset>\n",
+       "Dimensions:  (along_track: 23, time: 2920)\n",
+       "Coordinates:\n",
+       "    lat      (along_track) float32 75.0 72.5 70.0 67.5 ... 22.5 20.0 17.5 15.0\n",
+       "    lon      (along_track) float32 200.0 205.0 212.5 217.5 ... 317.5 325.0 330.0\n",
+       "  * time     (time) datetime64[ns] 2013-01-01 ... 2014-12-31T18:00:00\n",
+       "    lat2d    (along_track) float32 75.0 72.5 70.0 67.5 ... 22.5 20.0 17.5 15.0\n",
+       "    lon2d    (along_track) float32 200.0 205.0 212.5 217.5 ... 317.5 325.0 330.0\n",
+       "Dimensions without coordinates: along_track\n",
+       "Data variables:\n",
+       "    air      (time, along_track) float32 241.2 244.7 244.4 ... 296.7 295.7 295.7\n",
+       "Attributes:\n",
+       "    Conventions:  COARDS\n",
+       "    title:        4x daily NMC reanalysis (1948)\n",
+       "    description:  Data is from NMC initialized reanalysis\\n(4x/day).  These a...\n",
+       "    platform:     Model\n",
+       "    references:   http://www.esrl.noaa.gov/psd/data/gridded/data.ncep.reanaly...
" + ], + "text/plain": [ + "\n", + "Dimensions: (along_track: 23, time: 2920)\n", + "Coordinates:\n", + " lat (along_track) float32 75.0 72.5 70.0 67.5 ... 22.5 20.0 17.5 15.0\n", + " lon (along_track) float32 200.0 205.0 212.5 217.5 ... 317.5 325.0 330.0\n", + " * time (time) datetime64[ns] 2013-01-01 ... 2014-12-31T18:00:00\n", + " lat2d (along_track) float32 75.0 72.5 70.0 67.5 ... 22.5 20.0 17.5 15.0\n", + " lon2d (along_track) float32 200.0 205.0 212.5 217.5 ... 317.5 325.0 330.0\n", + "Dimensions without coordinates: along_track\n", + "Data variables:\n", + " air (time, along_track) float32 241.2 244.7 244.4 ... 296.7 295.7 295.7\n", + "Attributes:\n", + " Conventions: COARDS\n", + " title: 4x daily NMC reanalysis (1948)\n", + " description: Data is from NMC initialized reanalysis\\n(4x/day). These a...\n", + " platform: Model\n", + " references: http://www.esrl.noaa.gov/psd/data/gridded/data.ncep.reanaly..." + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "# ds.xoak.sel() does the following two steps\n", + "\n", + "indices = ds.xoak.query(\n", + " lat2d=ds_trajectory.latitude,\n", + " lon2d=ds_trajectory.longitude,\n", + ")\n", + "selected = ds.isel(indices)\n", + "\n", + "display(selected)" + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "metadata": {}, + "outputs": [], + "source": [ + "selected.coords[\"along_track\"] = ds_trajectory.coords[\"along_track\"]" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## We can now easily calculated the mismatch" + ] + }, + { + "cell_type": "code", + "execution_count": 7, + "metadata": {}, + "outputs": [], + "source": [ + "selected[\"diff_lat\"] = ds_trajectory[\"latitude\"] - selected[\"lat\"]\n", + "selected[\"diff_lon\"] = ds_trajectory[\"longitude\"] - selected[\"lon\"]" + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "metadata": {}, + "outputs": [ + { + "data": { + "text/html": [ + "
<xarray.Dataset>\n",
+       "Dimensions:      (along_track: 23, time: 2920)\n",
+       "Coordinates:\n",
+       "    lat          (along_track) float32 75.0 72.5 70.0 67.5 ... 20.0 17.5 15.0\n",
+       "    lon          (along_track) float32 200.0 205.0 212.5 ... 317.5 325.0 330.0\n",
+       "  * time         (time) datetime64[ns] 2013-01-01 ... 2014-12-31T18:00:00\n",
+       "    lat2d        (along_track) float32 75.0 72.5 70.0 67.5 ... 20.0 17.5 15.0\n",
+       "    lon2d        (along_track) float32 200.0 205.0 212.5 ... 317.5 325.0 330.0\n",
+       "  * along_track  (along_track) int64 0 1 2 3 4 5 6 7 ... 15 16 17 18 19 20 21 22\n",
+       "Data variables:\n",
+       "    air          (time, along_track) float32 241.2 244.7 244.4 ... 295.7 295.7\n",
+       "    diff_lat     (along_track) float64 0.0 -0.2273 -0.4545 ... 0.4545 0.2273 0.0\n",
+       "    diff_lon     (along_track) float64 0.0 0.9091 -0.6818 ... 0.6818 -0.9091 0.0\n",
+       "Attributes:\n",
+       "    Conventions:  COARDS\n",
+       "    title:        4x daily NMC reanalysis (1948)\n",
+       "    description:  Data is from NMC initialized reanalysis\\n(4x/day).  These a...\n",
+       "    platform:     Model\n",
+       "    references:   http://www.esrl.noaa.gov/psd/data/gridded/data.ncep.reanaly...
" + ], + "text/plain": [ + "\n", + "Dimensions: (along_track: 23, time: 2920)\n", + "Coordinates:\n", + " lat (along_track) float32 75.0 72.5 70.0 67.5 ... 20.0 17.5 15.0\n", + " lon (along_track) float32 200.0 205.0 212.5 ... 317.5 325.0 330.0\n", + " * time (time) datetime64[ns] 2013-01-01 ... 2014-12-31T18:00:00\n", + " lat2d (along_track) float32 75.0 72.5 70.0 67.5 ... 20.0 17.5 15.0\n", + " lon2d (along_track) float32 200.0 205.0 212.5 ... 317.5 325.0 330.0\n", + " * along_track (along_track) int64 0 1 2 3 4 5 6 7 ... 15 16 17 18 19 20 21 22\n", + "Data variables:\n", + " air (time, along_track) float32 241.2 244.7 244.4 ... 295.7 295.7\n", + " diff_lat (along_track) float64 0.0 -0.2273 -0.4545 ... 0.4545 0.2273 0.0\n", + " diff_lon (along_track) float64 0.0 0.9091 -0.6818 ... 0.6818 -0.9091 0.0\n", + "Attributes:\n", + " Conventions: COARDS\n", + " title: 4x daily NMC reanalysis (1948)\n", + " description: Data is from NMC initialized reanalysis\\n(4x/day). These a...\n", + " platform: Model\n", + " references: http://www.esrl.noaa.gov/psd/data/gridded/data.ncep.reanaly..." + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "display(selected)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Visualize" + ] + }, + { + "cell_type": "code", + "execution_count": 9, + "metadata": {}, + "outputs": [], + "source": [ + "from matplotlib import pyplot as plt" + ] + }, + { + "cell_type": "code", + "execution_count": 10, + "metadata": {}, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjQAAADQCAYAAAAQwfu+AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAsTAAALEwEAmpwYAABH5ElEQVR4nO3deXhcZfn/8fdn0iTdW0jLWkrZZC1LKUtYA0VAQAoWEEVZtaCgoILCF0UUEZSfiooKVcEii6AVWQRZAsNmWFpAyiqLRQoVaGhLS2kmmbl/f5wn7TSd5WSdmeR+Xde5Zs5+z8mZk2eeVWaGc84551wlS5Q6AOecc8657vIEjXPOOecqnidonHPOOVfxPEHjnHPOuYrnCRrnnHPOVTxP0DjnnHOu4nmCZgCTNE/SAeV+zC7E0CBpfoH1JmnzPo7pSknf6ctzOheXpD9I+oGkvSW9nLV8S0lPS1oq6auShki6XdISSX8uZczOdeQJmiyl/mdc6vOXi/aHa6nj6CpJJ0p6JHuZmZ1mZheVKibn4jCzh81sy6xF3wSSZjbCzH4BHAWsC9SZ2dElCdK5PDxB04MkVZU6Buec60EbA893mP+3mbV19kCSBvVYVM7l4AmaQNIfgfHA7ZKWSfpmWP5nSf8LWawPSdo2a58/SPqNpDslfQjsJ2lSVhbtnyXdlJ3bIOkwSc9IWizpn5K2L3T+DjGuJekOSe9JWhTej8tan5R0kaRHw/nvkTQma/3nJb0hqVnS+UWuxyGSXgjHeUvS2cU+Q45jJCSdK+m1cM6bJa2dtX6vsP9iSW+GnI3pwHHAN8N1uD1su4GkWeGz/0fSV7OOMyT8LRZJegHYpdBn6xDjKEnXhuO+IenbkhJZ678o6cVwHV6QNCksb/9c7cuPDMu3Bq4E6kP8i8Py1XKdwnFflfS+pNskbZC1ziSdJumV8Jl+JUlh3eaSHgz340JJN8X9rM61k7STpKfC/XsTMDgsX1lcK+l+YD/ginAv3whcAHw6zJ8Stjs5fEcWSbpb0sZZ5zFJp0t6BXglLMv7/FCUS322pGfDPX6TpMFZ66eGfT8I37+Dw/JRkn4vaUF4Xv1A4Qemf2cGEDPzKUzAPOCADstOBkYAtcDlwDNZ6/4ALAH2JEocjgTeAM4EqoFPASngB2H7ScC7wG5AFXBCOGdtvvN3iKUOmAYMDTH9Gfhb1vok8BrwMWBImL80rNsGWAbsEz7LT4G2fOcDFgB7h/drAZM6+xmAs4DHgHHhnFcBN4Z144GlwGfCtaoDdsy6rj/IiiUBzCF6mNYAmwKvAweF9ZcCDwNrAxsBzwHzC1xHAzYP768Fbg3XcwLwb+CUsO5o4C2iBJKAzYGNs9ZtEGL7NPAhsH5YdyLwSIdzrvxMwP7AwnAta4FfAg91iO8OYHS4Tu8BB4d1NwLnh/MOBvYq9ffGp8qawnfoDeBr4bt3FNAK/ABoyP7uED1DvpA1fyFwXdb8EcCrwNbAIODbwD+z1htwb/huDiHe8+OJ8N1aG3gROC2s25XoefvxcP9vCGwV1v2N6PkyDFgnHOPUsM6/MwNkKnkA5TRRPEExOnxBR4X5PwDXZq3fJ/wDVNayR7L+kf0GuKjDMV8G9o1z/hzx7AgsyppPAt/Omv8y8I/w/gLgT1nrhhEltvIlaP4LnAqM7LA89mcID6MpWdutT/TgHAScB9yS59x/YPUEzW7Afztscx5wTXj/OuEffpifTowETXigtgDbZK07lajOAMDdwJkx/xbPAFPD+xMpnKD5PfDjrHXDw3WZkBXfXlnrbwbODe+vBWYA40r1PfGpsqfwnHq7w3Pqn3QtQXMX4QdAmE8Ay1mV8Ddg/6z1cZ4fn8ta92PgyvD+KuBnOT7PuuF7PCRr2WeAB8J7/84MkMmLnAqQVCXp0pC1+QHRlw1gTNZmb2a93wB4y8K3KMf6jYFvhKzWxaE4YqOwX5x4hkq6KhSNfAA8BIzW6nV3/pf1fjnRP8v22FbGYmYfAs0FTjcNOAR4I2TX1nfhM2wM3JK13YtAmugBtBFRblIcGwMbdDjn/4XjrPHZiH59xjGGVb9Ws/fdMLzPG6Ok47OyzRcD27H6fVHIBtnnNLNlRH+LDbO2yfd3/CZRbtETkp6XdHLMczrXLtdzKu53pqONgZ9nfQ/eJ7o/s+/lzj4D8937+b6PGxPlNC3IOuZVRDk14N+ZAcMraa2u49DjnwWmAgcQJWZGAYuIvhy59lkAbChJWQ+L7C/hm8DFZnZxzPN39A1gS2A3M/ufpB2BpzvEk88ComxhIEocERXz5A7E7ElgqqRq4AyiXIKNYnyGbG8CJ5vZox1XSHqTKAs55+lzHOc/ZrZFnu0XhNjaKy+OjxEbRMU+rUQPxBey9n0r67ybddwp1BH4LTAFaDKztKRnWPV3KPZ3fDucs/14w4j+Fm/l3aP9wGb/A74Y9tsLuE/SQ2b2arF9nQtyPafGE/8HRrb258H1Bbbp+AMv7vMj17nW+D6G5S3AGMtRWdm/MwOH59Cs7h2i+hntRhB9UZqJ6q38sMj+TUQ5EGdIGiRpKqv/0/4tcJqk3RQZJulQSSPynL+jEcBHwGJFlWu/G/eDAX8BDlNUEbcG+D55/v6SaiQdJ2mUmbUCH4TPFeczZLsSuLi9kqCkseGaAFwPHCDpmHCt6kICDda8Dk8AH0j6lqIKwFWStpPUXvn3ZuA8RZWmxwFfiXNBzCwd9r1Y0ogQ59eB68ImvwPOlrRz+Kybh22GET2k3wuf6ySiHJp27wDjwnXO5QbgJEk7Sqoluq8eN7N5xWKWdLRWVQRfFOJIF9jFuY6aiOrPfTV89z5F/h8XxVxJ9N3bFlZWzi3UnLszz4+Ofk/0vZmiqMHBhpK2MrMFwD3ATySNDOs2k7RviMm/MwOEJ2hWdwnw7ZBteTZR2esbRL+cXyCq4JqXmaWIKgKfAiwGPkdUubMlrJ9N9EvhCqIv1qtE9S3ynb+jy4kq1i0Msfwj7gczs+eB04n+mS4I58/b+RzweWBeKNo6LXyWOJ8h28+B24B7JC0NMe8WjvNfoiKtbxBlUz8D7BD2+z2wTbgOfwsJj08S1Rn6T/j8vyPKMQP4HtHf6T9ED7Y/Fr8iK32FqELv60T1nW4Arg4x/hm4OCxbSlTxcG0zewH4CdE/hneAiUB2LtT9RLlF/5O0sOMJzawR+A4wi+hvsRlwbMx4dwEel7SM6NqeaWb/if9x3UCX9Zw6keg7/Gngr1081i3Aj4A/hWfFc8AnCmzfmedHx32fAE4CfkZUOfhBVuV0Hk9UfPxCOO5fiOrsgX9nBgytXozqepqkx4kqtV1T6licc865/spzaHqYpH0lrReyck8AtqcTOSnOOeec6zyvFNzztiSqlzGcqJLdUaGM1znnnHO9xIucnHPOOVfxvMjJOeeccxWvYoqcxowZYxMmTCh1GM4NWHPmzFloZmNLHUdv8ueM60868509aL9h1vx+/tbsc55tudvMDu6x4HpBxSRoJkyYwOzZs0sdhnMDlqSu9iZbMfw54/qTznxnF76f5vG7x+VdX73+a3F7Qi+ZiknQOOecc653GEarVXZ/g56gcc45N2C8uOB17n3lIWqqajhiu4NYb1TeEWAGFE/QOOeccxXiosaf89aI2VALmHjmpX9wQM1xHLXTQaUOreQMaCVT6jC6pV8laObOmEHzrFnUTZvGxOnTSx2O68daW1uZP38+K1asKHUoPW7w4MGMGzeO6urqUofiXI9pfOlx3ho+h8Sg9q5Kotf7UtczZXk9aw0dWbrgykSm6Li65a3fJGjmzpjBZqeeytZA6p57mAueqHG9Zv78+YwYMYIJEyYgxRnsvDKYGc3NzcyfP59NNtmk1OE412PufeMeVLdmDoRl4La5jZyw25EliKp8GNBa4f3S9Zt+aJpnzaKGKIVWHead6y0rVqygrq6uXyVmACRRV1fXL3Oe3MCWJg25vq6CNmvr83jKjZmRKjBVgn6ToKmbNo0U0BqmumnTShyR6+/6W2KmXX/9XG5g23O9vbHWNf/lKWF8YquGvg+ozBiitcBUCfpNkdPE6dOZC16Hxjnn3BoOn7gfD9+b5MPRb6BBGSwDZMT2rQexweh+3V9kLAZkKiMjJq9+k6CBUGfGEzJuAFi8eDE33HADX/7ylzu13yGHHMINN9zA6NGjeycw58pUIpHgso9/j3tfauLhtx6lRjVM3fIT7LDRlqUOrSwYkKrwQpt+laApqqmJN669lgeBLY4/nvr6+lJH5FyXLF68mF//+tdrJGjS6TRVVVV597vzzjt7OzTnylYikeCgbfbkoG32LHUoZSeqFNz1BI2kjYBrgfWADDDDzH4uaW3gJmACMA84xswWdTfeXCo7OdYZTU2k99uPDa+8kqOuvJLzGhpoamoqdVTOdcm5557La6+9xo477sguu+zCfvvtx2c/+1kmTpwIwBFHHMHOO+/Mtttuy4wZM1buN2HCBBYuXFiqsHuMpKslvSvpuTzrJekXkl6V9KykSX0do3OVxBCtVpV3iqEN+IaZbQ3sDpwuaRvgXKDRzLYAGsN8rxg4OTTJJEqlqCJKie7Z2koymfRcGtdnmpqaSCaTNDQ0dPu+u/TSS3nuued45plnSCaTHHrooTz33HMrm1pfffXVrL322nz00UfssssuTJs2jbq6ftUj6h+AK4h+EebyCWCLMO0G/Ca8OudyMCDdjcq/ZrYAWBDeL5X0IrAhMBVoCJvNBJLAt7oRal4FEzSSPhXjGCvMrPzzsRsasJoaWltaaAUera7mkoaGUkflBoimpiamTJlCKpWipqaGxsbGHk1M77rrrqv1G/OLX/yCW265BYA333yTV155pV8laMzsIUkTCmwyFbjWzAx4TNJoSeuHh65zroMoh6ZgkmCMpOyRW2eY2YxcG4bv5k7A48C67d87M1sgaZ0eCnkNxXJofgvcSu7W++32AXImaCRtSVR21m5T4AKiX1V9Uqa2Un09VQ88wPxQh+YSr0Pj+lAymSSVSpFOp0mlUj2eOzhs2LDVznXffffR1NTE0KFDaWhoGIj9ymwIvJk1Pz8sWyNBI2k6MB1g/PjxfRKcc+XGTKQKFy0tNLPJxY4jaTgwCzjLzD7oy24giiVo7jKzkwttIOm6fOvM7GVgx7BdFfAWcAurytQulXRumO+VLKjV1NezcX09x8fdvqkJkkloaABP/LhuaGhooKamZmUOTUM3cwdHjBjB0qVLc65bsmQJa621FkOHDuWll17iscce69a5KlSup2jORqnhV+YMgMmTJ+dtuLpo+RIWL1/C+LU2LFjx2rlKFI3l1L37WlI1UWLmejP7a1j8TnvuqKT1gXe7F2l+BRM0Zva5YgeIs00wBXjNzN6Q1Gdlal01d8YMtjrjDAal06i2FhobPVHjuqy+vp7GxsYeq0NTV1fHnnvuyXbbbceQIUNYd911V647+OCDufLKK9l+++3Zcsst2X333Vfbd4B0nDcf2ChrfhzwdlcOtHj5En7+xAVUrf0qINrmVbPdoOM4asejeiJO58qESHevlZOA3wMvmtlPs1bdBpwAXBpeb+1OlIXErhQsaQ+iIqKV+5hZvgp5uRwL3BjexypTK1VWcFNTE3edfjoXtLUhwFpaUDLpCRrXLfX19T1azHTDDTfkXF5bW8tdd921xvJ0Os3SpUsZOXJADMJ3G3CGpD8RVQZe0tX6Mz+b/S2q155P1SADjKpBLbzYOpOHX9uAvTfboydjdq5kombb3cqh2RP4PDBX0jNh2f8RJWRulnQK8F/g6O6cpJBYCRpJfwQ2A54B0mGxkb+FQcf9a4DDgfM6E1zcrOCelkwmuT+T4VxCHnVVFdVegdhVuG233ZYvfOEL/WIUbUk3EuXyjpE0H/gu0TBumNmVRPX6DgFeBZYDJ3XlPK+/N4/q0fOpGrT6oIZVVWnun3+DJ2hcvxGjUnDh/c0eIX992yldPnAnxI1+MrBNaDHQFZ8AnjKzd8J8n5WpdUVDQwMX1dZyYEsL+ycSHH3FFUz03BlX4V566aVSh9BjzOwzRdYbcHp3zzN/8Xwy6cQaCRolIF39fncP71zZyFC0UnDZi5ugeY6o97+uNnn8DKuKm6APy9S6omN9B0/MODcwbbP+Njz0emaN5em0GNK6SY49nKtcmQrva7dYPzS3E5W6jABekPQE0NK+3swOL3YCSUOBjwOnZi3uszK1ropV38FbQTnXr40ZvjZDl9Tz4cjHGFQdlbZnMpBpq+LYbb5Q4uic6zlm6m4dmpIrlkPz/7p7AjNbDtR1WNZMH5Wp9ZowlIJSKaymhqoHHvBEjXP90Fl7nsf1c67jleV3kqhuIbFsPEdt+WU2H+s5NK7/MCDVjTo05aBYs+0HAST9yMxWa1Yt6UfAg70YW1l749pr2bClhSqgtaWF+ddey8aeoHGu30kkEnx+l+Mhfg9WzlWc9rGcKlncArOP51j2iZ4MpNI8CKSA1jAN2JSdK3tXXnkl1167ZoPEefPmsd1225UgIudcuTEgY4m8UyUoVofmS8CXgc0kPZu1agTwaG8GVu62OP54Drn6avZsbY3GhTref7258tPW1sZpp51W6jCcc2WuP+TQFCswuwG4C7iE1Yf8XmpmA7rNYn19PZckkySTSS7J1/OrVxp2veyiiy7i+uuvZ6ONNmLMmDHsvPPO3HHHHeyxxx48+uijHH744SxdupThw4dz9tlnM2fOHE4++WSGDh3KXnvtVerwnXNlot9XCjazJZKWAhPN7I0+iqliFGwJ1dQEU6ZAKgU1NT50guvxBO7s2bOZNWsWTz/9NG1tbUyaNImdd94ZgMWLF/Pgg1FB6IUXXrhyn5NOOolf/vKX7LvvvpxzzjndjsE51z/0QE/BJVe0YMzMMsC/JPkwtJ2RTEaJmXSa9IoVvJGjDoMbQNoTuN/5TvTa1NTtQz7yyCNMnTqVIUOGMGLECD75yU+uXPfpT396je2XLFnC4sWL2XfffQH4/Oc/3+0YnHP9RTSWU76pEsRto7U+8Hzoh+bD9oVx+qEZsBoaSA8aRCadptWME66+mkuOP75Hx/JxFSQrgUsqFc13814o1HH3sGHDcm4/QAamdL3g7SVvcfsLP8Nq55JuG8qEYcdw6LYFO2x2FWRA5NAE3wMOA74P/CRrcvnU13P9SSdxocQU4JF0mmQyWeqoXKk0NERFj1VV0WsPjA221157cfvtt7NixQqWLVvG3//+94Lbjx49mlGjRvHII48AcP3113c7BjcwvPPBAv7+xnEMWftRRo5czFprv83CxC+55vELSx2a6yGGaLOqvFMliJVDY2YPSloX2CUsesLMymr8pXK0xfHHc9rMmaRSKWpqamjwAS4Hrvr6qB5VD9ah2WWXXTj88MPZYYcd2HjjjZk8eTKjRo0quM8111yzslLwQQcd1O0Y3MBwx4u/pnZUC4OqVg0DUV2dpm3UPTR/eDp1w8aWMDrXE8ygNVMZRUv5KM54k5KOAS4DkkSjae4NnGNmf+nV6LJMnjzZZs+e3Ven6zFNTU0rx4QqWNzkLaIqyosvvsjWW29d6jBYtmwZw4cPZ/ny5eyzzz7MmDGDSZMmdfu4uT6fpDlmNrnbBy9jlfqc6W2/enIqa631vzWWt7RUs2nt+ey7+YDulqxsdeY7u842dfbp6w7Ou/6KnW8o++9/3Do05wO7tOfKSBoL3Af0WYKmUsUeE8pbRLkumD59Oi+88AIrVqzghBNO6JHEjHMdWdvaZDL/I9HhB3wikWHssA1LE5TrUVEdmsrOoYmboEl0KGJqJn79G1dMMom1tKBMJnrtgQqjbmC44YYbSh2CGwAmr3ciL684j0QivXJZOi2WLRvLNltvX8LIXM9RxdSVySduouQfku6WdKKkE4G/A3f2XlgDy9y6Oj7KZGgFPspkmFtXV3QfV3pximsrUX/9XK7rdt9kX0amTqGlpZpUahBtbQmWLB7HtC1/U+rQXA+J6tBU5Z0qQdxKwedImgbsSVSHZoaZ3RJnX0mjgd8B2xHlap0MvAzcBEwA5gHHmNmiTsbeb9zR3MzfEwn2zmR4OJHg0OZmJpY6KFfQ4MGDaW5upq6url81hTYzmpubGTx4cKlDcWXmyB1OoaXtc7z0zlxGD65j4618tPH+xBAZq+xnWeyxws1sFjCrC+f4OfAPMztKUg0wFPg/oNHMLpV0LtGwCt8qdJD+rKGhgYtqa3kstIa6zFtDlb1x48Yxf/583nvvvVKH0uMGDx7MuHHjSh2GK0O1g2rZYcOyrhfqusiAtoFQh0bSp4AfAesQ5dAIMDMbWWS/kcA+wIlEO6SAlKSpQEPYbCZR66kBm6Cpr6+nsbExXmsoVxaqq6vZZBP/heqc6ydMtFVI0VI+cXNofgx80sxe7OTxNwXeA66RtAMwBzgTWNfMFgCY2QJJ6+TaWdJ0YDrA+PH9e+SFWK2hnHPOuV7QH3Jo4kb/ThcSMxAlmCYBvzGznYiGTTi38C6rmNkMM5tsZpPHjvWOm5qamrjkkkto6oFxgJyrdJIOlvSypFdD0XXH9Q2Slkh6JkwXlCJO5yqBARlT3qkSxM2hmS3pJuBvQEv7QjP7a5H95gPzzezxMP8XogTNO5LWD7kz6wPe63ARTU1NTJkyZWWvw42NjZ6j4wYsSVXAr4CPEz1nnpR0m5m90GHTh83ssD4P0LkKY4i2Cu8pOG6CZiSwHDgwa5kBBRM0ZvY/SW9K2tLMXgamAC+E6QTg0vB6a2cDH2iSySSpVIp0Ok0qlSKZTHqCxg1kuwKvmtnrAJL+BEwlerY412sWL1/Io6/9kRVtC5mw9r7stNEBJDr2OFiJrPKLnOI22z6p0HpJ55nZJXlWfwW4PrRweh04iaio62ZJpwD/BY6OH/LAdFhdHR9J3J9I8JSPC+XchsCbWfPzgd1ybFcv6V/A28DZZvZ8xw0GUl091z1z/nsf76/4KrWDMwyrSrMoNYubntqSo3b8M9WDakodXrcYdCuHRtLVRINYv2tm24VlFwJfJKpLC/B/ZtZrfdj1VHIsb4LEzJ4J9WC2N7MjzGyRmTWb2RQz2yK8vt9DcfRPTU1MPOssvpfJ8EBVFY9ffrnnzriBLlehfsceAZ8CNjazHYBfEhWZr7mT19VzMaTTaf637GwG16SorW6jKmHUVrcxdtTL3PNS5Xcw2N4PTTfq0PwByDUY1M/MbMcw9WqHvD2VoKmMGkOVKpmEVAplMlRnMkxsbi51RM6V2nxgo6z5cUS5MCuZ2Qdmtiy8vxOoljSm70J0/cnzC/5JTXXLGstrq9v4sLV/1JpIWyLvVIyZPQSUNHOipxI03ld6b2poiAatrKqKXgsUN3lLKDdAPAlsIWmTUJx9LHBb9gaS1lPoxlnSrkTPO/814Loo/7859YN/gWaQziTyTsAYSbOzpukxD32GpGclXS1prV78CPF7Ci7Cc2h6U319NAJ3MhklZvIUN3lLKDdQmFmbpDOAu4Eq4Goze17SaWH9lcBRwJcktQEfAceaD1Tlumjb9fdk/r9rGFzdutryltZBDBnUHxrSqT3hks9CM+tsN9G/AS4iSg1eBPyEaPij/FFI3zezC7Lmq4Brzey4YifrqQTNn3voOC6f+vqiI3Bnt4Sa1NJCy4UXwoUX+sjdrl8KxUh3dlh2Zdb7K4Ar+jou1z9VVVWxzpAfs6j1LBIYgwa10dpWzcIPNmXaDl8pdXjd1t4PTY8e0+yd9veSfgvcEWO38e0NjSTVEqUvnopzvlhFTpJ+LGmkpGpJjZIWSvpcVtA/jHMc17saGhqoqalhz0SCezIZ9r3vPpgyBbz4yTnnum3yhIPYaf0kyz86iXcXHcmIQZdxzE63V3wLJwAM0qa8U1eEfubaHQk8F2O3k4CJks4DbgceMLML45wvbg7NgWb2TUlHElXGOxp4ALgu5v6uD7SPCdVy4YUMue8+lMlAKhUVVXkujXPOdduY4evyye3739CDVrzIqSBJNxKN0ThG0nzgu0CDpB2JMoDmAacW2H9S1uzPgauAR4EHJU0ys6K5NHETNNXh9RDgRjN7P9S1c2Wmvr4+KmZ6+OEoMVOkErFzzjkHkMl0/f+6mX0mx+Lfd+IQP+kwvwjYJiw3YP9iB4iboLld0ktEFeu+LGkssKITgbq+FKMScVNTk4/u7ZxzDohaOVkJx2wys/26e4y4PQWfK+lHwAdmlpa0nKibcVeuClQi9tZQzjnnOkp3I4emp0j6IfBjM1sc5tcCvmFm3y62b9xKwUOB04maYAFsAHS2+ZYrE8lkkkktLZwTWkMlk8lSh+Scc66EDJHJJPJOfegT7YkZADNbRFTdpai4RU7XAHOAPcL8fKKmVHGaYLkyc1hdHWdmMtQAqUyG1+rqSh2Sc/1eS+uHzF/0InXDNmL0sHVLHY5zq7Oeb7bdRVWSas2sBUDSEKA2zo5xEzSbmdmnJX0GwMw+ktcKrlgTm5uxRAJlMlQlEj6UgnO9rPH5Cxk37BoAmlekefK1ndlrq+sYUjOixJE5l6U8up28DmiUdA1RRCcDM+PsGDdBkwqpJAOQtBmw5qAWrjI0NKDa2mh8KG8F5Vyvanr1asYPv5rBg1b1MLvJqDk8/NIJHLj9X0sYmXOr604rp55iZj+WNBeYQjQKwUVmdnecfeMmaL4L/APYSNL1wJ7AiXF2lDQPWAqkgTYzmyxpbeAmYAJR2/RjQjmZ6wsxh1JwrrskbWJm/ym2rD9Lp37D4CGrd5dfOyjNZqOfYOlHCxkxxMfLdKVnBta3dWXyMrO7gLs6u1+s6M3sXuBTRImYG4HJZpbsxHn2C0OHt1ckPhdoNLMtgMYw7/pSfT2cd17BxIwPdOl6wKwcy/7S51GU0IiaJTmXZyzB4o/eybnOuVKwTP6pr0jaXdKTkpZJSklKS/ogzr6xcmhCfZlPAJua2fcljZe0q5k90cWYpxL1KAhR2VgS6H9dL1awuTNmcNfpp3N/JsNFtbXetNt1iqStgG2BUZI+lbVqJDC4NFGVxrvLt2RU7RyqEqtXUGhNVzFh1BYlisq5jlTSfmiyXAEcS9TwaDJwPLB5nB3j5i/9GqgH2nsCXAr8Kua+BtwjaU7WcOPrmtkCgPC6Tsxjub7Q1MRWZ5zBBW1t3JPJeNNu1xVbAocBo4FPZk2TgC+WLqy+t/n636clXU1bVv2EFW3VvJs6g+qqfjAGkOsfDCyjvFOfhmL2KlBlZmkzuwaI1ele3Do0u5nZJElPh5MtkhT3m7inmb0taR3g3tDjcCwhATQdYPz48XF3c92VTDIonUaE/qYTCRq84rDrBDO7FbhVUr2ZDegyy03G7sQ8/Y1/v30RYwa/xNLU2gwb9lUatjqm1KE5t7ryyKFZHtIXz0j6MbAAGBZnx7gJmlZJVaxq5TQWiFWqZmZvh9d3Jd0C7Aq8I2l9M1sQRuN8N8++M4AZAJMnTy6PBmUDQWgFZS0tUFXF0VdcwUQvbnJd87Sk04mKn1YWNZnZyaULqe9NGLMTE8Z4iyZX5vqwrkwBnycqPToD+BqwETAtzo5xi5x+AdwCrCPpYuAR4IfFdpI0TNKI9vfAgUTDh98GnBA2OwG4NWYcri+EVlD6wQ+ofvBBJk6fXnwf53L7I7AecBDwIDCOqMjaOVdOjCiHJt/UB0LGycVmtsLMPjCz75nZ10MRVFFFc2gkJYD/AN9kVbvwI8zsxRjHXxe4JfTBNwi4wcz+IelJ4GZJpwD/BY6OE6zrQwXGgsrWPsjlYXV1UQd93gzcrW5zMzta0lQzmynpBiBWnxLOub7Vl62Zcp4/GityrKQaM0t1dv+iCRozy0j6iZnVA7Hrv4R9Xwd2yLG8mShx5CpY+yCXk1paODOTiXofrq2N+rjxRI2LtHfAsljSdsD/iPqfcs6VGZVBx3pEfdM9Kuk24MP2hWb202I7xi1yukfSNB/uwGVLJpOkUin2DuNCKZOBVCrqsM+5yIwwWu53iIqaXwB+XNqQnHNrMEGmwNR33iYaJzIBjAjT8Dg7xq0U/HWiWsZtklYQFTuZmY3sfKyuv2hoaKCmpoaHW1pIhXGhfCgFl83MfhfePghsWspYnHNFlEfTmxfM7M/ZCyTFqpYSK0FjZj6CmltDfX09jY2NJJNJXvM6NC6LpK8XWh8n+zjGOQ4Gfg5UAb8zs0s7rFdYfwiwHDjRzJ7q7nld16RSK7j9vm+x/YR/MqS6jafnbc4OW1/M+HEfK3Vorl15tHI6j6hTvWLL1hC3p+BJORYvAd4ws7Y4x3D9U319vfcg7HLp1R9BoTXEr4CPA/OBJyXdZmYvZG32CWCLMO0G/Ca8uhJINk3joO1fZ+igNABjtpnLe8s/zcL372LM2uuVODqHlbYOjaRPEP342FDSL7JWjQRipTPiFjn9mqiHz7lhfiLwL6BO0mlmdk/M47gBaO6MGTTPmkXdtGneBHyAMLPvxdlO0nlmdkkXTrEr8GpoeICkPxENqZKdoJkKXGtmBjwmaXR7/1ddOJ/rhmeeS7LXZq8zOCRmAGqqMqw1uIW7my5l2qGXly44t0ppi5zeBmYDhwNzspYvJeqPpqi4CZp5wClm9jyApG2Ac4CLgL8CnqBxOc2dMYPNTj2VrYHUPfcwFzxR47IdDXQlQbMh8GbW/HzWzH3Jtc2GRD2Puj40b/4DbDoyAaRXWz60uo21RnWq8azrRSphgsbM/gX8S9INZtaabztJs8wsZ0d7cVs5bdWemAknfgHYqf3XkXP5NM+aRQ1Ryrk6zDuXpat53Ln26/g4jrMNkqZLmi1p9nvvvdfFcFwhI4ZNQDn+W6bSCZZ+WFeCiNwajLJo5VQoMRPkbVwQN0HzsqTfSNo3TL8G/i2pllX9TDi3hrpp00gR3SStYd65LF39TTifqEv0duOIsqw7uw1mNsPMJpvZ5LFjx3YxHFfIfnscz7vLhq02QCdAWybBZhufXqKoXEfK5J/KSN5nRtwEzYnAq8BZRGVZr4dlrcQcBdMNTBOnT+e1q67i0QMP5LWrrvLiJtdRV3/6PQlsIWmTMJDdsUT93GS7DThekd2BJV5/pjQSVVVY7e95dsG6tKQTfNRWxVtLh/HP17/OdlvtUerwXDsrMFWAuM22Pwq5MneY2csdVi/r+bBcfzJx+nQolJBpaoo64/Mm3/2GpB+Z2bckHd2xT4kOijbFzMXM2iSdQTSMQhVwtZk9L+m0sP5K4E6iVhOvEjXbPqkr53I9Y/NNtwce4o03X2L58iVsuflkNtqiqtRhuUDWvZwYSVcDhwHvmtl2YdnawE1EvYPPA44xs0XdDTXfilg5NJIOB54B/hHmdwzdEjvXPU1NMGUKfOc7pPfbj2u/9CWamppKHZXrvkMkVRP1H5GXmRUd5LbAvnea2cfMbDMzuzgsuzIkZrDI6WH9RDOb3dVzuZ6z8UZbsfWWu5Go8sRM2eleHZo/AAd3WHYu0GhmWwCNYb67vpVvRdwip+8SNZNcDGBmz+DjsbiekExGwyWk02RaWnj5qquYMmWKJ2oq3z+AhcD2kj7ImpZK+qDUwTnn1tSdOjRm9hDwfofFU4GZ4f1M4IiiMUhzJT3bYXpY0s8k1RXqJiZugqbNzJbE3Na5+BoaoKaGtEQrcL8ZqVSKpI8HVem+bWajgL+b2cisaYQPmeJcGbJQ7JRnAsa0twYMU5wKkeu211sLr+vE2Ocu4O/AcWG6HXiIaGDbPxTaMW4/NM9J+ixQJWkL4KvAP2Pu61x+9fXQ2Mj8a6/lhKuv5sl0mpqaGhp8PKhK10TUGafnxjhXKQrnxCw0s8l9EMWeZrZn1vxcSY+a2Z6SPldox7gJmq8A5wMtwI1EFfEuihtd6KZ8NvCWmR3WSxWFXKWqr2fj+nouOf54kskkDQ0NPpxC5auRdAKwh6RPdVxpZn8tQUzOuQJ6oWO9d9p755a0PvBujH2GS9rNzB4HkLQrq0bbLjgEQtxWTsuJEjTnx9k+hzOBF4nGZIBVFYUulXRumM9b0ccNDLHHhfJWUZXgNKLs4tHAJzusM6Iexp0bcFpSrVx5+/08/t832GjkaL429eOsN2atUocVxnLq8aPeBpwAXBpeb42xzxeAqyUNJ2rR9AHwBUnDKNKreMEEjaTbKdAC3cwOLxaZpHHAocDFQPsIvFOBhvB+JpDEEzQuhrkzZrDVGWcwKJ1GtbXQ2OiJmjJkZo8Aj0iabWa/L3U8zpWD997/gP1+92s+HCYyIxIkWhbx5z/N4Le7HcqUXbYrdXjd6m9G0o1E/9fHSJpP1JjoUuBmSacA/yUa6qRwCGZPAhMljQJkZouzVt9caN9iOTT/L7x+ClgPuC7Mf4aoqCiOy4Fvsvrou6tVFJKUs6JQqHQ0HWD8+PExT+f6q6amJu46/XQuaGtDgLW0oGTSEzRlSNL+ZnY/sMiLnJyLnDnzJpaNEFYdtcfJ1CQgY3w1eQfPlzhBI7qXQ2Nmn8mzakqn4ohGIJhGVCVlkKT243+/2L4FEzRm9mA4wUVmtk/WqtslPRQjsPZOduZIaii2fY7zzwBmAEyePLlC+ip0vSWZTHJ/JsO5hB8SVVVUe+XhcrUPcD9RcZMRPS+zXz1B4wac2daMVXfofychPhwuXpr3NltN2KA0gUFvFTl1xa3AEqIRt1s6s2PcSsFjJW3aPhilpE2AOIOe7AkcLukQYDAwUtJ1dK2ikBvgGhoauKi2lgNbWtg/keDoK65goufOlKulkr4OPMeqhAxUTCfqzvW8vJVuJaoHlUFHg+WRoBlnZh076IslboLma0BSUvvo2hMIRUGFmNl5hJ5CQw7N2Wb2OUmX0fmKQm6Aq6+vp7GxcWVLqLyJGa80XA7aWyVsCexC9B0XUY5N0dxd5/qjPWvW5f7Wd1cWOQGQMUYuMTYbt27pAgt6oZVTV/xT0kQzm9vZHeO2cvpH6H9mq7DoJTPrVFZQB52uKOQcxGgJ1T6UQioFNTVeabhEzOx7AJLuASaZ2dIwfyFdHL/JuUr385OPZb9fXkHziAyZKpFIG1Vtxm8PParUoUV5p+WRQ7MXcKKk/xAVOYloJJPti+1YrJXTJDN7iuhoLcC/Cm1TiJkliVozYWbNdLKikHOxJJNRZeFMxisNl4fxQCprPoUPm+IGqBHDhvDEN7/Bzfc/zqOvvMb4urU47fD9GTF0SKlDA8qmDs0nurpjsRyaa0JRUaGRqX4P7NTVAJzrSXPr6tgsk6EaaM1keK2ujomlDmpg+yPwhKRbiH4DHsmqsV2cG3ASiQTHHlDPsQeU3w+tckjQmNkbknYA9g6LHjazNTJTcimWoBlFVNO4UILmvTgncq4v3NHczN8TCfbOZHg4keDQ5mZP0JSQmV0s6S5WPZxOMrOnSxmTcy4Hoyyq7Es6E/giq1pCXidphpn9sti+xZptT+h+eM71nfaWUI+lUtTU1HCZN+suuVAkXbRY2jlXOqJsKgWfAuxmZh8CSPoR0dhw3UvQOFdpOraEylmB2FtBOefcGsqhyIkobZXOmk9TuJRoJU/QuH6nYEsobwXlnHNrKp+O9a4BHg/17gCOIKqrW1Si+CbO9SPJZJSYSaexlhaSF15IU1NTqaNyrqhMJkPyXy/ywDMv0JZOF9/Buc6yAlNfhWD2U+Ak4H1gEVG9u8vj7Bsrh0bRYArHAZua2fcljQfWM7MnuhaycyXS0AA1NVhLCx9lMnz7vvt46uGHaWxsjDfSt3MlcOcT/+IbT99B6+BoftCTcMm2B3PkHjuXNjDXr5Qyh0bS2lmz88gaL1LS2mb2frFjxM2h+TVQTzQoJcBS4Fcx93WufNTXQ2MjDx5wAAcmEjyayZBKpUgmk6WOzLmc3lv8AWc+fwctI0WmJppSw8U5//4Hby0s+ox3Lp5Q5JRv6gNzgNnhtf397Kz3RcWtQ7ObmU2S9DSAmS2SVNP5eJ0rA/X11F54IU89/DBVoTVUQ6HWUF6J2JXQr+95AMvx09MEv7i7kR8d5x2tu+7r7mjb3WVmm6yMJcqt2YJoDMjY4iZoWiVVEUrSJI2lXDpJdq4LYrWGAq9EXIbCw+4moh6H5wHHmNmiHNvNI8pNTgNtZja576LsOe8s+4BMjp+PVg3vfbis7wNy/VcZNNuW9AXgTGAc8AywO/BPYowuELfI6RfALcA6ki4GHgF+2JVgnSsX9fX1nHfeeYXrzmRVIiaViuZdqZ0LNJrZFkBjmM9nPzPbsVITMwBTPrYViRx1gBOtsO+mW/R9QK5/MlDG8k596EyiAW3fMLP9iEYiWBhnx1gJGjO7HvgmcAmwADjCzHyAOdfvza2rozWRwBKJKIfGO+orB1NZNXzCTKJmnf3WkXvszHpLalFq1T8VpYwxH1Tz+f32KGFkrr8pcR2adivMbAWApFozewnYMs6OBRM0ktZun4B3gRuBG4B3OtRIdq7faWpqYrezzmK/dJrvJhLMvfxyL24qD+ua2QKA8LpOnu0MuEfSHEnT8x1M0nRJsyXNfu+98hvJJZFI0HjamXy6ektGNScY1SymVW1B4xfPJJHwnjdczymTBM18SaOBvwH3SroVeDvOjsXq0MwheiiIaNTcReH9aOC/wCZ59wQkDQYeAmrDuf5iZt+NWwbuXCklk0lSqRSPZjI8JjHEx4XqM5LuA9bLser8ThxmTzN7W9I6RA/Gl8zsoY4bmdkMYAbA5MmTy6AWwZpqa6r54WeP8nJ+13usPIY+MLMjw9sLJT1ANKbkP+LsW2wsp00AJF0J3GZmd4b5TwAHxDh+C7C/mS2TVA08Egaq+xRRGfilks4lKgP/VpyAnesrDQ0N1NTUkIrTEsr1KDPL+3yR9I6k9c1sgaT1iXKPcx3j7fD6buh1dFeiH1jOuQ5K3copFzN7sDPbx82v3KU9MRNOchewb4xgzMzaq+FXh8kYYGXgrjK1t4S66KKLvOO98nIbcEJ4fwJwa8cNJA2TNKL9PXAg8FyfRehcBSqTSsFdFrfZ9kJJ3wauI0qQfA5ojrNjaO49B9gc+JWZPS5ptTLwkCWca9/pwHSA8ePHxwzVuZ5TcFwoVyqXAjdLOoWo6PtoAEkbAL8zs0OAdYFbok7OGQTcYGaxsq2dG5AMVOEjasRN0HwG+C5R022Ism0/k3/zVcwsDewYKvncImm7uMFVQtm2c01NTcX7s3E9xsyaydEnRShiOiS8fx3YoY9Dc66yVfh/2VgJmjCGwpndOZGZLZaUBA4maiVVtAzcuXLX1NTElClTVtaz8aIp51xFCv3QVLK4g1M+QI60m5ntX2S/sUBrSMwMIapI/CNWlYFfSp4ycOcqQTKZZFJLC3tnMjzc0kIymfQEjXO95IZH/skV/06ydGiK4curOW2zfThh371LHVa/0d1KwaXunTtukdPZWe8HA9OAthj7rQ/MDPVoEsDNZnaHpCZylIE7V2kOq6vjzEyGGiCVyfBaXV2pQ3KuX7om+SA/WfgAFr5iHwxp5bL3GvnovhSnHVC0V3xXhKzHKv/uZ2axevbtaXGLnOZ0WPSopKLNqczsWaJuizsuz1kG7lylmdjcjCUSKJOhKpFgYnOsuvLOuU761RsPYaNXX2bVMOPtRznN/530iHLoh6Y74hY5ZfcKnAB2JnenV84NLA0NqLYWUinkQyM412s+GtlG1FvK6lpGpslkMt5rcg/ogX5o2nvnNuCq0LCnz8QtcsruMbgN+A9wSm8F5VzFqK+PRuBOJqPETIH6M94ayrmuG7Q8QdvwNbMQqj6SJ2Z6ggHpglk0YyTNzpqfkSPBEqt37t4SN0GzdftgUe0k1fZCPM5Vnvr6omM8eWso57rnk4O34W+p57GaVcuUggMTHytdUP1MkTo0C4tV8i1179xxk7X/zLGsqScDca4/ax8Xapd0mq+tWMEr115b6pCcqyjfP3IaB6S3ILEC1AqJFtg3tQk/PvrYUofWb8jyT0X3LYPeuQvm0EhaD9gQGCJpJ1YVYI4EhvZybM71Gw0NDexVVcWd6TQ1Zuiaa+D44330budiSiQSXH7scbSkWnn7/fdZb621GFJbU3xHF4u63w9NyXvnLlbkdBBwIjAO+GnW8qXA//VSTM71O/X19cw8+WQGX3UVCTNoa4vq3XiCxrlOqa2pZpP11i11GP2SCtehKagceucuNtr2TKJ+ZKaZ2aw+ism5fmnj44+HmTMhlQJvEeWcKydm0J97Cpb0OTO7Dpgg6esd15vZT3Ps5pzLJUaLKG8J5Zwrlf7eD82w8Do8x7oK/+jOlUCBFlHeEso5VzLWvSKnclCsyOmq8PY+M3s0e52kPXstKucGIB8XyjlXUv25yCnLL4FJMZY557rIx4VyzsVx17+eYtbLTRjGpz5Wz6E77twjx1Wm+10Fl1KxOjT1wB7A2A51aEYCVb0ZmHMDjY8L5Zwr5vS/XMW/hr0Go6P5Sxb8l9v+/QRXHfOl7h3YgMpOzxTtWK+GqP7MIGBE1vQBcFTvhubcANM+LlRVVfTqraCcc1maXnk5SsxUE/33TgDV8PyI//DQS89369jCUCaTd6oExerQPAg8KOkPZvZGZw8uaSPgWqKBLDNEYz/8PAx2eRMwAZgHHGNmizp7fOf6lZjjQnlLKOcGpr89/zgMzrEiAbe+8AT7bLVt1w9efCynshe3Ds1ySZcB25J1Oc1s/yL7tQHfMLOnQpfIcyTdS9RZX6OZXSrpXOBc4Fudjt65/qbIuFBzZ8zgrtNP5/5Mhotqa70llHMDSE1iUO72xQY1VXH/nedXKTkx+cQdy+l64CVgE+B7RLkqTxbbycwWmNlT4f1S4EWioRSmAjPDZjOBIzoTtHMDUlMTW51xBhe0tXFPJsOk0BLKOTcwfG7yvrlXGBw3Kc+62CzqXC/fVAHiJmjqzOz3QKuZPWhmJwO7d+ZEkiYAOwGPA+ua2QKIEj3AOnn2mS5ptqTZ7733XmdO51z/k0wyKJ1mEFER+v6JBA0DsJ6NpKMlPS8pIynv6L+SDpb0sqRXQ06wcxVtyw025FOD9oJWIAWkBK1wuOrZbqPx3Tt4e5FTvqkCxM2jag2vCyQdCrxNNL5TLJKGA7OAs8zsgzB4VVFmNgOYATB58uTKuKLO9ZZQadhaWqCqiqOvuIKJA7O46TngU8BV+TaQVAX8Cvg4MB94UtJtZvZC34ToXO/4+scP5+jmPfjT7EcA49OT9mL82LE9cuxKL3KKm6D5gaRRwDeI+p8ZCZwVZ0dJ1USJmevN7K9h8TuS1jezBZLWB97tXNjODUCh0rCSSaobGgZqYgYzexGgyA+jXYFXw4B5SPoTUVG3J2hcxduobgznHHREzx7UDNIDIEFjZneEt0uA/QAknVVsP0VPnN8DL3YY9+k24ATg0vB6a/yQnRvAilQabuctodgQeDNrfj6wW64NJU0HpgOMH9/NbHvnKlmF1JXJpzvVor8OXF5kmz2BzwNzJT0Tlv0fUULmZkmnAP8Fju5GHM65LNljQu1VVcXMk0+ORvquoISNpPuIunvo6Hwzi/MDKFf2Tc6ntRdtO0eoQzMAcmjyKFoRxsweKbDdlG6c2zmXRzKZJJVKsUs6zZ3pNIOvugpmzoz6uKmQRI2ZHdDNQ8wHNsqaH0dU9885l5NBJl3qILolbiunXPyXjHNlqKGhgZqaGvaXqAESZpBKRR32DRxPAltI2kRSDXAsUVG3cy6X9hyafFMFKJigkbRU0gc5pqXABn0Uo3OuE+rr62lsbGTLU09dOZQCNTX9ZigFSUdKmg/UA3+XdHdYvoGkOwHMrA04A7ibqP+rm82se33D9wPPvvk6f33qYd5s9nYYLocK74em2NAHI/oqEOdcz6mvr48qAx9/fNGhFCqNmd0C3JJj+dvAIVnzdwJ39mFoZat56Qd89YEf8dGoxZAW1y0zxi3elMsP+xqJRHcy6l2/YQbpgVvk5Jwrd/X1cN55eRMzc2fMIHnQQcydMaOPA3N96Zz7L+ej0YtQtaHBGTTImD/ydS67/0+lDs2Vk0wm/1QBuj/4g3OuIs2dMYPNTj2VrYHUPfcwF5g4fXqpw3I9bNGHS1m81v9Qh6e9aownP3oS+GxJ4nJlxgzzHBrnXCVqnjWLGlg5lELzrFkljsj1hiXLPyRfY9NMTVvfBuPKW8byTxXAEzTODVB106aRIhrXpDXMu/5nfN066KM1M+MtA6OX5hxGzw1E7XVo8k0VwIucnBugJk6fzlyinJm6adO8uKmfSiQSHL32VG5ePguqDCXA2oC2BF/b+bhSh+fKSKUXOXmCxrkBbOL06eAJmX7v2F32Z9y/xzDzxdv5ILGEDdiA03c9ms3X3bDUoblyMVDGcnLOOVfZ9vrY9uz1se1LHYYrZ9b1BI2kg4GfA1XA78zs0p4KKy5P0DjnnHMDnHWjlZOkKuBXwMeJhh15UtJtZtano9t7gsY555xz3alDsyvwqpm9DiDpT8BUoE8TNLIK6dJY0nvAGzE2HQMs7OVwuqqcY4Pyjq+cY4Pyjq+nYtvYzMb2wHHKVieeM8WU2/3g8RTWX+OJ/Z2V9I9w3nwGAyuy5meEkeqRdBRwsJl9Icx/HtjNzM7oWthdUzE5NJ34o8w2s8m9HU9XlHNsUN7xlXNsUN7xlXNs5aanEmzlds09nsI8HjCzg7uxe66Ojvo8t8T7oXHOOedcd8wHNsqaHwe83ddBeILGOeecc93xJLCFpE0k1QDHArf1dRAVU+TUCeU8yl45xwblHV85xwblHV85x9Zflds193gK83i6wczaJJ0B3E3UbPtqM3u+r+OomErBzjnnnHP5eJGTc8455yqeJ2icc845V/EqKkEjaSNJD0h6UdLzks4My9eWdK+kV8LrWln7nCfpVUkvSzqoRPFdJuklSc9KukXS6LB8gqSPJD0TpitLENuFkt7KiuGQrH365NoViO2mrLjmSXomLO+z6xbON1jSE5L+FeL7Xlhe8vuuQGwlv+cGEklHh+ufkZS3ua2kg8M98aqkc3sxnrz3Zoft5kmaG+6F2b0QR8HPq8gvwvpnJU3q6Rg6GU+DpCVZ348LejGWqyW9K+m5POv79Nr0C2ZWMROwPjApvB8B/BvYBvgxcG5Yfi7wo/B+G+BfQC2wCfAaUFWC+A4EBoXlP8qKbwLwXImv3YXA2Tm277Nrly+2Dtv8BLigr69bOJ+A4eF9NfA4sHs53HcFYiv5PTeQJmBrYEsgCUzOs01VuBc2BWrCPbJNL8WT897Msd08YEwvxVD08wKHAHeF+3h34PFe/BvFiacBuKOP7pl9gEn5vo99eW36y1RROTRmtsDMngrvlwIvAhsSdbE8M2w2EzgivJ8K/MnMWszsP8CrRF0092l8ZnaPmbWFzR4jaqPfpwpcu3z67NoVi02SgGOAG3vj/DHiMzNbFmarw2SUwX2XL7ZyuOcGEjN70cxeLrLZyu7hzSwFtHcP3xvy3Zt9Kc7nnQpcG+7jx4DRktYvYTx9xsweAt4vsElfXpt+oaISNNkkTQB2IvpFuq6ZLYDonyOwTthsQ+DNrN3mU/ifeG/Fl+1kolR3u00kPS3pQUl7lyi2M0KW5tVZWdMluXZ5rtvewDtm9krWsj69bpKqQpHXu8C9ZlY2912e2LKV/J5zQN/eF/nuzY4MuEfSHEnTeziGOJ+3L69J3HPVhyLcuyRt20uxxFGy/1+VqiL7oZE0HJgFnGVmH0Q/4HNvmmNZr7dT7xhf1vLzgTbg+rBoATDezJol7Qz8TdK22fv0dmySfgNcRHRdLiIq2jmZEly7fNcN+Ayr5870+XUzszSwY6iLcouk7Qps3qfXLldsZvYclMc9119Iug9YL8eq883s1jiHyLGsy/dFoXg6cZg9zextSesA90p6KeQc9IQ4n7cvvytxzvUU0fhHyxTVJ/wbsEUvxVNMWQwnUEkqLkEjqZron971ZvbXsPgdSeub2YKQJfduWN7n3THniQ9JJwCHAVPMogJSM2sBWsL7OZJeAz4G9HjlvHyxmdk7Wet/C9wRZvv02hW4boOATwE7ty/r6+uWzcwWS0oCB1NG912O2J4rh3uuPzGzA7p5iB69LwrFIynfvdnxGG+H13cl3UJULNNTCZo4n7cvvytFz5WdsDezOyX9WtIYMyvFwJVlMZxAJamoIqdQl+L3wItm9tOsVbcBJ4T3JwC3Zi0/VlKtpE2IUtpP9HV8kg4GvgUcbmbLs5aPlVQV3m8a4nu9j2PLLpM9Emivcd9n167A3xXgAOAlM5uftX2fXbes840O74e0x0QZ3Hf5YiuHe86toS+7h893b64kaZikEe3viSqS52xx00VxPu9twPGK7A4saS8q6wVF45G0XngeIWlXov+Rzb0UTzF9eW36h67UJC7VBOxFlOX2LPBMmA4B6oBG4JXwunbWPucT1Wx/GfhEieJ7lagstH3ZlWH7acDzRLXtnwI+WYLY/gjMDctvA9bv62uXL7aw7g/AaR2277PrFs63PfB0iO85VrW2Kvl9VyC2kt9zA2ki+jEwnyj36x3g7rB8A+DOrO0OIWrF9xpRUVVvxZPz3syOh6i1z7/C9HxvxJPr8wKntX+niYpVfhXWzyVPC7E+jOeMrO/HY8AevRjLjURFwK3h3jmllNemP0w+9IFzzjnnKl5FFTk555xzzuXiCRrnnHPOVTxP0DjnnHOu4nmCxjnnnHMVzxM0zjnnnKt4nqBxzjnnXMXzBE0fk7Ss+FbdOv6dkkaH6ctd2L9B0h3Ft1xt+yWS7syz/g+SjupsHJUoXIs9sua/Jum/kq4oZVzO9RZJF0o6W9L3JR0Qlu0t6XlJz0gaIumyMH9ZoWP0beSuP6q4oQ9cYWZ2CKwc5PHLwK/74LQPm9lhvXkCSYNs1ejR5aoBWAb8E8DMfiZpETC5lEE519vM7IKs2eOA/2dm1wBIOhUYa9GwG871Gs+hKQOSdpT0mKIRr29RGPFaUlLSjyQ9IenfCiMjSxoq6eaw/U2SHpc0OaybJ2kMcCmwWfiVdFnHnBdJV0g6Mbw/WNJLkh4hGjepfZthikbgflLR6MxTY3wWhWO/IOnvZI3yK2lnRSM8z5F0d/uwC5J2CZ+lKcTaPrDiiZL+LOl2ohGBc8ajaLTpy8LyZ8MDFEnrS3ooXIPnVGBkaUkHhvM/Fc45PCy/IBz3OUkzsrpF/2r4jM9K+lNIQJ4GfC2cz0exdv2SpPMlvaxocMwtw7I/SDpK0heAY4ALJF0v6TZgGPC4pE/HOHannoXOZfMETXm4FviWmW1P1MX1d7PWDTKzXYGzspZ/GVgUtr+IrIEbs5wLvGZmO5rZOflOLGkw8Fvgk8DerD567/nA/Wa2C7AfcJmiMV8KOZLoITcR+CKwRzhPNfBL4Cgz2xm4Grg47HMNUXff9UC6w/HqgRPMbP8C8ZxCNM7JLsAuwBcVjaH0WaIu6HcEdiAaAiDXNRgDfBs4wMwmEQ3U+PWw+goz28XMtgOGEA32CNH13Sn8DU4zs3nAlcDPwjV/uMh1cq7iKBqh/VhgJ6IfP7tkrzez3xENoXKOmR1nZocDH4XvxE0xTtHZZ6FzK3mRU4lJGgWMNrMHw6KZwJ+zNmkfeXoOMCG83wv4OYCZPSfp2W6EsBXwHzN7JcRzHTA9rDsQODyrfHswMB54scDx9gFuNLM08Lak+8PyLYHtgHtDJkcVsEDRwIojzOyfYbsbWJVoALjXzN4vEs+BwPZaVVdnFNGgi08CV4fE1N/M7Jk8Me8ObAM8GmKrAZrCuv0kfRMYCqxNNM7L7URjJ10v6W/A3wpcD+f6k72BWywMeBpyYHpEF5+Fzq3kCZry117unGbV30tdOE4bq+fIDc56n29ALwHTzOzlTp4r1/EEPB9yYVYtDFnKBXxYLJ5QDPQVM7t7jZNK+wCHAn+UdJmZXZsntnvN7DMd9h1MVAdpspm9KelCVl23Q4kSb4cD35G0bZHP4Vx/UaoBAHM9C51byYucSszMlgCLssqEPw88WGAXgEeIyqmRtA1R8U5HS4ERWfNvANtIqg2/hKaE5S8Bm0jaLMxn/1O/G/hKVr2RnWJ8pIeAY0O9lvWJioYgGnV6rKT6cKxqSdua2SJgqaTdw3bHFjh2vnjuBr4UcmKQ9LFQ32Zj4F0z+y3we2BSnuM+BuwpafOw/1BJH2NV4mVhqFNzVFifADYysweAbwKjgeGsec2d628eAo5U1HppBFFRdY/o4rPQuZU8ldv3hkqanzX/U+AE4EpJQ4HXgZOKHOPXwMxQ1PQ0UfHHkuwNzKxZ0qOhgu1dZnaOpJvDtq+E/TCzFZKmA3+XtJAosbRdOMxFwOXAsyERMY/Vi4NyuQXYn6j8+9+EB5KZpUKR0C9CgmpQOPbzRHVgfivpQyDZ8bNkyRfP74iyoJ8Ky98DjiBqdXSOpFai1kfH5zqomb2nqIL0jZJqw+Jvm9m/Jf02fJZ5REVYEBWXXRc+h4jqzSxWVHn5L4oqK3/F69G4/sbMnpJ0E1F9tDeAnr7HO/ssdG4lmZUq99B1laQqoDokRjYDGoGPmVmqBLE0AGd3p9m2pOFmtiy8PxdY38zO7JkISysklCab2RmljsU55/ozz6GpTEOBB0IRi4AvlSIxE6SA7STd2d4HThccKuk8ovvxDeDEngqulCR9jagp96xSx+Kcc/2d59C4AUXS40Bth8WfN7O5pYjHuYFC0vnA0R0W/9nMLs61vXOd5Qka55xzzlU8b+XknHPOuYrnCRrnnHPOVTxP0DjnnHOu4nmCxjnnnHMV7/8D7uqsMLTH2ecAAAAASUVORK5CYII=\n", + "text/plain": [ + "
" + ] + }, + "metadata": { + "needs_background": "light" + }, + "output_type": "display_data" + } + ], + "source": [ + "fig, ax = plt.subplots(1, 2, figsize=(8, 3))\n", + "\n", + "ds_trajectory.plot.scatter(\n", + " x=\"longitude\", y=\"latitude\", c=\"k\",\n", + " s=9, ax=ax[0], label=\"traj\"\n", + ")\n", + "selected.plot.scatter(\n", + " x=\"lon2d\", y=\"lat2d\", c=\"r\",\n", + " s=9, ax=ax[0], label=\"grid\"\n", + ")\n", + "ax[0].legend()\n", + "ax[0].set_title(\"target and selected locations\")\n", + "\n", + "selected.plot.scatter(\n", + " x=\"diff_lon\", y=\"diff_lat\", hue=\"along_track\",\n", + " ax=ax[1]\n", + ")\n", + "ax[1].set_title(\"differences\")\n", + "\n", + "fig.tight_layout();" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.8.8" + } + }, + "nbformat": 4, + "nbformat_minor": 4 +} From 629432c4a9e230d6e8716f96d89ab79d6d46045b Mon Sep 17 00:00:00 2001 From: Willi Rath Date: Sun, 1 Aug 2021 13:04:26 +0200 Subject: [PATCH 11/11] Link correct query method --- doc/api.rst | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/doc/api.rst b/doc/api.rst index 25a1799..9e5306e 100644 --- a/doc/api.rst +++ b/doc/api.rst @@ -36,7 +36,7 @@ properties listed below. Proper use of this accessor should be like: Dataset.xoak.set_index Dataset.xoak.sel - Dataset.xoak.get_indexers + Dataset.xoak.query DataArray.xoak -------------- @@ -59,7 +59,7 @@ The accessor above is also registered for :py:class:`xarray.DataArray`. DataArray.xoak.set_index DataArray.xoak.sel - DataArray.xoak.get_indexers + DataArray.xoak.query Indexes -------