From 38beef72f8ed3522cd7357f09b4eda44e1bc0f97 Mon Sep 17 00:00:00 2001 From: Sonja Stockhaus Date: Mon, 18 Mar 2024 22:30:39 -0400 Subject: [PATCH 01/23] prototype version with changed behavior regarding element order --- src/spatialdata_plot/pl/basic.py | 4 ++ src/spatialdata_plot/pl/render.py | 59 +++++++++++++++++++++++- src/spatialdata_plot/pl/render_params.py | 4 ++ 3 files changed, 65 insertions(+), 2 deletions(-) diff --git a/src/spatialdata_plot/pl/basic.py b/src/spatialdata_plot/pl/basic.py index 0f18c5d7..a13e1366 100644 --- a/src/spatialdata_plot/pl/basic.py +++ b/src/spatialdata_plot/pl/basic.py @@ -342,6 +342,7 @@ def render_shapes( outline_alpha=outline_alpha, fill_alpha=fill_alpha, transfunc=kwargs.get("transfunc", None), + zorder=n_steps, ) return sdata @@ -499,6 +500,7 @@ def render_points( alpha=alpha, transfunc=kwargs.get("transfunc", None), size=size, + zorder=n_steps, ) return sdata @@ -654,6 +656,7 @@ def render_images( alpha=alpha, quantiles_for_norm=quantiles_for_norm, scale=scale, + zorder=n_steps, ) return sdata @@ -815,6 +818,7 @@ def render_labels( fill_alpha=fill_alpha, transfunc=kwargs.get("transfunc", None), scale=scale, + zorder=n_steps, ) return sdata diff --git a/src/spatialdata_plot/pl/render.py b/src/spatialdata_plot/pl/render.py index f3c57af9..07897cd3 100644 --- a/src/spatialdata_plot/pl/render.py +++ b/src/spatialdata_plot/pl/render.py @@ -5,6 +5,7 @@ from typing import Union import dask +import datashader as ds import geopandas as gpd import matplotlib import matplotlib.transforms as mtransforms @@ -13,10 +14,12 @@ import scanpy as sc import spatialdata as sd from anndata import AnnData +from dask.array import from_array from matplotlib.colors import ListedColormap, Normalize from multiscale_spatial_image.multiscale_spatial_image import MultiscaleSpatialImage from pandas.api.types import is_categorical_dtype from scanpy._settings import settings as sc_settings +from spatial_image import to_spatial_image from spatialdata._core.data_extent import get_extent from spatialdata.models import ( PointsModel, @@ -24,6 +27,7 @@ from spatialdata.transformations import ( get_transformation, ) +from spatialdata.transformations.transformations import Identity from spatialdata_plot._logging import logger from spatialdata_plot.pl.render_params import ( @@ -132,7 +136,8 @@ def _render_shapes( cmap=render_params.cmap_params.cmap, norm=norm, fill_alpha=render_params.fill_alpha, - outline_alpha=render_params.outline_alpha + outline_alpha=render_params.outline_alpha, + zorder=render_params.zorder, # **kwargs, ) @@ -236,6 +241,47 @@ def _render_points( points = dask.dataframe.from_pandas(points, npartitions=1) sdata_filt.points[e] = PointsModel.parse(points, coordinates={"x": "x", "y": "y"}) + # TODO: maybe move this, add heuristic + if len(points) > 100: + extent = get_extent(sdata_filt.points[e], coordinate_system=coordinate_system) + x_ext = extent["x"][1] + y_ext = extent["y"][1] + previous_xlim = fig_params.ax.get_xlim() + previous_ylim = fig_params.ax.get_ylim() + x_range = [0, x_ext] + y_range = [0, y_ext] + # if sth else was plotted before: might need larger extent/size + # else: circles representing a point are "cut off at image borders" + if previous_ylim[0] > previous_ylim[1]: + if previous_xlim[0] < x_range[0]: + x_range[0] = previous_xlim[0] + if previous_xlim[1] > x_range[1]: + x_range[1] = previous_xlim[1] + # different order for x and y! + if previous_ylim[0] > y_range[1]: + y_range[1] = previous_ylim[0] + if previous_ylim[1] < y_range[0]: + y_range[0] = previous_ylim[1] + # round because we need integers + plot_width = int(np.round(x_range[1] - x_range[0])) + plot_height = int(np.round(y_range[1] - y_range[0])) + + # use datashader for the visualization of points + # TODO: incorporate color, size etc + cvs = ds.Canvas(plot_width=plot_width, plot_height=plot_height, x_range=x_range, y_range=y_range) + agg = cvs.points(sdata.points[e], "x", "y") + ds_result = ds.tf.shade(ds.tf.spread(agg, px=5), rescale_discrete_levels=True) + + # TODO: we might not even need this conversion (directly use numpy array) for direct rendering + dsimage = to_spatial_image( + from_array(np.transpose(ds_result.to_numpy().base, (2, 0, 1))), dims=["c", "y", "x"] + ) + dsimage.attrs["transform"] = {"global": Identity()} + # render image + stacked = np.stack([dsimage.data[c] for c in [0, 1, 2, 3]], axis=-1) + ax.imshow(stacked, zorder=render_params.zorder) + return # comment out to not show normal points + if render_params.col_for_color is not None: cols = sc.get.obs_df(adata, render_params.col_for_color) # maybe set color based on type @@ -284,7 +330,8 @@ def _render_points( cmap=render_params.cmap_params.cmap, norm=norm, alpha=render_params.alpha, - transform=trans + transform=trans, + zorder=render_params.zorder, # **kwargs, ) cax = ax.add_collection(_cax) @@ -419,6 +466,7 @@ def _render_images( im = ax.imshow( layer, cmap=cmap, + zorder=render_params.zorder, ) im.set_transform(trans_data) @@ -467,6 +515,7 @@ def _render_images( im = ax.imshow( stacked, alpha=render_params.alpha, + zorder=render_params.zorder, ) im.set_transform(trans_data) @@ -489,6 +538,7 @@ def _render_images( im = ax.imshow( colored, alpha=render_params.alpha, + zorder=render_params.zorder, ) im.set_transform(trans_data) @@ -508,6 +558,7 @@ def _render_images( im = ax.imshow( colored, alpha=render_params.alpha, + zorder=render_params.zorder, ) im.set_transform(trans_data) @@ -523,6 +574,7 @@ def _render_images( im = ax.imshow( colored, alpha=render_params.alpha, + zorder=render_params.zorder, ) im.set_transform(trans_data) @@ -642,6 +694,7 @@ def _render_labels( norm=None if categorical else render_params.cmap_params.norm, alpha=render_params.fill_alpha, origin="lower", + zorder=render_params.zorder, ) _cax.set_transform(trans_data) cax = ax.add_image(_cax) @@ -665,6 +718,7 @@ def _render_labels( norm=None if categorical else render_params.cmap_params.norm, alpha=render_params.outline_alpha, origin="lower", + zorder=render_params.zorder, ) else: # Default: no alpha, contour = infill @@ -686,6 +740,7 @@ def _render_labels( norm=None if categorical else render_params.cmap_params.norm, alpha=render_params.fill_alpha, origin="lower", + zorder=render_params.zorder, ) _cax.set_transform(trans_data) cax = ax.add_image(_cax) diff --git a/src/spatialdata_plot/pl/render_params.py b/src/spatialdata_plot/pl/render_params.py index 8eca04bf..f2e7a722 100644 --- a/src/spatialdata_plot/pl/render_params.py +++ b/src/spatialdata_plot/pl/render_params.py @@ -81,6 +81,7 @@ class ShapesRenderParams: fill_alpha: float = 0.3 scale: float = 1.0 transfunc: Callable[[float], float] | None = None + zorder: int | None = None @dataclass @@ -96,6 +97,7 @@ class PointsRenderParams: alpha: float = 1.0 size: float = 1.0 transfunc: Callable[[float], float] | None = None + zorder: int | None = None @dataclass @@ -109,6 +111,7 @@ class ImageRenderParams: alpha: float = 1.0 quantiles_for_norm: tuple[float | None, float | None] = (None, None) scale: str | list[str] | None = None + zorder: int | None = None @dataclass @@ -127,3 +130,4 @@ class LabelsRenderParams: fill_alpha: float = 0.4 transfunc: Callable[[float], float] | None = None scale: str | list[str] | None = None + zorder: int | None = None From 9f01e910ef35ad3a967c22ba7eaf6552457ec4c6 Mon Sep 17 00:00:00 2001 From: Sonja Stockhaus Date: Tue, 19 Mar 2024 22:22:04 -0400 Subject: [PATCH 02/23] add alpha & size control --- src/spatialdata_plot/pl/render.py | 19 +++++++------------ 1 file changed, 7 insertions(+), 12 deletions(-) diff --git a/src/spatialdata_plot/pl/render.py b/src/spatialdata_plot/pl/render.py index 07897cd3..3a56796b 100644 --- a/src/spatialdata_plot/pl/render.py +++ b/src/spatialdata_plot/pl/render.py @@ -14,12 +14,10 @@ import scanpy as sc import spatialdata as sd from anndata import AnnData -from dask.array import from_array from matplotlib.colors import ListedColormap, Normalize from multiscale_spatial_image.multiscale_spatial_image import MultiscaleSpatialImage from pandas.api.types import is_categorical_dtype from scanpy._settings import settings as sc_settings -from spatial_image import to_spatial_image from spatialdata._core.data_extent import get_extent from spatialdata.models import ( PointsModel, @@ -27,7 +25,6 @@ from spatialdata.transformations import ( get_transformation, ) -from spatialdata.transformations.transformations import Identity from spatialdata_plot._logging import logger from spatialdata_plot.pl.render_params import ( @@ -267,19 +264,17 @@ def _render_points( plot_height = int(np.round(y_range[1] - y_range[0])) # use datashader for the visualization of points - # TODO: incorporate color, size etc + # TODO: parameters of render_params: color, col_for_color, groups, palette, cmap_params cvs = ds.Canvas(plot_width=plot_width, plot_height=plot_height, x_range=x_range, y_range=y_range) agg = cvs.points(sdata.points[e], "x", "y") - ds_result = ds.tf.shade(ds.tf.spread(agg, px=5), rescale_discrete_levels=True) - - # TODO: we might not even need this conversion (directly use numpy array) for direct rendering - dsimage = to_spatial_image( - from_array(np.transpose(ds_result.to_numpy().base, (2, 0, 1))), dims=["c", "y", "x"] + px = int(np.round(render_params.size)) + ds_result = ds.tf.shade( + ds.tf.spread(agg, px=px), rescale_discrete_levels=True, alpha=render_params.alpha * 255 ) - dsimage.attrs["transform"] = {"global": Identity()} + # render image - stacked = np.stack([dsimage.data[c] for c in [0, 1, 2, 3]], axis=-1) - ax.imshow(stacked, zorder=render_params.zorder) + rbga_image = np.transpose(ds_result.to_numpy().base, (0, 1, 2)) + ax.imshow(rbga_image, zorder=render_params.zorder) return # comment out to not show normal points if render_params.col_for_color is not None: From 32f0e200bd4a94efef93175a7b87709ace1102ea Mon Sep 17 00:00:00 2001 From: Sonja Stockhaus Date: Thu, 21 Mar 2024 22:31:48 -0400 Subject: [PATCH 03/23] simple color control --- src/spatialdata_plot/pl/render.py | 85 +++++++++++++++++-------------- 1 file changed, 46 insertions(+), 39 deletions(-) diff --git a/src/spatialdata_plot/pl/render.py b/src/spatialdata_plot/pl/render.py index 3a56796b..75119e15 100644 --- a/src/spatialdata_plot/pl/render.py +++ b/src/spatialdata_plot/pl/render.py @@ -238,45 +238,6 @@ def _render_points( points = dask.dataframe.from_pandas(points, npartitions=1) sdata_filt.points[e] = PointsModel.parse(points, coordinates={"x": "x", "y": "y"}) - # TODO: maybe move this, add heuristic - if len(points) > 100: - extent = get_extent(sdata_filt.points[e], coordinate_system=coordinate_system) - x_ext = extent["x"][1] - y_ext = extent["y"][1] - previous_xlim = fig_params.ax.get_xlim() - previous_ylim = fig_params.ax.get_ylim() - x_range = [0, x_ext] - y_range = [0, y_ext] - # if sth else was plotted before: might need larger extent/size - # else: circles representing a point are "cut off at image borders" - if previous_ylim[0] > previous_ylim[1]: - if previous_xlim[0] < x_range[0]: - x_range[0] = previous_xlim[0] - if previous_xlim[1] > x_range[1]: - x_range[1] = previous_xlim[1] - # different order for x and y! - if previous_ylim[0] > y_range[1]: - y_range[1] = previous_ylim[0] - if previous_ylim[1] < y_range[0]: - y_range[0] = previous_ylim[1] - # round because we need integers - plot_width = int(np.round(x_range[1] - x_range[0])) - plot_height = int(np.round(y_range[1] - y_range[0])) - - # use datashader for the visualization of points - # TODO: parameters of render_params: color, col_for_color, groups, palette, cmap_params - cvs = ds.Canvas(plot_width=plot_width, plot_height=plot_height, x_range=x_range, y_range=y_range) - agg = cvs.points(sdata.points[e], "x", "y") - px = int(np.round(render_params.size)) - ds_result = ds.tf.shade( - ds.tf.spread(agg, px=px), rescale_discrete_levels=True, alpha=render_params.alpha * 255 - ) - - # render image - rbga_image = np.transpose(ds_result.to_numpy().base, (0, 1, 2)) - ax.imshow(rbga_image, zorder=render_params.zorder) - return # comment out to not show normal points - if render_params.col_for_color is not None: cols = sc.get.obs_df(adata, render_params.col_for_color) # maybe set color based on type @@ -315,6 +276,52 @@ def _render_points( affine_trans = trans.to_affine_matrix(input_axes=("x", "y"), output_axes=("x", "y")) trans = mtransforms.Affine2D(matrix=affine_trans) + ax.transData + # optionally render points using datashader + # TODO: maybe move this, add heuristic + if len(points) > 100: + extent = get_extent(sdata_filt.points[e], coordinate_system=coordinate_system) + x_ext = extent["x"][1] + y_ext = extent["y"][1] + previous_xlim = fig_params.ax.get_xlim() + previous_ylim = fig_params.ax.get_ylim() + x_range = [0, x_ext] + y_range = [0, y_ext] + # if sth else was plotted before: might need larger extent/size + # else: circles representing a point are "cut off at image borders" + if previous_ylim[0] > previous_ylim[1]: + if previous_xlim[0] < x_range[0]: + x_range[0] = previous_xlim[0] + if previous_xlim[1] > x_range[1]: + x_range[1] = previous_xlim[1] + # different order for x and y! + if previous_ylim[0] > y_range[1]: + y_range[1] = previous_ylim[0] + if previous_ylim[1] < y_range[0]: + y_range[0] = previous_ylim[1] + # round because we need integers + plot_width = int(np.round(x_range[1] - x_range[0])) + plot_height = int(np.round(y_range[1] - y_range[0])) + + # use datashader for the visualization of points + # TODO: parameters of render_params: color, col_for_color, groups, palette, cmap_params + # TODO: what about trans at this point? + cvs = ds.Canvas(plot_width=plot_width, plot_height=plot_height, x_range=x_range, y_range=y_range) + agg = cvs.points(sdata_filt.points[e], "x", "y") + px = int(np.round(render_params.size)) + if render_params.color is None: + render_params.color = "black" + ds_result = ds.tf.shade( + ds.tf.spread(agg, px=px), + rescale_discrete_levels=True, + alpha=render_params.alpha * 255, + cmap=render_params.color, # TODO: can we use color_vector somehow? + ) + + # render image + rbga_image = np.transpose(ds_result.to_numpy().base, (0, 1, 2)) + ax.imshow(rbga_image, zorder=render_params.zorder) + return # comment out to not show normal points + norm = copy(render_params.cmap_params.norm) _cax = ax.scatter( adata[:, 0].X.flatten(), From 0685c357fefc939c53e637e557fe73f135f088f0 Mon Sep 17 00:00:00 2001 From: Sonja Stockhaus Date: Thu, 28 Mar 2024 21:53:50 -0400 Subject: [PATCH 04/23] update color handling --- src/spatialdata_plot/pl/render.py | 61 ++++++++++++++++++------------- 1 file changed, 36 insertions(+), 25 deletions(-) diff --git a/src/spatialdata_plot/pl/render.py b/src/spatialdata_plot/pl/render.py index 75119e15..1fa493b6 100644 --- a/src/spatialdata_plot/pl/render.py +++ b/src/spatialdata_plot/pl/render.py @@ -276,9 +276,11 @@ def _render_points( affine_trans = trans.to_affine_matrix(input_axes=("x", "y"), output_axes=("x", "y")) trans = mtransforms.Affine2D(matrix=affine_trans) + ax.transData + norm = copy(render_params.cmap_params.norm) + # optionally render points using datashader # TODO: maybe move this, add heuristic - if len(points) > 100: + if len(points) > 50: extent = get_extent(sdata_filt.points[e], coordinate_system=coordinate_system) x_ext = extent["x"][1] y_ext = extent["y"][1] @@ -303,40 +305,49 @@ def _render_points( plot_height = int(np.round(y_range[1] - y_range[0])) # use datashader for the visualization of points - # TODO: parameters of render_params: color, col_for_color, groups, palette, cmap_params - # TODO: what about trans at this point? + # TODO: parameters of render_params: color, col_for_color, palette, cmap_params + # TODO: what about trans/norm at this point? cvs = ds.Canvas(plot_width=plot_width, plot_height=plot_height, x_range=x_range, y_range=y_range) - agg = cvs.points(sdata_filt.points[e], "x", "y") + # TODO: design choice: do we want count as aggregate? for many overlapping points it could be cool + # currently it would be enough to just go full opacity when there is a point + if col_for_color is not None: + agg = cvs.points(sdata_filt.points[e], "x", "y", agg=ds.by(col_for_color, ds.count())) + else: + agg = cvs.points(sdata_filt.points[e], "x", "y") px = int(np.round(render_params.size)) - if render_params.color is None: - render_params.color = "black" + + color_key = ( + [x[:-2] for x in color_vector.categories.values] + if type(color_vector) == pd.core.arrays.categorical.Categorical + else None + ) # TODO: is order always correct? ds_result = ds.tf.shade( ds.tf.spread(agg, px=px), rescale_discrete_levels=True, alpha=render_params.alpha * 255, - cmap=render_params.color, # TODO: can we use color_vector somehow? + cmap=color_vector[0][:-2], + color_key=color_key, ) - # render image rbga_image = np.transpose(ds_result.to_numpy().base, (0, 1, 2)) ax.imshow(rbga_image, zorder=render_params.zorder) - return # comment out to not show normal points - - norm = copy(render_params.cmap_params.norm) - _cax = ax.scatter( - adata[:, 0].X.flatten(), - adata[:, 1].X.flatten(), - s=render_params.size, - c=color_vector, - rasterized=sc_settings._vector_friendly, - cmap=render_params.cmap_params.cmap, - norm=norm, - alpha=render_params.alpha, - transform=trans, - zorder=render_params.zorder, - # **kwargs, - ) - cax = ax.add_collection(_cax) + cax = None + else: + # original way of plotting points + _cax = ax.scatter( + adata[:, 0].X.flatten(), + adata[:, 1].X.flatten(), + s=render_params.size, + c=color_vector, + rasterized=sc_settings._vector_friendly, + cmap=render_params.cmap_params.cmap, + norm=norm, + alpha=render_params.alpha, + transform=trans, + zorder=render_params.zorder, + # **kwargs, + ) + cax = ax.add_collection(_cax) if len(set(color_vector)) != 1 or list(set(color_vector))[0] != to_hex(render_params.cmap_params.na_color): if color_source_vector is None: From defb86b70aa5be9467454763ef468ec7b3a7db7d Mon Sep 17 00:00:00 2001 From: Sonja Stockhaus Date: Fri, 29 Mar 2024 22:02:57 -0400 Subject: [PATCH 05/23] first finished draft --- src/spatialdata_plot/pl/render.py | 19 ++++++++++--------- 1 file changed, 10 insertions(+), 9 deletions(-) diff --git a/src/spatialdata_plot/pl/render.py b/src/spatialdata_plot/pl/render.py index 1fa493b6..55bffd1f 100644 --- a/src/spatialdata_plot/pl/render.py +++ b/src/spatialdata_plot/pl/render.py @@ -305,22 +305,23 @@ def _render_points( plot_height = int(np.round(y_range[1] - y_range[0])) # use datashader for the visualization of points - # TODO: parameters of render_params: color, col_for_color, palette, cmap_params # TODO: what about trans/norm at this point? cvs = ds.Canvas(plot_width=plot_width, plot_height=plot_height, x_range=x_range, y_range=y_range) - # TODO: design choice: do we want count as aggregate? for many overlapping points it could be cool - # currently it would be enough to just go full opacity when there is a point - if col_for_color is not None: + # TODO: design choice: do we want count as reduction? + if col_for_color is not None and (render_params.groups is None or len(render_params.groups) > 1): agg = cvs.points(sdata_filt.points[e], "x", "y", agg=ds.by(col_for_color, ds.count())) + # agg = cvs.points(sdata_filt.points[e], "x", "y", agg=ds.by(col_for_color, ds.any())) else: - agg = cvs.points(sdata_filt.points[e], "x", "y") - px = int(np.round(render_params.size)) - + agg = cvs.points(sdata_filt.points[e], "x", "y", agg=ds.count()) + # agg = cvs.points(sdata_filt.points[e], "x", "y", agg=ds.any()) + # NOTE: s in matplotlib is in units of points**2 + px = int(np.round(np.sqrt(render_params.size))) color_key = ( [x[:-2] for x in color_vector.categories.values] - if type(color_vector) == pd.core.arrays.categorical.Categorical + if (type(color_vector) == pd.core.arrays.categorical.Categorical) + and (len(color_vector.categories.values) > 1) else None - ) # TODO: is order always correct? + ) ds_result = ds.tf.shade( ds.tf.spread(agg, px=px), rescale_discrete_levels=True, From 3bfab39233809725ed65cc8e8afc90bd6045934c Mon Sep 17 00:00:00 2001 From: Sonja Stockhaus Date: Mon, 1 Apr 2024 22:05:41 -0400 Subject: [PATCH 06/23] add min_alpha to avoid undersaturation --- src/spatialdata_plot/pl/render.py | 1 + 1 file changed, 1 insertion(+) diff --git a/src/spatialdata_plot/pl/render.py b/src/spatialdata_plot/pl/render.py index 55bffd1f..0e1fef21 100644 --- a/src/spatialdata_plot/pl/render.py +++ b/src/spatialdata_plot/pl/render.py @@ -328,6 +328,7 @@ def _render_points( alpha=render_params.alpha * 255, cmap=color_vector[0][:-2], color_key=color_key, + min_alpha=200, # TODO: choose some value (any gives wierd background) ) # render image rbga_image = np.transpose(ds_result.to_numpy().base, (0, 1, 2)) From bcecb0eca6ebfbd10e1bf6543a211a558e2b43e3 Mon Sep 17 00:00:00 2001 From: ckmah Date: Wed, 3 Apr 2024 08:37:59 -0700 Subject: [PATCH 07/23] basic shapes datashader rendering --- src/spatialdata_plot/pl/render.py | 92 +++++++++++++++++++++---------- 1 file changed, 63 insertions(+), 29 deletions(-) diff --git a/src/spatialdata_plot/pl/render.py b/src/spatialdata_plot/pl/render.py index 0e1fef21..29f00805 100644 --- a/src/spatialdata_plot/pl/render.py +++ b/src/spatialdata_plot/pl/render.py @@ -122,35 +122,6 @@ def _render_shapes( shapes = shapes.reset_index() color_source_vector = color_source_vector[mask] color_vector = color_vector[mask] - shapes = gpd.GeoDataFrame(shapes, geometry="geometry") - - _cax = _get_collection_shape( - shapes=shapes, - s=render_params.scale, - c=color_vector, - render_params=render_params, - rasterized=sc_settings._vector_friendly, - cmap=render_params.cmap_params.cmap, - norm=norm, - fill_alpha=render_params.fill_alpha, - outline_alpha=render_params.outline_alpha, - zorder=render_params.zorder, - # **kwargs, - ) - - # Sets the limits of the colorbar to the values instead of [0, 1] - if not norm and not values_are_categorical: - _cax.set_clim(min(color_vector), max(color_vector)) - - cax = ax.add_collection(_cax) - - # Apply the transformation to the PatchCollection's paths - trans = get_transformation(sdata_filt.shapes[e], get_all=True)[coordinate_system] - affine_trans = trans.to_affine_matrix(input_axes=("x", "y"), output_axes=("x", "y")) - trans = mtransforms.Affine2D(matrix=affine_trans) - - for path in _cax.get_paths(): - path.vertices = trans.transform(path.vertices) # Using dict.fromkeys here since set returns in arbitrary order # remove the color of NaN values, else it might be assigned to a category @@ -160,6 +131,69 @@ def _render_shapes( else: palette = ListedColormap(dict.fromkeys(color_vector[~pd.Categorical(color_source_vector).isnull()])) + # Apply the transformation to the PatchCollection's paths + trans = get_transformation(sdata_filt.shapes[e], get_all=True)[coordinate_system] + affine_trans = trans.to_affine_matrix(input_axes=("x", "y"), output_axes=("x", "y")) + trans = mtransforms.Affine2D(matrix=affine_trans) + ax.transData + + shapes = gpd.GeoDataFrame(shapes, geometry="geometry") + + if len(shapes) < 1: + logger.info("Using matplotlib") + _cax = _get_collection_shape( + shapes=shapes, + s=render_params.scale, + c=color_vector, + render_params=render_params, + rasterized=sc_settings._vector_friendly, + cmap=render_params.cmap_params.cmap, + norm=norm, + fill_alpha=render_params.fill_alpha, + outline_alpha=render_params.outline_alpha, + zorder=render_params.zorder, + # **kwargs, + ) + cax = ax.add_collection(_cax) + + # Transform the paths in PatchCollection + for path in _cax.get_paths(): + path.vertices = trans.transform(path.vertices) + cax = ax.add_collection(_cax) + else: + logger.info("Using datashader") + extent = get_extent(sdata.shapes[e]) + x_ext = extent["x"][1] + y_ext = extent["y"][1] + # previous_xlim = fig_params.ax.get_xlim() + # previous_ylim = fig_params.ax.get_ylim() + x_range = [0, x_ext] + y_range = [0, y_ext] + # round because we need integers + plot_width = int(np.round(x_range[1] - x_range[0])) + plot_height = int(np.round(y_range[1] - y_range[0])) + + cvs = ds.Canvas(plot_width=plot_width, plot_height=plot_height, x_range=x_range, y_range=y_range) + + _geometry = shapes["geometry"] + is_point = _geometry.type == "Point" + + # Handle circles encoded as points with radius + if is_point.any(): # TODO + shapes.loc[is_point, "geometry"] = _geometry[is_point].buffer(shapes[is_point]["radius"]) + + agg = cvs.polygons(shapes, geometry="geometry", agg=ds.count()) + ds_result = ds.tf.shade(agg) + + # Render image + rgba_image = np.transpose(ds_result.to_numpy().base, (0, 1, 2)) + _cax = ax.imshow(rgba_image, cmap=palette, zorder=render_params.zorder) + _cax.set_transform(trans) + cax = ax.add_image(_cax) + + # Sets the limits of the colorbar to the values instead of [0, 1] + if not norm and not values_are_categorical: + _cax.set_clim(min(color_vector), max(color_vector)) + if not ( len(set(color_vector)) == 1 and list(set(color_vector))[0] == to_hex(render_params.cmap_params.na_color) ): From cc02f62a00b301d9261c02dbfbd030f3c165a302 Mon Sep 17 00:00:00 2001 From: ckmah Date: Wed, 3 Apr 2024 12:06:09 -0700 Subject: [PATCH 08/23] handle render shape colors; add parameter for render method --- src/spatialdata_plot/pl/basic.py | 22 +++++++++ src/spatialdata_plot/pl/render.py | 57 +++++++++++++++++++----- src/spatialdata_plot/pl/render_params.py | 2 + 3 files changed, 70 insertions(+), 11 deletions(-) diff --git a/src/spatialdata_plot/pl/basic.py b/src/spatialdata_plot/pl/basic.py index a13e1366..19f4eb3a 100644 --- a/src/spatialdata_plot/pl/basic.py +++ b/src/spatialdata_plot/pl/basic.py @@ -160,6 +160,7 @@ def render_shapes( cmap: Colormap | str | None = None, norm: bool | Normalize = False, scale: float | int = 1.0, + method: str | None = None, **kwargs: Any, ) -> sd.SpatialData: """ @@ -204,6 +205,9 @@ def render_shapes( Colormap normalization for continuous annotations. scale : float | int, default 1.0 Value to scale circles, if present. + method : str | None, optional + Whether to use 'matplotlib' and 'datashader'. When None, the method is + chosen based on the size of the data. **kwargs : Any Additional arguments to be passed to cmap and norm. @@ -317,6 +321,12 @@ def render_shapes( if scale < 0: raise ValueError("Parameter 'scale' must be a positive number.") + if method is not None: + if not isinstance(method, str): + raise TypeError("Parameter 'method' must be a string.") + if method not in ["matplotlib", "datashader"]: + raise ValueError("Parameter 'method' must be either 'matplotlib' or 'datashader'.") + sdata = self._copy() sdata = _verify_plotting_tree(sdata) n_steps = len(sdata.plotting_tree.keys()) @@ -343,6 +353,7 @@ def render_shapes( fill_alpha=fill_alpha, transfunc=kwargs.get("transfunc", None), zorder=n_steps, + method=method, ) return sdata @@ -358,6 +369,7 @@ def render_points( cmap: Colormap | str | None = None, norm: None | Normalize = None, size: float | int = 1.0, + method: str | None = None, **kwargs: Any, ) -> sd.SpatialData: """ @@ -392,6 +404,9 @@ def render_points( Colormap normalization for continuous annotations. size : float | int, default 1.0 Size of the points + method : str | None, optional + Whether to use 'matplotlib' and 'datashader'. When None, the method is + chosen based on the size of the data. kwargs Additional arguments to be passed to cmap and norm. @@ -479,6 +494,12 @@ def render_points( if size < 0: raise ValueError("Parameter 'size' must be a positive number.") + if method is not None: + if not isinstance(method, str): + raise TypeError("Parameter 'method' must be a string.") + if method not in ["matplotlib", "datashader"]: + raise ValueError("Parameter 'method' must be either 'matplotlib' or 'datashader'.") + sdata = self._copy() sdata = _verify_plotting_tree(sdata) n_steps = len(sdata.plotting_tree.keys()) @@ -501,6 +522,7 @@ def render_points( transfunc=kwargs.get("transfunc", None), size=size, zorder=n_steps, + method=method, ) return sdata diff --git a/src/spatialdata_plot/pl/render.py b/src/spatialdata_plot/pl/render.py index 29f00805..4ed4f4e5 100644 --- a/src/spatialdata_plot/pl/render.py +++ b/src/spatialdata_plot/pl/render.py @@ -134,12 +134,21 @@ def _render_shapes( # Apply the transformation to the PatchCollection's paths trans = get_transformation(sdata_filt.shapes[e], get_all=True)[coordinate_system] affine_trans = trans.to_affine_matrix(input_axes=("x", "y"), output_axes=("x", "y")) - trans = mtransforms.Affine2D(matrix=affine_trans) + ax.transData + trans = mtransforms.Affine2D(matrix=affine_trans) shapes = gpd.GeoDataFrame(shapes, geometry="geometry") - if len(shapes) < 1: - logger.info("Using matplotlib") + # Determine which method to use for rendering. Default is matplotlib for under 100 shapes and datashader for more + # User can also specify the method to use + method = render_params.method + + if method is None: + method = "datashader" if len(shapes) > 100 else "matplotlib" + elif method not in ["matplotlib", "datashader"]: + raise ValueError("Method must be either 'matplotlib' or 'datashader'.") + + if method == "matplotlib": + logger.info(f"Using {method}") _cax = _get_collection_shape( shapes=shapes, s=render_params.scale, @@ -159,8 +168,12 @@ def _render_shapes( for path in _cax.get_paths(): path.vertices = trans.transform(path.vertices) cax = ax.add_collection(_cax) - else: - logger.info("Using datashader") + elif method == "datashader": + logger.info(f"Using {method}") + + # Where to put this + trans = mtransforms.Affine2D(matrix=affine_trans) + ax.transData + extent = get_extent(sdata.shapes[e]) x_ext = extent["x"][1] y_ext = extent["y"][1] @@ -179,10 +192,28 @@ def _render_shapes( # Handle circles encoded as points with radius if is_point.any(): # TODO - shapes.loc[is_point, "geometry"] = _geometry[is_point].buffer(shapes[is_point]["radius"]) + scale = shapes[is_point]["radius"] * render_params.scale + shapes.loc[is_point, "geometry"] = _geometry[is_point].buffer(scale) agg = cvs.polygons(shapes, geometry="geometry", agg=ds.count()) - ds_result = ds.tf.shade(agg) + + # Render shapes with datashader + if render_params.col_for_color is not None and ( + render_params.groups is None or len(render_params.groups) > 1 + ): + agg = cvs.polygons(shapes, geometry="geometry", agg=ds.by(render_params.col_for_color, ds.count())) + else: + agg = cvs.polygons(shapes, geometry="geometry", agg=ds.count()) + + color_key = ( + [x[:-2] for x in color_vector.categories.values] + if (type(color_vector) == pd.core.arrays.categorical.Categorical) + and (len(color_vector.categories.values) > 1) + else None + ) + ds_result = ds.tf.shade( + agg, cmap=color_vector[0][:-2], alpha=render_params.fill_alpha * 255, color_key=color_key, min_alpha=200 + ) # Render image rgba_image = np.transpose(ds_result.to_numpy().base, (0, 1, 2)) @@ -312,9 +343,13 @@ def _render_points( norm = copy(render_params.cmap_params.norm) - # optionally render points using datashader - # TODO: maybe move this, add heuristic - if len(points) > 50: + method = render_params.method + if method is None: + method = "datashader" if len(points.shape[0]) > 10000 else "matplotlib" + elif method not in ["matplotlib", "datashader"]: + raise ValueError("Method must be either 'matplotlib' or 'datashader'.") + + if method == "datashader": extent = get_extent(sdata_filt.points[e], coordinate_system=coordinate_system) x_ext = extent["x"][1] y_ext = extent["y"][1] @@ -368,7 +403,7 @@ def _render_points( rbga_image = np.transpose(ds_result.to_numpy().base, (0, 1, 2)) ax.imshow(rbga_image, zorder=render_params.zorder) cax = None - else: + elif method == "matplotlib": # original way of plotting points _cax = ax.scatter( adata[:, 0].X.flatten(), diff --git a/src/spatialdata_plot/pl/render_params.py b/src/spatialdata_plot/pl/render_params.py index f2e7a722..39439887 100644 --- a/src/spatialdata_plot/pl/render_params.py +++ b/src/spatialdata_plot/pl/render_params.py @@ -81,6 +81,7 @@ class ShapesRenderParams: fill_alpha: float = 0.3 scale: float = 1.0 transfunc: Callable[[float], float] | None = None + method: str | None = None zorder: int | None = None @@ -97,6 +98,7 @@ class PointsRenderParams: alpha: float = 1.0 size: float = 1.0 transfunc: Callable[[float], float] | None = None + method: str | None = None zorder: int | None = None From 8e809f53a17273228110e699ce13a729d66bea4a Mon Sep 17 00:00:00 2001 From: ckmah Date: Wed, 3 Apr 2024 12:09:14 -0700 Subject: [PATCH 09/23] cleanup --- src/spatialdata_plot/pl/render.py | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/src/spatialdata_plot/pl/render.py b/src/spatialdata_plot/pl/render.py index 4ed4f4e5..6ee60291 100644 --- a/src/spatialdata_plot/pl/render.py +++ b/src/spatialdata_plot/pl/render.py @@ -138,10 +138,8 @@ def _render_shapes( shapes = gpd.GeoDataFrame(shapes, geometry="geometry") - # Determine which method to use for rendering. Default is matplotlib for under 100 shapes and datashader for more - # User can also specify the method to use + # Determine which method to use for rendering method = render_params.method - if method is None: method = "datashader" if len(shapes) > 100 else "matplotlib" elif method not in ["matplotlib", "datashader"]: From 9b39852046c5fae6127ef8a011d23d3eaa2adef2 Mon Sep 17 00:00:00 2001 From: Sonja Stockhaus Date: Tue, 30 Apr 2024 20:51:51 -0400 Subject: [PATCH 10/23] extent fixes, sum aggregation for coloring by numerical variable --- src/spatialdata_plot/pl/render.py | 104 ++++++++++++++++++------------ 1 file changed, 62 insertions(+), 42 deletions(-) diff --git a/src/spatialdata_plot/pl/render.py b/src/spatialdata_plot/pl/render.py index 6ee60291..145e1307 100644 --- a/src/spatialdata_plot/pl/render.py +++ b/src/spatialdata_plot/pl/render.py @@ -14,6 +14,7 @@ import scanpy as sc import spatialdata as sd from anndata import AnnData +from matplotlib.cm import ScalarMappable from matplotlib.colors import ListedColormap, Normalize from multiscale_spatial_image.multiscale_spatial_image import MultiscaleSpatialImage from pandas.api.types import is_categorical_dtype @@ -43,6 +44,7 @@ _get_linear_colormap, _map_color_seg, _maybe_set_colors, + _mpl_ax_contains_elements, _multiscale_to_spatial_image, _normalize, _rasterize_if_necessary, @@ -141,12 +143,12 @@ def _render_shapes( # Determine which method to use for rendering method = render_params.method if method is None: - method = "datashader" if len(shapes) > 100 else "matplotlib" + method = "datashader" if len(shapes) > 10000 else "matplotlib" elif method not in ["matplotlib", "datashader"]: raise ValueError("Method must be either 'matplotlib' or 'datashader'.") + logger.info(f"Using {method}") if method == "matplotlib": - logger.info(f"Using {method}") _cax = _get_collection_shape( shapes=shapes, s=render_params.scale, @@ -167,8 +169,6 @@ def _render_shapes( path.vertices = trans.transform(path.vertices) cax = ax.add_collection(_cax) elif method == "datashader": - logger.info(f"Using {method}") - # Where to put this trans = mtransforms.Affine2D(matrix=affine_trans) + ax.transData @@ -343,65 +343,80 @@ def _render_points( method = render_params.method if method is None: - method = "datashader" if len(points.shape[0]) > 10000 else "matplotlib" + method = "datashader" if len(points) > 10000 else "matplotlib" elif method not in ["matplotlib", "datashader"]: raise ValueError("Method must be either 'matplotlib' or 'datashader'.") if method == "datashader": + # NOTE: s in matplotlib is in units of points**2 + px = int(np.round(np.sqrt(render_params.size))) + extent = get_extent(sdata_filt.points[e], coordinate_system=coordinate_system) - x_ext = extent["x"][1] - y_ext = extent["y"][1] - previous_xlim = fig_params.ax.get_xlim() - previous_ylim = fig_params.ax.get_ylim() - x_range = [0, x_ext] - y_range = [0, y_ext] - # if sth else was plotted before: might need larger extent/size - # else: circles representing a point are "cut off at image borders" - if previous_ylim[0] > previous_ylim[1]: - if previous_xlim[0] < x_range[0]: - x_range[0] = previous_xlim[0] - if previous_xlim[1] > x_range[1]: - x_range[1] = previous_xlim[1] - # different order for x and y! - if previous_ylim[0] > y_range[1]: - y_range[1] = previous_ylim[0] - if previous_ylim[1] < y_range[0]: - y_range[0] = previous_ylim[1] + x_ext = [min(0, extent["x"][0]), extent["x"][1]] + y_ext = [min(0, extent["y"][0]), extent["y"][1]] + previous_xlim = ax.get_xlim() + previous_ylim = ax.get_ylim() + # increase range if sth larger was rendered before + if _mpl_ax_contains_elements(ax): + x_ext = [min(x_ext[0], previous_xlim[0]), max(x_ext[1], previous_xlim[1])] + if ax.yaxis_inverted(): # case for e.g. images + y_ext = [min(y_ext[0], previous_ylim[1]), max(y_ext[1], previous_ylim[0])] + else: # case for e.g. labels + y_ext = [min(y_ext[0], previous_ylim[0]), max(y_ext[1], previous_ylim[1])] # round because we need integers - plot_width = int(np.round(x_range[1] - x_range[0])) - plot_height = int(np.round(y_range[1] - y_range[0])) + plot_width = int(np.round(x_ext[1] - x_ext[0])) + plot_height = int(np.round(y_ext[1] - y_ext[0])) # use datashader for the visualization of points # TODO: what about trans/norm at this point? - cvs = ds.Canvas(plot_width=plot_width, plot_height=plot_height, x_range=x_range, y_range=y_range) - # TODO: design choice: do we want count as reduction? + cvs = ds.Canvas(plot_width=plot_width, plot_height=plot_height, x_range=x_ext, y_range=y_ext) + + color_by_categorical = col_for_color is not None and points[col_for_color].values.dtype == object + aggregate_with_sum = None if col_for_color is not None and (render_params.groups is None or len(render_params.groups) > 1): - agg = cvs.points(sdata_filt.points[e], "x", "y", agg=ds.by(col_for_color, ds.count())) - # agg = cvs.points(sdata_filt.points[e], "x", "y", agg=ds.by(col_for_color, ds.any())) + if color_by_categorical: + agg = cvs.points(sdata_filt.points[e], "x", "y", agg=ds.by(col_for_color, ds.count())) + else: + # numerical + agg = cvs.points(sdata_filt.points[e], "x", "y", agg=ds.sum(column=col_for_color)) + # save min and max values for drawing the colorbar + aggregate_with_sum = (agg.min(), agg.max()) else: agg = cvs.points(sdata_filt.points[e], "x", "y", agg=ds.count()) - # agg = cvs.points(sdata_filt.points[e], "x", "y", agg=ds.any()) - # NOTE: s in matplotlib is in units of points**2 - px = int(np.round(np.sqrt(render_params.size))) + color_key = ( [x[:-2] for x in color_vector.categories.values] if (type(color_vector) == pd.core.arrays.categorical.Categorical) and (len(color_vector.categories.values) > 1) else None ) - ds_result = ds.tf.shade( - ds.tf.spread(agg, px=px), - rescale_discrete_levels=True, - alpha=render_params.alpha * 255, - cmap=color_vector[0][:-2], - color_key=color_key, - min_alpha=200, # TODO: choose some value (any gives wierd background) - ) + if color_by_categorical or col_for_color is None: + ds_result = ds.tf.shade( + ds.tf.spread(agg, px=px), + rescale_discrete_levels=True, + cmap=color_vector[0][:-2], + color_key=color_key, + min_alpha=np.min([150, render_params.alpha * 255]), + ) # TODO: choose other value than 150 for min_alpha (here and below)? + else: + ds_result = ds.tf.shade( + ds.tf.spread(agg, px=px), + rescale_discrete_levels=True, + cmap=render_params.cmap_params.cmap, + # color_key=color_key, + ) # render image rbga_image = np.transpose(ds_result.to_numpy().base, (0, 1, 2)) - ax.imshow(rbga_image, zorder=render_params.zorder) - cax = None + cax = ax.imshow(rbga_image, zorder=render_params.zorder, alpha=render_params.alpha) + if aggregate_with_sum is not None: + cax = ScalarMappable( + norm=matplotlib.colors.Normalize(vmin=aggregate_with_sum[0], vmax=aggregate_with_sum[1]), + cmap=render_params.cmap_params.cmap, + ) + elif method == "matplotlib": + # update axis limits if plot was empty before (necessary if datashader comes after) + update_parameters = not _mpl_ax_contains_elements(ax) # original way of plotting points _cax = ax.scatter( adata[:, 0].X.flatten(), @@ -417,6 +432,11 @@ def _render_points( # **kwargs, ) cax = ax.add_collection(_cax) + if update_parameters: + # necessary if points are plotted with mpl first and then with datashader + extent = get_extent(sdata_filt.points[e], coordinate_system=coordinate_system) + ax.set_xbound(extent["x"]) + ax.set_ybound(extent["y"]) if len(set(color_vector)) != 1 or list(set(color_vector))[0] != to_hex(render_params.cmap_params.na_color): if color_source_vector is None: From 31cee7de95db06565b55635879fc3d1a28e4ed6a Mon Sep 17 00:00:00 2001 From: Sonja Stockhaus Date: Thu, 20 Jun 2024 18:25:30 +0200 Subject: [PATCH 11/23] post merge with main adaptations --- src/spatialdata_plot/pl/basic.py | 18 +- src/spatialdata_plot/pl/render.py | 342 ++++++++++++++---------------- src/spatialdata_plot/pl/utils.py | 3 +- 3 files changed, 167 insertions(+), 196 deletions(-) diff --git a/src/spatialdata_plot/pl/basic.py b/src/spatialdata_plot/pl/basic.py index ab648cb4..ca64d6cd 100644 --- a/src/spatialdata_plot/pl/basic.py +++ b/src/spatialdata_plot/pl/basic.py @@ -21,7 +21,7 @@ from multiscale_spatial_image.multiscale_spatial_image import MultiscaleSpatialImage from spatial_image import SpatialImage from spatialdata._core.data_extent import get_extent -from spatialdata._utils import deprecation_alias +from spatialdata._utils import _deprecation_alias from spatialdata_plot._accessor import register_spatial_data_accessor from spatialdata_plot.pl.render import ( @@ -150,7 +150,7 @@ def _copy( return sdata - @deprecation_alias(elements="element", version="0.3.0") + @_deprecation_alias(elements="element", version="0.3.0") def render_shapes( self, element: str | None = None, @@ -291,12 +291,14 @@ def render_shapes( fill_alpha=param_values["fill_alpha"], transfunc=kwargs.get("transfunc", None), table_name=param_values["table_name"], + zorder=n_steps, + method=method, ) n_steps += 1 return sdata - @deprecation_alias(elements="element", version="0.3.0") + @_deprecation_alias(elements="element", version="0.3.0") def render_points( self, element: str | None = None, @@ -411,12 +413,14 @@ def render_points( transfunc=kwargs.get("transfunc", None), size=param_values["size"], table_name=param_values["table_name"], + zorder=n_steps, + method=method, ) n_steps += 1 return sdata - @deprecation_alias(elements="element", quantiles_for_norm="percentiles_for_norm", version="version 0.3.0") + @_deprecation_alias(elements="element", quantiles_for_norm="percentiles_for_norm", version="version 0.3.0") def render_images( self, element: str | None = None, @@ -524,12 +528,13 @@ def render_images( alpha=param_values["alpha"], percentiles_for_norm=param_values["percentiles_for_norm"], scale=param_values["scale"], + zorder=n_steps, ) n_steps += 1 return sdata - @deprecation_alias(elements="element", version="0.3.0") + @_deprecation_alias(elements="element", version="0.3.0") def render_labels( self, element: str | None = None, @@ -542,7 +547,7 @@ def render_labels( norm: Normalize | None = None, na_color: ColorLike | None = (0.0, 0.0, 0.0, 0.0), outline_alpha: float | int = 1.0, - fill_alpha: float | int = 0.3, + fill_alpha: float | int = 0.4, scale: str | None = None, table_name: str | None = None, **kwargs: Any, @@ -646,6 +651,7 @@ def render_labels( transfunc=kwargs.get("transfunc", None), scale=param_values["scale"], table_name=param_values["table_name"], + zorder=n_steps, ) n_steps += 1 return sdata diff --git a/src/spatialdata_plot/pl/render.py b/src/spatialdata_plot/pl/render.py index 9911bdd7..9564ad98 100644 --- a/src/spatialdata_plot/pl/render.py +++ b/src/spatialdata_plot/pl/render.py @@ -15,6 +15,7 @@ import scanpy as sc import spatialdata as sd from anndata import AnnData +from datatree.datatree import DataTree from matplotlib.cm import ScalarMappable from matplotlib.colors import ListedColormap, Normalize from multiscale_spatial_image.multiscale_spatial_image import MultiscaleSpatialImage @@ -193,7 +194,7 @@ def _render_shapes( colorbar = False if col_for_color is None else legend_params.colorbar # Apply the transformation to the PatchCollection's paths - trans = get_transformation(sdata_filt.shapes[e], get_all=True)[coordinate_system] + trans = get_transformation(sdata_filt.shapes[element], get_all=True)[coordinate_system] affine_trans = trans.to_affine_matrix(input_axes=("x", "y"), output_axes=("x", "y")) trans = mtransforms.Affine2D(matrix=affine_trans) @@ -231,7 +232,7 @@ def _render_shapes( # Where to put this trans = mtransforms.Affine2D(matrix=affine_trans) + ax.transData - extent = get_extent(sdata.shapes[e]) + extent = get_extent(sdata.shapes[element]) x_ext = extent["x"][1] y_ext = extent["y"][1] # previous_xlim = fig_params.ax.get_xlim() @@ -297,7 +298,7 @@ def _render_shapes( cax=cax, fig_params=fig_params, adata=table, - value_to_plot=render_params.col_for_color, + value_to_plot=col_for_color, color_source_vector=color_source_vector, palette=palette, alpha=render_params.fill_alpha, @@ -409,108 +410,88 @@ def _render_points( affine_trans = trans.to_affine_matrix(input_axes=("x", "y"), output_axes=("x", "y")) trans = mtransforms.Affine2D(matrix=affine_trans) + ax.transData - norm = copy(render_params.cmap_params.norm) - _cax = ax.scatter( - adata[:, 0].X.flatten(), - adata[:, 1].X.flatten(), - s=render_params.size, - c=color_vector, - rasterized=sc_settings._vector_friendly, - cmap=render_params.cmap_params.cmap, - norm=norm, - alpha=render_params.alpha, - transform=trans, - ) - cax = ax.add_collection(_cax) - - if len(set(color_vector)) != 1 or list(set(color_vector))[0] != to_hex(render_params.cmap_params.na_color): - if color_source_vector is None: - palette = ListedColormap(dict.fromkeys(color_vector)) - else: - palette = ListedColormap(dict.fromkeys(color_vector[~pd.Categorical(color_source_vector).isnull()])) - - # color_source_vector is None when the values aren't categorical - if color_source_vector is None and render_params.transfunc is not None: - color_vector = render_params.transfunc(color_vector) - - trans = get_transformation(sdata.points[e], get_all=True)[coordinate_system] - affine_trans = trans.to_affine_matrix(input_axes=("x", "y"), output_axes=("x", "y")) - trans = mtransforms.Affine2D(matrix=affine_trans) + ax.transData + # color_source_vector is None when the values aren't categorical + if color_source_vector is None and render_params.transfunc is not None: + color_vector = render_params.transfunc(color_vector) - norm = copy(render_params.cmap_params.norm) + trans = get_transformation(sdata.points[element], get_all=True)[coordinate_system] + affine_trans = trans.to_affine_matrix(input_axes=("x", "y"), output_axes=("x", "y")) + trans = mtransforms.Affine2D(matrix=affine_trans) + ax.transData - method = render_params.method - if method is None: - method = "datashader" if len(points) > 10000 else "matplotlib" - elif method not in ["matplotlib", "datashader"]: - raise ValueError("Method must be either 'matplotlib' or 'datashader'.") + norm = copy(render_params.cmap_params.norm) - if method == "datashader": - # NOTE: s in matplotlib is in units of points**2 - px = int(np.round(np.sqrt(render_params.size))) - - extent = get_extent(sdata_filt.points[e], coordinate_system=coordinate_system) - x_ext = [min(0, extent["x"][0]), extent["x"][1]] - y_ext = [min(0, extent["y"][0]), extent["y"][1]] - previous_xlim = ax.get_xlim() - previous_ylim = ax.get_ylim() - # increase range if sth larger was rendered before - if _mpl_ax_contains_elements(ax): - x_ext = [min(x_ext[0], previous_xlim[0]), max(x_ext[1], previous_xlim[1])] - if ax.yaxis_inverted(): # case for e.g. images - y_ext = [min(y_ext[0], previous_ylim[1]), max(y_ext[1], previous_ylim[0])] - else: # case for e.g. labels - y_ext = [min(y_ext[0], previous_ylim[0]), max(y_ext[1], previous_ylim[1])] - # round because we need integers - plot_width = int(np.round(x_ext[1] - x_ext[0])) - plot_height = int(np.round(y_ext[1] - y_ext[0])) - - # use datashader for the visualization of points - # TODO: what about trans/norm at this point? - cvs = ds.Canvas(plot_width=plot_width, plot_height=plot_height, x_range=x_ext, y_range=y_ext) - - color_by_categorical = col_for_color is not None and points[col_for_color].values.dtype == object - aggregate_with_sum = None - if col_for_color is not None and (render_params.groups is None or len(render_params.groups) > 1): - if color_by_categorical: - agg = cvs.points(sdata_filt.points[e], "x", "y", agg=ds.by(col_for_color, ds.count())) - else: - # numerical - agg = cvs.points(sdata_filt.points[e], "x", "y", agg=ds.sum(column=col_for_color)) - # save min and max values for drawing the colorbar - aggregate_with_sum = (agg.min(), agg.max()) + method = render_params.method + if method is None: + method = "datashader" if len(points) > 10000 else "matplotlib" + elif method not in ["matplotlib", "datashader"]: + raise ValueError("Method must be either 'matplotlib' or 'datashader'.") + + if method == "datashader": + # NOTE: s in matplotlib is in units of points**2 + px = int(np.round(np.sqrt(render_params.size))) + + extent = get_extent(sdata_filt.points[element], coordinate_system=coordinate_system) + x_ext = [min(0, extent["x"][0]), extent["x"][1]] + y_ext = [min(0, extent["y"][0]), extent["y"][1]] + previous_xlim = ax.get_xlim() + previous_ylim = ax.get_ylim() + # increase range if sth larger was rendered before + if _mpl_ax_contains_elements(ax): + x_ext = [min(x_ext[0], previous_xlim[0]), max(x_ext[1], previous_xlim[1])] + if ax.yaxis_inverted(): # case for e.g. images + y_ext = [min(y_ext[0], previous_ylim[1]), max(y_ext[1], previous_ylim[0])] + else: # case for e.g. labels + y_ext = [min(y_ext[0], previous_ylim[0]), max(y_ext[1], previous_ylim[1])] + # round because we need integers + plot_width = int(np.round(x_ext[1] - x_ext[0])) + plot_height = int(np.round(y_ext[1] - y_ext[0])) + + # use datashader for the visualization of points + # TODO: what about trans/norm at this point? + cvs = ds.Canvas(plot_width=plot_width, plot_height=plot_height, x_range=x_ext, y_range=y_ext) + + color_by_categorical = col_for_color is not None and points[col_for_color].values.dtype == object + aggregate_with_sum = None + if col_for_color is not None and (render_params.groups is None or len(render_params.groups) > 1): + if color_by_categorical: + agg = cvs.points(sdata_filt.points[element], "x", "y", agg=ds.by(col_for_color, ds.count())) else: - agg = cvs.points(sdata_filt.points[e], "x", "y", agg=ds.count()) + # numerical + agg = cvs.points(sdata_filt.points[element], "x", "y", agg=ds.sum(column=col_for_color)) + # save min and max values for drawing the colorbar + aggregate_with_sum = (agg.min(), agg.max()) + else: + agg = cvs.points(sdata_filt.points[element], "x", "y", agg=ds.count()) - color_key = ( - [x[:-2] for x in color_vector.categories.values] - if (type(color_vector) == pd.core.arrays.categorical.Categorical) - and (len(color_vector.categories.values) > 1) - else None + color_key = ( + [x[:-2] for x in color_vector.categories.values] + if (type(color_vector) == pd.core.arrays.categorical.Categorical) + and (len(color_vector.categories.values) > 1) + else None + ) + if color_by_categorical or col_for_color is None: + ds_result = ds.tf.shade( + ds.tf.spread(agg, px=px), + rescale_discrete_levels=True, + cmap=color_vector[0][:-2], + color_key=color_key, + min_alpha=np.min([150, render_params.alpha * 255]), + ) # TODO: choose other value than 150 for min_alpha (here and below)? + else: + ds_result = ds.tf.shade( + ds.tf.spread(agg, px=px), + rescale_discrete_levels=True, + cmap=render_params.cmap_params.cmap, + # color_key=color_key, + ) + # render image + rbga_image = np.transpose(ds_result.to_numpy().base, (0, 1, 2)) + cax = ax.imshow(rbga_image, zorder=render_params.zorder, alpha=render_params.alpha) + if aggregate_with_sum is not None: + cax = ScalarMappable( + norm=matplotlib.colors.Normalize(vmin=aggregate_with_sum[0], vmax=aggregate_with_sum[1]), + cmap=render_params.cmap_params.cmap, ) - if color_by_categorical or col_for_color is None: - ds_result = ds.tf.shade( - ds.tf.spread(agg, px=px), - rescale_discrete_levels=True, - cmap=color_vector[0][:-2], - color_key=color_key, - min_alpha=np.min([150, render_params.alpha * 255]), - ) # TODO: choose other value than 150 for min_alpha (here and below)? - else: - ds_result = ds.tf.shade( - ds.tf.spread(agg, px=px), - rescale_discrete_levels=True, - cmap=render_params.cmap_params.cmap, - # color_key=color_key, - ) - # render image - rbga_image = np.transpose(ds_result.to_numpy().base, (0, 1, 2)) - cax = ax.imshow(rbga_image, zorder=render_params.zorder, alpha=render_params.alpha) - if aggregate_with_sum is not None: - cax = ScalarMappable( - norm=matplotlib.colors.Normalize(vmin=aggregate_with_sum[0], vmax=aggregate_with_sum[1]), - cmap=render_params.cmap_params.cmap, - ) elif method == "matplotlib": # update axis limits if plot was empty before (necessary if datashader comes after) @@ -532,7 +513,7 @@ def _render_points( cax = ax.add_collection(_cax) if update_parameters: # necessary if points are plotted with mpl first and then with datashader - extent = get_extent(sdata_filt.points[e], coordinate_system=coordinate_system) + extent = get_extent(sdata_filt.points[element], coordinate_system=coordinate_system) ax.set_xbound(extent["x"]) ax.set_ybound(extent["y"]) @@ -547,7 +528,7 @@ def _render_points( cax=cax, fig_params=fig_params, adata=adata, - value_to_plot=render_params.col_for_color, + value_to_plot=col_for_color, color_source_vector=color_source_vector, palette=palette, alpha=render_params.alpha, @@ -648,65 +629,39 @@ def _render_images( # render_params.cmap_params.cmap if render_params.palette is None else _get_linear_colormap(palette, "k")[0] ) - im = ax.imshow( - layer, - cmap=cmap, - zorder=render_params.zorder, - ) - im.set_transform(trans_data) + # Overwrite alpha in cmap: https://stackoverflow.com/a/10127675 + cmap._init() + cmap._lut[:, -1] = render_params.alpha - # 2) Image has any number of channels but 1 - else: - layers = {} - for i, c in enumerate(channels): - layers[c] = img.sel(c=c).copy(deep=True).squeeze() - - if render_params.quantiles_for_norm != (None, None): - layers[c] = _normalize( - layers[c], - pmin=render_params.quantiles_for_norm[0], - pmax=render_params.quantiles_for_norm[1], - clip=True, - ) - - if not isinstance(render_params.cmap_params, list): - if render_params.cmap_params.norm is not None: - layers[c] = render_params.cmap_params.norm(layers[c]) - else: - if render_params.cmap_params[i].norm is not None: - layers[c] = render_params.cmap_params[i].norm(layers[c]) - - # 2A) Image has 3 channels, no palette info, and no/only one cmap was given - if n_channels == 3 and render_params.palette is None and not isinstance(render_params.cmap_params, list): - if render_params.cmap_params.is_default: # -> use RGB - stacked = np.stack([layers[c] for c in channels], axis=-1) - else: # -> use given cmap for each channel - channel_cmaps = [render_params.cmap_params.cmap] * n_channels - # Apply cmaps to each channel, add up and normalize to [0, 1] - stacked = ( - np.stack([channel_cmaps[i](layers[c]) for i, c in enumerate(channels)], 0).sum(0) / n_channels - ) - # Remove alpha channel so we can overwrite it from render_params.alpha - stacked = stacked[:, :, :3] - logger.warning( - "One cmap was given for multiple channels and is now used for each channel. " - "You're blending multiple cmaps. " - "If the plot doesn't look like you expect, it might be because your " - "cmaps go from a given color to 'white', and not to 'transparent'. " - "Therefore, the 'white' of higher layers will overlay the lower layers. " - "Consider using 'palette' instead." - ) - - im = ax.imshow( - stacked, - alpha=render_params.alpha, - zorder=render_params.zorder, + _ax_show_and_transform(layer, trans_data, ax, cmap=cmap) + # im = ax.imshow( + # layer, + # cmap=cmap, + # zorder=render_params.zorder, + # ) + # im.set_transform(trans_data) + + # 2) Image has any number of channels but 1 + else: + layers = {} + for ch_index, c in enumerate(channels): + layers[c] = img.sel(c=c).copy(deep=True).squeeze() + + if render_params.percentiles_for_norm != (None, None): + layers[c] = _normalize( + layers[c], + pmin=render_params.percentiles_for_norm[0], + pmax=render_params.percentiles_for_norm[1], + clip=True, ) - if not isinstance(render_params.cmap_params, list) and render_params.cmap_params.norm: - layers[c] = render_params.cmap_params.norm(layers[c]) - elif isinstance(render_params.cmap_params, list) and render_params.cmap_params[ch_index].norm: - layers[c] = render_params.cmap_params[ch_index].norm(layers[c]) + # TODO: can't be list anymore??? + if not isinstance(render_params.cmap_params, list): + if render_params.cmap_params.norm is not None: + layers[c] = render_params.cmap_params.norm(layers[c]) + else: + if render_params.cmap_params[ch_index].norm is not None: + layers[c] = render_params.cmap_params[ch_index].norm(layers[c]) # 2A) Image has 3 channels, no palette info, and no/only one cmap was given if palette is None and n_channels == 3 and not isinstance(render_params.cmap_params, list): @@ -729,7 +684,16 @@ def _render_images( "Consider using 'palette' instead." ) + # im = ax.imshow( + # stacked, + # alpha=render_params.alpha, + # zorder=render_params.zorder, + # ) _ax_show_and_transform(stacked, trans_data, ax, render_params.alpha) + # if not isinstance(render_params.cmap_params, list) and render_params.cmap_params.norm: + # layers[c] = render_params.cmap_params.norm(layers[c]) + # elif isinstance(render_params.cmap_params, list) and render_params.cmap_params[ch_index].norm: + # layers[c] = render_params.cmap_params[ch_index].norm(layers[c]) # 2B) Image has n channels, no palette/cmap info -> sample n categorical colors elif palette is None and not got_multiple_cmaps: @@ -741,12 +705,12 @@ def _render_images( channel_cmaps = [_get_linear_colormap([c], "k")[0] for c in seed_colors] - im = ax.imshow( - colored, - alpha=render_params.alpha, - zorder=render_params.zorder, - ) - im.set_transform(trans_data) + # im = ax.imshow( + # colored, + # alpha=render_params.alpha, + # zorder=render_params.zorder, + # ) + # im.set_transform(trans_data) # Apply cmaps to each channel and add up colored = np.stack([channel_cmaps[ind](layers[ch]) for ind, ch in enumerate(channels)], 0).sum(0) @@ -762,12 +726,12 @@ def _render_images( channel_cmaps = [_get_linear_colormap([c], "k")[0] for c in palette if isinstance(c, str)] - im = ax.imshow( - colored, - alpha=render_params.alpha, - zorder=render_params.zorder, - ) - im.set_transform(trans_data) + # im = ax.imshow( + # colored, + # alpha=render_params.alpha, + # zorder=render_params.zorder, + # ) + # im.set_transform(trans_data) # Apply cmaps to each channel and add up colored = np.stack([channel_cmaps[i](layers[c]) for i, c in enumerate(channels)], 0).sum(0) @@ -779,12 +743,12 @@ def _render_images( elif palette is None and got_multiple_cmaps: channel_cmaps = [cp.cmap for cp in render_params.cmap_params] # type: ignore[union-attr] - im = ax.imshow( - colored, - alpha=render_params.alpha, - zorder=render_params.zorder, - ) - im.set_transform(trans_data) + # im = ax.imshow( + # colored, + # alpha=render_params.alpha, + # zorder=render_params.zorder, + # ) + # im.set_transform(trans_data) # Apply cmaps to each channel, add up and normalize to [0, 1] colored = ( np.stack([channel_cmaps[ind](layers[ch]) for ind, ch in enumerate(channels)], 0).sum(0) / n_channels @@ -828,7 +792,8 @@ def _render_labels( extent = get_extent(label, coordinate_system=coordinate_system) # get best scale out of multiscale label - if isinstance(label, MultiscaleSpatialImage): + # if isinstance(label, MultiscaleSpatialImage): + if isinstance(label, DataTree): label = _multiscale_to_spatial_image( multiscale_image=label, dpi=fig_params.fig.dpi, @@ -837,6 +802,7 @@ def _render_labels( scale=scale, is_label=True, ) + # rasterize spatial image if necessary to speed up performance if rasterize: label = _rasterize_if_necessary( @@ -887,6 +853,17 @@ def _render_labels( seg_boundaries=render_params.outline, na_color=render_params.cmap_params.na_color, ) + _cax = ax.imshow( + labels_infill, + rasterized=True, + cmap=None if categorical else render_params.cmap_params.cmap, + norm=None if categorical else render_params.cmap_params.norm, + alpha=render_params.fill_alpha, + origin="lower", + zorder=render_params.zorder, + ) + _cax.set_transform(trans_data) + cax = ax.add_image(_cax) # Then overlay the contour labels_contour = _map_color_seg( @@ -899,7 +876,6 @@ def _render_labels( seg_boundaries=render_params.outline, na_color=render_params.cmap_params.na_color, ) - _cax = ax.imshow( labels_contour, rasterized=True, @@ -909,17 +885,7 @@ def _render_labels( origin="lower", zorder=render_params.zorder, ) - _cax = ax.imshow( - labels_infill, - rasterized=True, - cmap=None if categorical else render_params.cmap_params.cmap, - norm=None if categorical else render_params.cmap_params.norm, - alpha=render_params.fill_alpha, - origin="lower", - zorder=render_params.zorder, - ) - _cax.set_transform(trans_data) - cax = ax.add_image(_cax) + else: # Default: no alpha, contour = infill label = _map_color_seg( diff --git a/src/spatialdata_plot/pl/utils.py b/src/spatialdata_plot/pl/utils.py index 25246625..f958dbcf 100644 --- a/src/spatialdata_plot/pl/utils.py +++ b/src/spatialdata_plot/pl/utils.py @@ -46,6 +46,7 @@ from pandas.api.types import CategoricalDtype from scanpy import settings from scanpy.plotting._tools.scatterplots import _add_categorical_legend +from scanpy.plotting._utils import add_colors_for_categorical_sample_annotation from scanpy.plotting.palettes import default_20, default_28, default_102 from shapely.geometry import LineString, Polygon from skimage.color import label2rgb @@ -847,8 +848,6 @@ def _get_palette( def _maybe_set_colors( source: AnnData, target: AnnData, key: str, palette: str | ListedColormap | Cycler | Sequence[Any] | None = None ) -> None: - from scanpy.plotting._utils import add_colors_for_categorical_sample_annotation - color_key = f"{key}_colors" try: if palette is not None: From 0e4440378a359ceac7563fb5cf9a50ea47323997 Mon Sep 17 00:00:00 2001 From: Sonja Stockhaus Date: Mon, 24 Jun 2024 20:56:21 +0200 Subject: [PATCH 12/23] bugfixes --- src/spatialdata_plot/pl/basic.py | 2 +- src/spatialdata_plot/pl/render.py | 101 +++++++++++++++--------------- src/spatialdata_plot/pl/utils.py | 2 +- 3 files changed, 54 insertions(+), 51 deletions(-) diff --git a/src/spatialdata_plot/pl/basic.py b/src/spatialdata_plot/pl/basic.py index ca64d6cd..340d3c3f 100644 --- a/src/spatialdata_plot/pl/basic.py +++ b/src/spatialdata_plot/pl/basic.py @@ -540,7 +540,7 @@ def render_labels( element: str | None = None, color: str | None = None, groups: list[str] | str | None = None, - contour_px: int = 3, + contour_px: int | None = 3, outline: bool = False, palette: list[str] | str | None = None, cmap: Colormap | str | None = None, diff --git a/src/spatialdata_plot/pl/render.py b/src/spatialdata_plot/pl/render.py index 9564ad98..783ce402 100644 --- a/src/spatialdata_plot/pl/render.py +++ b/src/spatialdata_plot/pl/render.py @@ -358,6 +358,9 @@ def _render_points( if groups is not None and col_for_color is not None: points = points[points[col_for_color].isin(groups)] + # in case no rows are left: + if len(points) <= 0: + raise ValueError(f"None of the groups {groups} could be found in the column '{col_for_color}'.") # we construct an anndata to hack the plotting functions if table_name is None: @@ -493,55 +496,55 @@ def _render_points( cmap=render_params.cmap_params.cmap, ) - elif method == "matplotlib": - # update axis limits if plot was empty before (necessary if datashader comes after) - update_parameters = not _mpl_ax_contains_elements(ax) - # original way of plotting points - _cax = ax.scatter( - adata[:, 0].X.flatten(), - adata[:, 1].X.flatten(), - s=render_params.size, - c=color_vector, - rasterized=sc_settings._vector_friendly, - cmap=render_params.cmap_params.cmap, - norm=norm, - alpha=render_params.alpha, - transform=trans, - zorder=render_params.zorder, - # **kwargs, - ) - cax = ax.add_collection(_cax) - if update_parameters: - # necessary if points are plotted with mpl first and then with datashader - extent = get_extent(sdata_filt.points[element], coordinate_system=coordinate_system) - ax.set_xbound(extent["x"]) - ax.set_ybound(extent["y"]) - - if len(set(color_vector)) != 1 or list(set(color_vector))[0] != to_hex(render_params.cmap_params.na_color): - if color_source_vector is None: - palette = ListedColormap(dict.fromkeys(color_vector)) - else: - palette = ListedColormap(dict.fromkeys(color_vector[~pd.Categorical(color_source_vector).isnull()])) - - _ = _decorate_axs( - ax=ax, - cax=cax, - fig_params=fig_params, - adata=adata, - value_to_plot=col_for_color, - color_source_vector=color_source_vector, - palette=palette, - alpha=render_params.alpha, - na_color=render_params.cmap_params.na_color, - legend_fontsize=legend_params.legend_fontsize, - legend_fontweight=legend_params.legend_fontweight, - legend_loc=legend_params.legend_loc, - legend_fontoutline=legend_params.legend_fontoutline, - na_in_legend=legend_params.na_in_legend, - colorbar=legend_params.colorbar, - scalebar_dx=scalebar_params.scalebar_dx, - scalebar_units=scalebar_params.scalebar_units, - ) + elif method == "matplotlib": + # update axis limits if plot was empty before (necessary if datashader comes after) + update_parameters = not _mpl_ax_contains_elements(ax) + # original way of plotting points + _cax = ax.scatter( + adata[:, 0].X.flatten(), + adata[:, 1].X.flatten(), + s=render_params.size, + c=color_vector, + rasterized=sc_settings._vector_friendly, + cmap=render_params.cmap_params.cmap, + norm=norm, + alpha=render_params.alpha, + transform=trans, + zorder=render_params.zorder, + # **kwargs, + ) + cax = ax.add_collection(_cax) + if update_parameters: + # necessary if points are plotted with mpl first and then with datashader + extent = get_extent(sdata_filt.points[element], coordinate_system=coordinate_system) + ax.set_xbound(extent["x"]) + ax.set_ybound(extent["y"]) + + if len(set(color_vector)) != 1 or list(set(color_vector))[0] != to_hex(render_params.cmap_params.na_color): + if color_source_vector is None: + palette = ListedColormap(dict.fromkeys(color_vector)) + else: + palette = ListedColormap(dict.fromkeys(color_vector[~pd.Categorical(color_source_vector).isnull()])) + + _ = _decorate_axs( + ax=ax, + cax=cax, + fig_params=fig_params, + adata=adata, + value_to_plot=col_for_color, + color_source_vector=color_source_vector, + palette=palette, + alpha=render_params.alpha, + na_color=render_params.cmap_params.na_color, + legend_fontsize=legend_params.legend_fontsize, + legend_fontweight=legend_params.legend_fontweight, + legend_loc=legend_params.legend_loc, + legend_fontoutline=legend_params.legend_fontoutline, + na_in_legend=legend_params.na_in_legend, + colorbar=legend_params.colorbar, + scalebar_dx=scalebar_params.scalebar_dx, + scalebar_units=scalebar_params.scalebar_units, + ) def _render_images( diff --git a/src/spatialdata_plot/pl/utils.py b/src/spatialdata_plot/pl/utils.py index f958dbcf..865c82a6 100644 --- a/src/spatialdata_plot/pl/utils.py +++ b/src/spatialdata_plot/pl/utils.py @@ -1646,7 +1646,7 @@ def _validate_label_render_params( cmap: list[Colormap | str] | Colormap | str | None, color: str | None, fill_alpha: float | int, - contour_px: int, + contour_px: int | None, outline: bool, groups: list[str] | str | None, palette: list[str] | str | None, From 75f9a3def5840f7c3e73905a3c818efa624ff148 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Tue, 25 Jun 2024 12:17:48 +0000 Subject: [PATCH 13/23] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- src/spatialdata_plot/pl/basic.py | 7 +++---- src/spatialdata_plot/pl/render.py | 5 +++-- 2 files changed, 6 insertions(+), 6 deletions(-) diff --git a/src/spatialdata_plot/pl/basic.py b/src/spatialdata_plot/pl/basic.py index da06fb6f..6668f3bb 100644 --- a/src/spatialdata_plot/pl/basic.py +++ b/src/spatialdata_plot/pl/basic.py @@ -19,11 +19,10 @@ from matplotlib.axes import Axes from matplotlib.colors import Colormap, Normalize from matplotlib.figure import Figure -from multiscale_spatial_image.multiscale_spatial_image import MultiscaleSpatialImage -from spatial_image import SpatialImage -#from spatialdata._core.data_extent import get_extent -from spatialdata._utils import _deprecation_alias from spatialdata import get_extent + +# from spatialdata._core.data_extent import get_extent +from spatialdata._utils import _deprecation_alias from xarray import DataArray from spatialdata_plot._accessor import register_spatial_data_accessor diff --git a/src/spatialdata_plot/pl/render.py b/src/spatialdata_plot/pl/render.py index 446f112f..c37f3ff9 100644 --- a/src/spatialdata_plot/pl/render.py +++ b/src/spatialdata_plot/pl/render.py @@ -15,9 +15,10 @@ import scanpy as sc import spatialdata as sd from anndata import AnnData -#from datatree.datatree import DataTree -from matplotlib.cm import ScalarMappable from datatree import DataTree + +# from datatree.datatree import DataTree +from matplotlib.cm import ScalarMappable from matplotlib.colors import ListedColormap, Normalize from scanpy._settings import settings as sc_settings from spatialdata import get_extent From d1e2709bf4f17870c4230bce65aef27018794bdb Mon Sep 17 00:00:00 2001 From: Sonja Stockhaus Date: Wed, 3 Jul 2024 17:58:40 +0200 Subject: [PATCH 14/23] bugfix point transformation info lost --- src/spatialdata_plot/pl/render.py | 61 ++++++++++++++++++++++++++----- 1 file changed, 52 insertions(+), 9 deletions(-) diff --git a/src/spatialdata_plot/pl/render.py b/src/spatialdata_plot/pl/render.py index c37f3ff9..04171a92 100644 --- a/src/spatialdata_plot/pl/render.py +++ b/src/spatialdata_plot/pl/render.py @@ -25,6 +25,7 @@ from spatialdata.models import PointsModel, get_table_keys from spatialdata.transformations import ( get_transformation, + set_transformation, ) from spatialdata_plot._logging import logger @@ -254,15 +255,34 @@ def _render_shapes( scale = shapes[is_point]["radius"] * render_params.scale shapes.loc[is_point, "geometry"] = _geometry[is_point].buffer(scale) - agg = cvs.polygons(shapes, geometry="geometry", agg=ds.count()) + # if isinstance(shapes, geopandas.geodataframe.GeoDataFrame): + # shapes = spatialpandas.GeoDataFrame(shapes) + + # in case we are coloring by a column in table + if col_for_color is not None and col_for_color not in sdata_filt.shapes[element].columns: + # numerical + if color_source_vector is None: + sdata_filt.shapes[element][col_for_color] = color_vector + else: # categorical + sdata_filt.shapes[element][col_for_color] = color_source_vector # Render shapes with datashader - if render_params.col_for_color is not None and ( - render_params.groups is None or len(render_params.groups) > 1 - ): - agg = cvs.polygons(shapes, geometry="geometry", agg=ds.by(render_params.col_for_color, ds.count())) + color_by_categorical = col_for_color is not None and color_source_vector is not None + aggregate_with_sum = None + if col_for_color is not None and (render_params.groups is None or len(render_params.groups) > 1): + if color_by_categorical: + agg = cvs.polygons( + sdata_filt.shapes[element], geometry="geometry", agg=ds.by(col_for_color, ds.count()) + ) + else: + # numerical + agg = cvs.polygons( + sdata_filt.shapes[element], geometry="geometry", agg=ds.sum(column=col_for_color) + ) + # save min and max values for drawing the colorbar + aggregate_with_sum = (agg.min(), agg.max()) else: - agg = cvs.polygons(shapes, geometry="geometry", agg=ds.count()) + agg = cvs.polygons(sdata_filt.shapes[element], geometry="geometry", agg=ds.count()) color_key = ( [x[:-2] for x in color_vector.categories.values] @@ -270,15 +290,33 @@ def _render_shapes( and (len(color_vector.categories.values) > 1) else None ) - ds_result = ds.tf.shade( - agg, cmap=color_vector[0][:-2], alpha=render_params.fill_alpha * 255, color_key=color_key, min_alpha=200 - ) + # ds_result = ds.tf.shade( + # agg, cmap=color_vector[0][:-2], alpha=render_params.fill_alpha * 255, color_key=color_key, + # min_alpha=200 + # ) + if color_by_categorical or col_for_color is None: + ds_result = ds.tf.shade( + agg, + cmap=color_vector[0][:-2], + color_key=color_key, + min_alpha=np.min([150, render_params.fill_alpha * 255]), + ) # TODO: choose other value than 150 for min_alpha (here and below)? + else: + ds_result = ds.tf.shade( + agg, + cmap=render_params.cmap_params.cmap, + ) # Render image rgba_image = np.transpose(ds_result.to_numpy().base, (0, 1, 2)) _cax = ax.imshow(rgba_image, cmap=palette, zorder=render_params.zorder) _cax.set_transform(trans) cax = ax.add_image(_cax) + if aggregate_with_sum is not None: + cax = ScalarMappable( + norm=matplotlib.colors.Normalize(vmin=aggregate_with_sum[0], vmax=aggregate_with_sum[1]), + cmap=render_params.cmap_params.cmap, + ) # Sets the limits of the colorbar to the values instead of [0, 1] if not norm and not values_are_categorical: @@ -380,8 +418,13 @@ def _render_points( # we can do this because of dealing with a copy # Convert back to dask dataframe to modify sdata + transformation_in_cs = sdata_filt.points[element].attrs["transform"][coordinate_system] points = dask.dataframe.from_pandas(points, npartitions=1) sdata_filt.points[element] = PointsModel.parse(points, coordinates={"x": "x", "y": "y"}) + # restore transformation in coordinate system of interest + set_transformation( + element=sdata_filt.points[element], transformation=transformation_in_cs, to_coordinate_system=coordinate_system + ) if col_for_color is not None: cols = sc.get.obs_df(adata, col_for_color) From 6a61693a6d9fc0eb0b922da0fb5050f0b3263013 Mon Sep 17 00:00:00 2001 From: Sonja Stockhaus Date: Tue, 9 Jul 2024 21:09:49 +0200 Subject: [PATCH 15/23] render_shapes cleanup --- src/spatialdata_plot/pl/render.py | 369 ++++++++++++++---------------- 1 file changed, 178 insertions(+), 191 deletions(-) diff --git a/src/spatialdata_plot/pl/render.py b/src/spatialdata_plot/pl/render.py index 04171a92..36c125eb 100644 --- a/src/spatialdata_plot/pl/render.py +++ b/src/spatialdata_plot/pl/render.py @@ -76,7 +76,6 @@ def _render_shapes( filter_tables=bool(render_params.table_name), ) - # for index, e in enumerate(elements): shapes = sdata[element] if (table_name := render_params.table_name) is None: @@ -125,7 +124,6 @@ def _render_shapes( color_vector = [render_params.cmap_params.na_color] # filter by `groups` - if isinstance(groups, list) and color_source_vector is not None: mask = color_source_vector.isin(groups) shapes = shapes[mask] @@ -134,50 +132,51 @@ def _render_shapes( color_vector = color_vector[mask] shapes = gpd.GeoDataFrame(shapes, geometry="geometry") - _cax = _get_collection_shape( - shapes=shapes, - s=render_params.scale, - c=color_vector, - render_params=render_params, - rasterized=sc_settings._vector_friendly, - cmap=render_params.cmap_params.cmap, - norm=norm, - fill_alpha=render_params.fill_alpha, - outline_alpha=render_params.outline_alpha, - # **kwargs, - ) - - # Sets the limits of the colorbar to the values instead of [0, 1] - if not norm and not values_are_categorical: - _cax.set_clim(min(color_vector), max(color_vector)) - - values_are_categorical = color_source_vector is not None - - # color_source_vector is None when the values aren't categorical - if values_are_categorical and render_params.transfunc is not None: - color_vector = render_params.transfunc(color_vector) - - norm = copy(render_params.cmap_params.norm) - - if len(color_vector) == 0: - color_vector = [render_params.cmap_params.na_color] - - # filter by `groups` - if render_params.groups is not None and color_source_vector is not None: - mask = color_source_vector.isin(render_params.groups) - shapes = shapes[mask] - shapes = shapes.reset_index() - color_source_vector = color_source_vector[mask] - color_vector = color_vector[mask] - cax = ax.add_collection(_cax) + # _cax = _get_collection_shape( + # shapes=shapes, + # s=render_params.scale, + # c=color_vector, + # render_params=render_params, + # rasterized=sc_settings._vector_friendly, + # cmap=render_params.cmap_params.cmap, + # norm=norm, + # fill_alpha=render_params.fill_alpha, + # outline_alpha=render_params.outline_alpha, + # zorder=render_params.zorder, + # # **kwargs, + # ) + + # Sets the limits of the colorbar to the right values instead of [0, 1] + # if not norm and not values_are_categorical: + # _cax.set_clim(min(color_vector), max(color_vector)) + + # values_are_categorical = color_source_vector is not None + + # # color_source_vector is None when the values aren't categorical + # if values_are_categorical and render_params.transfunc is not None: + # color_vector = render_params.transfunc(color_vector) + + # norm = copy(render_params.cmap_params.norm) + + # if len(color_vector) == 0: + # color_vector = [render_params.cmap_params.na_color] + + # # filter by `groups` + # if render_params.groups is not None and color_source_vector is not None: + # mask = color_source_vector.isin(render_params.groups) + # shapes = shapes[mask] + # shapes = shapes.reset_index() + # color_source_vector = color_source_vector[mask] + # color_vector = color_vector[mask] + # cax = ax.add_collection(_cax) # Apply the transformation to the PatchCollection's paths - trans = get_transformation(sdata_filt[element], get_all=True)[coordinate_system] - affine_trans = trans.to_affine_matrix(input_axes=("x", "y"), output_axes=("x", "y")) - trans = mtransforms.Affine2D(matrix=affine_trans) + # trans = get_transformation(sdata_filt[element], get_all=True)[coordinate_system] + # affine_trans = trans.to_affine_matrix(input_axes=("x", "y"), output_axes=("x", "y")) + # trans = mtransforms.Affine2D(matrix=affine_trans) - for path in _cax.get_paths(): - path.vertices = trans.transform(path.vertices) + # for path in _cax.get_paths(): + # path.vertices = trans.transform(path.vertices) # Using dict.fromkeys here since set returns in arbitrary order # remove the color of NaN values, else it might be assigned to a category @@ -195,163 +194,151 @@ def _render_shapes( # False if user specified color-like with 'color' parameter colorbar = False if col_for_color is None else legend_params.colorbar - # Apply the transformation to the PatchCollection's paths - trans = get_transformation(sdata_filt.shapes[element], get_all=True)[coordinate_system] - affine_trans = trans.to_affine_matrix(input_axes=("x", "y"), output_axes=("x", "y")) - trans = mtransforms.Affine2D(matrix=affine_trans) - - shapes = gpd.GeoDataFrame(shapes, geometry="geometry") - - # Determine which method to use for rendering - method = render_params.method - if method is None: - method = "datashader" if len(shapes) > 10000 else "matplotlib" - elif method not in ["matplotlib", "datashader"]: - raise ValueError("Method must be either 'matplotlib' or 'datashader'.") - logger.info(f"Using {method}") - - if method == "matplotlib": - _cax = _get_collection_shape( - shapes=shapes, - s=render_params.scale, - c=color_vector, - render_params=render_params, - rasterized=sc_settings._vector_friendly, - cmap=render_params.cmap_params.cmap, - norm=norm, - fill_alpha=render_params.fill_alpha, - outline_alpha=render_params.outline_alpha, - zorder=render_params.zorder, - # **kwargs, - ) + # Apply the transformation to the PatchCollection's paths + trans = get_transformation(sdata_filt.shapes[element], get_all=True)[coordinate_system] + affine_trans = trans.to_affine_matrix(input_axes=("x", "y"), output_axes=("x", "y")) + trans = mtransforms.Affine2D(matrix=affine_trans) + + shapes = gpd.GeoDataFrame(shapes, geometry="geometry") + + # Determine which method to use for rendering + method = render_params.method + if method is None: + method = "datashader" if len(shapes) > 10000 else "matplotlib" + elif method not in ["matplotlib", "datashader"]: + raise ValueError("Method must be either 'matplotlib' or 'datashader'.") + logger.info(f"Using {method}") + + if method == "matplotlib": + _cax = _get_collection_shape( + shapes=shapes, + s=render_params.scale, + c=color_vector, + render_params=render_params, + rasterized=sc_settings._vector_friendly, + cmap=render_params.cmap_params.cmap, + norm=norm, + fill_alpha=render_params.fill_alpha, + outline_alpha=render_params.outline_alpha, + zorder=render_params.zorder, + # **kwargs, + ) + cax = ax.add_collection(_cax) + + # Transform the paths in PatchCollection + for path in _cax.get_paths(): + path.vertices = trans.transform(path.vertices) cax = ax.add_collection(_cax) + elif method == "datashader": + # Where to put this + trans = mtransforms.Affine2D(matrix=affine_trans) + ax.transData + + extent = get_extent(sdata.shapes[element]) + x_ext = extent["x"][1] + y_ext = extent["y"][1] + x_range = [0, x_ext] + y_range = [0, y_ext] + # round because we need integers + plot_width = int(np.round(x_range[1] - x_range[0])) + plot_height = int(np.round(y_range[1] - y_range[0])) - # Transform the paths in PatchCollection - for path in _cax.get_paths(): - path.vertices = trans.transform(path.vertices) - cax = ax.add_collection(_cax) - elif method == "datashader": - # Where to put this - trans = mtransforms.Affine2D(matrix=affine_trans) + ax.transData - - extent = get_extent(sdata.shapes[element]) - x_ext = extent["x"][1] - y_ext = extent["y"][1] - # previous_xlim = fig_params.ax.get_xlim() - # previous_ylim = fig_params.ax.get_ylim() - x_range = [0, x_ext] - y_range = [0, y_ext] - # round because we need integers - plot_width = int(np.round(x_range[1] - x_range[0])) - plot_height = int(np.round(y_range[1] - y_range[0])) - - cvs = ds.Canvas(plot_width=plot_width, plot_height=plot_height, x_range=x_range, y_range=y_range) - - _geometry = shapes["geometry"] - is_point = _geometry.type == "Point" - - # Handle circles encoded as points with radius - if is_point.any(): # TODO - scale = shapes[is_point]["radius"] * render_params.scale - shapes.loc[is_point, "geometry"] = _geometry[is_point].buffer(scale) - - # if isinstance(shapes, geopandas.geodataframe.GeoDataFrame): - # shapes = spatialpandas.GeoDataFrame(shapes) - - # in case we are coloring by a column in table - if col_for_color is not None and col_for_color not in sdata_filt.shapes[element].columns: - # numerical - if color_source_vector is None: - sdata_filt.shapes[element][col_for_color] = color_vector - else: # categorical - sdata_filt.shapes[element][col_for_color] = color_source_vector - - # Render shapes with datashader - color_by_categorical = col_for_color is not None and color_source_vector is not None - aggregate_with_sum = None - if col_for_color is not None and (render_params.groups is None or len(render_params.groups) > 1): - if color_by_categorical: - agg = cvs.polygons( - sdata_filt.shapes[element], geometry="geometry", agg=ds.by(col_for_color, ds.count()) - ) - else: - # numerical - agg = cvs.polygons( - sdata_filt.shapes[element], geometry="geometry", agg=ds.sum(column=col_for_color) - ) - # save min and max values for drawing the colorbar - aggregate_with_sum = (agg.min(), agg.max()) - else: - agg = cvs.polygons(sdata_filt.shapes[element], geometry="geometry", agg=ds.count()) + cvs = ds.Canvas(plot_width=plot_width, plot_height=plot_height, x_range=x_range, y_range=y_range) - color_key = ( - [x[:-2] for x in color_vector.categories.values] - if (type(color_vector) == pd.core.arrays.categorical.Categorical) - and (len(color_vector.categories.values) > 1) - else None - ) - # ds_result = ds.tf.shade( - # agg, cmap=color_vector[0][:-2], alpha=render_params.fill_alpha * 255, color_key=color_key, - # min_alpha=200 - # ) - if color_by_categorical or col_for_color is None: - ds_result = ds.tf.shade( - agg, - cmap=color_vector[0][:-2], - color_key=color_key, - min_alpha=np.min([150, render_params.fill_alpha * 255]), - ) # TODO: choose other value than 150 for min_alpha (here and below)? - else: - ds_result = ds.tf.shade( - agg, - cmap=render_params.cmap_params.cmap, - ) + _geometry = shapes["geometry"] + is_point = _geometry.type == "Point" + + # Handle circles encoded as points with radius + if is_point.any(): # TODO + scale = shapes[is_point]["radius"] * render_params.scale + sdata_filt.shapes[element].loc[is_point, "geometry"] = _geometry[is_point].buffer(scale) - # Render image - rgba_image = np.transpose(ds_result.to_numpy().base, (0, 1, 2)) - _cax = ax.imshow(rgba_image, cmap=palette, zorder=render_params.zorder) - _cax.set_transform(trans) - cax = ax.add_image(_cax) - if aggregate_with_sum is not None: - cax = ScalarMappable( - norm=matplotlib.colors.Normalize(vmin=aggregate_with_sum[0], vmax=aggregate_with_sum[1]), - cmap=render_params.cmap_params.cmap, + # in case we are coloring by a column in table + if col_for_color is not None and col_for_color not in sdata_filt.shapes[element].columns: + # numerical + if color_source_vector is None: + sdata_filt.shapes[element][col_for_color] = color_vector + else: # categorical + sdata_filt.shapes[element][col_for_color] = color_source_vector + + # Render shapes with datashader + color_by_categorical = col_for_color is not None and color_source_vector is not None + aggregate_with_sum = None + if col_for_color is not None and (render_params.groups is None or len(render_params.groups) > 1): + if color_by_categorical: + agg = cvs.polygons( + sdata_filt.shapes[element], geometry="geometry", agg=ds.by(col_for_color, ds.count()) ) + else: + # numerical + agg = cvs.polygons(sdata_filt.shapes[element], geometry="geometry", agg=ds.sum(column=col_for_color)) + # save min and max values for drawing the colorbar + aggregate_with_sum = (agg.min(), agg.max()) + else: + agg = cvs.polygons(sdata_filt.shapes[element], geometry="geometry", agg=ds.count()) - # Sets the limits of the colorbar to the values instead of [0, 1] - if not norm and not values_are_categorical: - _cax.set_clim(min(color_vector), max(color_vector)) + color_key = ( + [x[:-2] for x in color_vector.categories.values] + if (type(color_vector) == pd.core.arrays.categorical.Categorical) + and (len(color_vector.categories.values) > 1) + else None + ) - if not ( - len(set(color_vector)) == 1 and list(set(color_vector))[0] == to_hex(render_params.cmap_params.na_color) - ): - # necessary in case different shapes elements are annotated with one table - if color_source_vector is not None and render_params.col_for_color is not None: - color_source_vector = color_source_vector.remove_unused_categories() - - # False if user specified color-like with 'color' parameter - colorbar = False if render_params.col_for_color is None else legend_params.colorbar - - _ = _decorate_axs( - ax=ax, - cax=cax, - fig_params=fig_params, - adata=table, - value_to_plot=col_for_color, - color_source_vector=color_source_vector, - palette=palette, - alpha=render_params.fill_alpha, - na_color=render_params.cmap_params.na_color, - legend_fontsize=legend_params.legend_fontsize, - legend_fontweight=legend_params.legend_fontweight, - legend_loc=legend_params.legend_loc, - legend_fontoutline=legend_params.legend_fontoutline, - na_in_legend=legend_params.na_in_legend, - colorbar=colorbar, - scalebar_dx=scalebar_params.scalebar_dx, - scalebar_units=scalebar_params.scalebar_units, + if color_by_categorical or col_for_color is None: + ds_result = ds.tf.shade( + agg, + cmap=color_vector[0][:-2], + color_key=color_key, + min_alpha=np.min([150, render_params.fill_alpha * 255]), + ) # TODO: choose other value than 150 for min_alpha (here and below)? + else: + ds_result = ds.tf.shade( + agg, + cmap=render_params.cmap_params.cmap, ) + # Render image + rgba_image = np.transpose(ds_result.to_numpy().base, (0, 1, 2)) + _cax = ax.imshow(rgba_image, cmap=palette, zorder=render_params.zorder) + _cax.set_transform(trans) + cax = ax.add_image(_cax) + if aggregate_with_sum is not None: + cax = ScalarMappable( + norm=matplotlib.colors.Normalize(vmin=aggregate_with_sum[0], vmax=aggregate_with_sum[1]), + cmap=render_params.cmap_params.cmap, + ) + + # Sets the limits of the colorbar to the values instead of [0, 1] + if not norm and not values_are_categorical: + _cax.set_clim(min(color_vector), max(color_vector)) + + if not (len(set(color_vector)) == 1 and list(set(color_vector))[0] == to_hex(render_params.cmap_params.na_color)): + # necessary in case different shapes elements are annotated with one table + if color_source_vector is not None and render_params.col_for_color is not None: + color_source_vector = color_source_vector.remove_unused_categories() + + # False if user specified color-like with 'color' parameter + colorbar = False if render_params.col_for_color is None else legend_params.colorbar + + _ = _decorate_axs( + ax=ax, + cax=cax, + fig_params=fig_params, + adata=table, + value_to_plot=col_for_color, + color_source_vector=color_source_vector, + palette=palette, + alpha=render_params.fill_alpha, + na_color=render_params.cmap_params.na_color, + legend_fontsize=legend_params.legend_fontsize, + legend_fontweight=legend_params.legend_fontweight, + legend_loc=legend_params.legend_loc, + legend_fontoutline=legend_params.legend_fontoutline, + na_in_legend=legend_params.na_in_legend, + colorbar=colorbar, + scalebar_dx=scalebar_params.scalebar_dx, + scalebar_units=scalebar_params.scalebar_units, + ) + def _render_points( sdata: sd.SpatialData, From e20811e1a862540e7f8d0f94780f5be8038d9dd6 Mon Sep 17 00:00:00 2001 From: Sonja Stockhaus Date: Tue, 9 Jul 2024 22:00:23 +0200 Subject: [PATCH 16/23] zorder bugfix and added tests for zorder and datashader stuff --- src/spatialdata_plot/pl/render.py | 96 ++---------------- src/spatialdata_plot/pl/render_params.py | 8 +- src/spatialdata_plot/pl/utils.py | 3 + tests/_images/Images_can_stick_to_zorder.png | Bin 0 -> 29437 bytes ...oints_datashader_can_color_by_category.png | Bin 0 -> 7509 bytes .../Points_datashader_continuous_color.png | Bin 0 -> 14453 bytes .../Points_datashader_matplotlib_stack.png | Bin 0 -> 14365 bytes ...hapes_datashader_can_color_by_category.png | Bin 0 -> 6103 bytes .../Shapes_datashader_can_color_by_value.png | Bin 0 -> 6815 bytes .../Shapes_datashader_can_render_shapes.png | Bin 0 -> 5776 bytes tests/pl/test_render_images.py | 3 + tests/pl/test_render_points.py | 15 +++ tests/pl/test_render_shapes.py | 21 ++++ 13 files changed, 53 insertions(+), 93 deletions(-) create mode 100644 tests/_images/Images_can_stick_to_zorder.png create mode 100644 tests/_images/Points_datashader_can_color_by_category.png create mode 100644 tests/_images/Points_datashader_continuous_color.png create mode 100644 tests/_images/Points_datashader_matplotlib_stack.png create mode 100644 tests/_images/Shapes_datashader_can_color_by_category.png create mode 100644 tests/_images/Shapes_datashader_can_color_by_value.png create mode 100644 tests/_images/Shapes_datashader_can_render_shapes.png diff --git a/src/spatialdata_plot/pl/render.py b/src/spatialdata_plot/pl/render.py index 36c125eb..b77cc40e 100644 --- a/src/spatialdata_plot/pl/render.py +++ b/src/spatialdata_plot/pl/render.py @@ -132,52 +132,6 @@ def _render_shapes( color_vector = color_vector[mask] shapes = gpd.GeoDataFrame(shapes, geometry="geometry") - # _cax = _get_collection_shape( - # shapes=shapes, - # s=render_params.scale, - # c=color_vector, - # render_params=render_params, - # rasterized=sc_settings._vector_friendly, - # cmap=render_params.cmap_params.cmap, - # norm=norm, - # fill_alpha=render_params.fill_alpha, - # outline_alpha=render_params.outline_alpha, - # zorder=render_params.zorder, - # # **kwargs, - # ) - - # Sets the limits of the colorbar to the right values instead of [0, 1] - # if not norm and not values_are_categorical: - # _cax.set_clim(min(color_vector), max(color_vector)) - - # values_are_categorical = color_source_vector is not None - - # # color_source_vector is None when the values aren't categorical - # if values_are_categorical and render_params.transfunc is not None: - # color_vector = render_params.transfunc(color_vector) - - # norm = copy(render_params.cmap_params.norm) - - # if len(color_vector) == 0: - # color_vector = [render_params.cmap_params.na_color] - - # # filter by `groups` - # if render_params.groups is not None and color_source_vector is not None: - # mask = color_source_vector.isin(render_params.groups) - # shapes = shapes[mask] - # shapes = shapes.reset_index() - # color_source_vector = color_source_vector[mask] - # color_vector = color_vector[mask] - # cax = ax.add_collection(_cax) - - # Apply the transformation to the PatchCollection's paths - # trans = get_transformation(sdata_filt[element], get_all=True)[coordinate_system] - # affine_trans = trans.to_affine_matrix(input_axes=("x", "y"), output_axes=("x", "y")) - # trans = mtransforms.Affine2D(matrix=affine_trans) - - # for path in _cax.get_paths(): - # path.vertices = trans.transform(path.vertices) - # Using dict.fromkeys here since set returns in arbitrary order # remove the color of NaN values, else it might be assigned to a category # order of color in the palette should agree to order of occurence @@ -230,7 +184,7 @@ def _render_shapes( path.vertices = trans.transform(path.vertices) cax = ax.add_collection(_cax) elif method == "datashader": - # Where to put this + # TODO: Where to put this trans = mtransforms.Affine2D(matrix=affine_trans) + ax.transData extent = get_extent(sdata.shapes[element]) @@ -248,7 +202,7 @@ def _render_shapes( is_point = _geometry.type == "Point" # Handle circles encoded as points with radius - if is_point.any(): # TODO + if is_point.any(): scale = shapes[is_point]["radius"] * render_params.scale sdata_filt.shapes[element].loc[is_point, "geometry"] = _geometry[is_point].buffer(scale) @@ -533,7 +487,6 @@ def _render_points( elif method == "matplotlib": # update axis limits if plot was empty before (necessary if datashader comes after) update_parameters = not _mpl_ax_contains_elements(ax) - # original way of plotting points _cax = ax.scatter( adata[:, 0].X.flatten(), adata[:, 1].X.flatten(), @@ -663,20 +616,13 @@ def _render_images( _get_linear_colormap(palette, "k")[0] if isinstance(palette, list) and all(isinstance(p, str) for p in palette) else render_params.cmap_params.cmap - # render_params.cmap_params.cmap if render_params.palette is None else _get_linear_colormap(palette, "k")[0] ) # Overwrite alpha in cmap: https://stackoverflow.com/a/10127675 cmap._init() cmap._lut[:, -1] = render_params.alpha - _ax_show_and_transform(layer, trans_data, ax, cmap=cmap) - # im = ax.imshow( - # layer, - # cmap=cmap, - # zorder=render_params.zorder, - # ) - # im.set_transform(trans_data) + _ax_show_and_transform(layer, trans_data, ax, cmap=cmap, zorder=render_params.zorder) # 2) Image has any number of channels but 1 else: @@ -721,16 +667,7 @@ def _render_images( "Consider using 'palette' instead." ) - # im = ax.imshow( - # stacked, - # alpha=render_params.alpha, - # zorder=render_params.zorder, - # ) - _ax_show_and_transform(stacked, trans_data, ax, render_params.alpha) - # if not isinstance(render_params.cmap_params, list) and render_params.cmap_params.norm: - # layers[c] = render_params.cmap_params.norm(layers[c]) - # elif isinstance(render_params.cmap_params, list) and render_params.cmap_params[ch_index].norm: - # layers[c] = render_params.cmap_params[ch_index].norm(layers[c]) + _ax_show_and_transform(stacked, trans_data, ax, render_params.alpha, zorder=render_params.zorder) # 2B) Image has n channels, no palette/cmap info -> sample n categorical colors elif palette is None and not got_multiple_cmaps: @@ -742,19 +679,13 @@ def _render_images( channel_cmaps = [_get_linear_colormap([c], "k")[0] for c in seed_colors] - # im = ax.imshow( - # colored, - # alpha=render_params.alpha, - # zorder=render_params.zorder, - # ) - # im.set_transform(trans_data) # Apply cmaps to each channel and add up colored = np.stack([channel_cmaps[ind](layers[ch]) for ind, ch in enumerate(channels)], 0).sum(0) # Remove alpha channel so we can overwrite it from render_params.alpha colored = colored[:, :, :3] - _ax_show_and_transform(colored, trans_data, ax, render_params.alpha) + _ax_show_and_transform(colored, trans_data, ax, render_params.alpha, zorder=render_params.zorder) # 2C) Image has n channels and palette info elif palette is not None and not got_multiple_cmaps: @@ -763,29 +694,17 @@ def _render_images( channel_cmaps = [_get_linear_colormap([c], "k")[0] for c in palette if isinstance(c, str)] - # im = ax.imshow( - # colored, - # alpha=render_params.alpha, - # zorder=render_params.zorder, - # ) - # im.set_transform(trans_data) # Apply cmaps to each channel and add up colored = np.stack([channel_cmaps[i](layers[c]) for i, c in enumerate(channels)], 0).sum(0) # Remove alpha channel so we can overwrite it from render_params.alpha colored = colored[:, :, :3] - _ax_show_and_transform(colored, trans_data, ax, render_params.alpha) + _ax_show_and_transform(colored, trans_data, ax, render_params.alpha, zorder=render_params.zorder) elif palette is None and got_multiple_cmaps: channel_cmaps = [cp.cmap for cp in render_params.cmap_params] # type: ignore[union-attr] - # im = ax.imshow( - # colored, - # alpha=render_params.alpha, - # zorder=render_params.zorder, - # ) - # im.set_transform(trans_data) # Apply cmaps to each channel, add up and normalize to [0, 1] colored = ( np.stack([channel_cmaps[ind](layers[ch]) for ind, ch in enumerate(channels)], 0).sum(0) / n_channels @@ -794,7 +713,7 @@ def _render_images( # Remove alpha channel so we can overwrite it from render_params.alpha colored = colored[:, :, :3] - _ax_show_and_transform(colored, trans_data, ax, render_params.alpha) + _ax_show_and_transform(colored, trans_data, ax, render_params.alpha, zorder=render_params.zorder) elif palette is not None and got_multiple_cmaps: raise ValueError("If 'palette' is provided, 'cmap' must be None.") @@ -829,7 +748,6 @@ def _render_labels( extent = get_extent(label, coordinate_system=coordinate_system) # get best scale out of multiscale label - # if isinstance(label, MultiscaleSpatialImage): if isinstance(label, DataTree): label = _multiscale_to_spatial_image( multiscale_image=label, diff --git a/src/spatialdata_plot/pl/render_params.py b/src/spatialdata_plot/pl/render_params.py index 2b0bf610..53b3a93e 100644 --- a/src/spatialdata_plot/pl/render_params.py +++ b/src/spatialdata_plot/pl/render_params.py @@ -81,7 +81,7 @@ class ShapesRenderParams: scale: float = 1.0 transfunc: Callable[[float], float] | None = None method: str | None = None - zorder: int | None = None + zorder: int = 0 table_name: str | None = None @@ -99,7 +99,7 @@ class PointsRenderParams: size: float = 1.0 transfunc: Callable[[float], float] | None = None method: str | None = None - zorder: int | None = None + zorder: int = 0 table_name: str | None = None @@ -114,7 +114,7 @@ class ImageRenderParams: alpha: float = 1.0 percentiles_for_norm: tuple[float | None, float | None] = (None, None) scale: str | None = None - zorder: int | None = None + zorder: int = 0 @dataclass @@ -133,4 +133,4 @@ class LabelsRenderParams: transfunc: Callable[[float], float] | None = None scale: str | None = None table_name: str | None = None - zorder: int | None = None + zorder: int = 0 diff --git a/src/spatialdata_plot/pl/utils.py b/src/spatialdata_plot/pl/utils.py index 7ed73e5c..8a9b1a3b 100644 --- a/src/spatialdata_plot/pl/utils.py +++ b/src/spatialdata_plot/pl/utils.py @@ -1939,16 +1939,19 @@ def _ax_show_and_transform( ax: Axes, alpha: float | None = None, cmap: ListedColormap | LinearSegmentedColormap | None = None, + zorder: int = 0, ) -> None: if not cmap and alpha is not None: im = ax.imshow( array, alpha=alpha, + zorder=zorder, ) im.set_transform(trans_data) else: im = ax.imshow( array, cmap=cmap, + zorder=zorder, ) im.set_transform(trans_data) diff --git a/tests/_images/Images_can_stick_to_zorder.png b/tests/_images/Images_can_stick_to_zorder.png new file mode 100644 index 0000000000000000000000000000000000000000..847266235415c5a114949860b64f0d99e3223986 GIT binary patch literal 29437 zcmV)aK&roqP)V>bE;%kREFfrfbZ~PzFE4FjbZ~5MbZlv2E^l&Y zFEj*5wg3PC32;bRa{vGaCjbBjCjm`V*v9|>a%xFLK~#90?ETHFZToiDjegqirT0F@ zoO7+c_jyg?dCXOwD~LV?#XR^Ah={?|5-hQop z+xPozVrItouJ2u+!xetFU%>aS|KGX37fj!~KDX<8!SucBbGyD5OrPQU)-$Tf6q_-)KBpT{=gsLzxRm;W-S(}_R%C;udWRX1{ZW43@B4kc|Ni^@;1B*F zfA-J*S^mmj`77Mt-}56s@*@EJ!Y}*+f9fozyCfz{nJ0qul?Gu@fZK%U*wPdu|M{^{gVEdzkXZ6^l$&|zvW;2 zi+_Qc{gz$K_pbl{yM9~8+|T~(&+_9x{^R_^fA|mi_P4*ySHAKUe*Wiwo^OBq+nmp5 z?(XissMqwlr^a=Gva|KK0w`@jGDKk@sYOU?Alzx>Pmzkc-pjj~dB$!osGukj83 zE5Dcj!VmCEyyRDbZ}Ey({3;jzD}`TS@aqgd#@K-d==>K2h!8S_74inBGd!H&ofY0q z@KWI1K?)pKc-)a!JMy{^FO8)Pwx-;7pt!Z6qkUK_ju9)TUGK%cL|z}u4u zvOqXPdI{^7;Qnjy?w8>!F9KhVk#~>G`q3TZ+wbuB{?~Z__4o7b^8I|smw1&o_EVv4 zAENocWA@*X`s;-KW8&pX+z0W9gmwr0J=A`H!S84A|73*UkEZWK*E{T?!+RgY_kEM_ z1OJWiegBR0?ys<;sDVCv5(f*noS4_XvDh z7QPhjc#+>DZQr22_yfH9%KyO!-~Wer{~P}gKKg+dJbtB89==6ezKMlL{y+cWKm5-H z)1Us+f0`fv@gK*`_{vwl0>B^t<9{50*C>NQq4CJhg+~e>qGSR)3CKK>_>jaSodcyI z^VYpNYh~WW1T4S_-CA_^=s<;ujtXR(y_ZMqyr0uyZ`g+zac(; zTfy{0KlDTV&=394r~cWuXE0SN7p}bKV-7NdPUMaZEL>Umkish>R|XA5-wLJMEsgV~ zt-u}FyTF&g-UE98UH#M09l|TH1%g1DfrMchtaq$~Hp2>NHMD}-{7`S-*5sCio?nHD zhq%q(5CNING@=L7o*_rDuB157(7A?#$8cc}uc=WO%YtqV&aYv)5IBgjk;Y&x%9@lV zK#~S>(B+Y8A0l!OoG{54g6d1?3+yl8;|6*V9-oLjlFB2?>kHdQ2is#KmOyWY4QOUm z-aIhOU>Y;P4&Y)<I@<9EDb&xL!Wy^3nSgIn!HVd_2(YV=O)T1N?vf#u)hApn?C+ zKoJ^+oeNJyt{8(IVV%KLXLfo`U`OZ(#+#gAJfAIOunm7Y-v#|DjITm}7wQR$LbyOY zK^GVU3k_KYmK0M^6-Mv$L#a_oatkD*b)hRLf`|~zKa60vtYjXqGpIy7Uc)V*0+a}s zVDz8_84RiZ?&XmqL|Uwjn6ca;KOxI&#ID3q$*HqOgM9dZxzz#yOmjLe&I_x3?{JOR%W-sYd_@TajA@7@gP z1C9e?pRu>VdWS{BB%^7e#|t76>B7iR7YN;tRf`eT47L)-!ZH%e2;{B=?Wk26g-WGc zM-AdomS(I4@(w9_Q9gE~9!jeXy1`IH8$lAW2~rG|4y_^DkaBCW6xuA{)MgQxjwyRr#HvB(NThR zgagw@)B-~XV_eX-VItUU68t`Wfx+~v^PmbC5r`;}0SpF_LgY%O5NWIoR+NkpW_99y zh5vZVP_y5#D~wlo9o1nh&<3O(;uE9_OYkExMr_&8v4SM9fSJ*QQ6q#!kq{7D5j&CW zOzd~8ZRM^8HW|cT&~c%1p;(~Wfes@zV=HjF7;A<&K!Z^axLy=Fn8BgmwOL@OFq$z? zEDU6L5ZP7f#c0jn%I)lNctYqw&?s0|Y;D;2fV)Rn#gCfu3T+2HKxIShT};n3D~MfD zo@hK#c%tEFhG8BQGoa+H`p_%hj#$OGVrHl{OrN|75SjPs_P~FE!Sn%s{)Xp3x>Z8{ z!W)4D;ULpEF}TMohB#l`ffsfc-lWm)fyBq40cc>aK(l9n*!(N?`TkZ28x}g)ilzvv zSZ$2a2|jd*R!{par3a*_x{{x?xSRVjf7NWQ2f4fixgrVKG<- z@EX&P(6lo`pvep1#XPzifhUxY7<@oO=;HOC%A3pt<{hie2E*Pu{N{;ZG|VnQ!bEse zirgOdFDRHUylL!ZhDew{H*8h|k;X#h%-{ve9Y&sidY%8$t)}w-0(S4y8K4)i)q$Qc z|7sNfwavOo7BItlf;dnPtnX+qhzw+fu&%`I0m)y+_!^;ng?M?%8WSfDR{aPxQ_YAp zl2xoXav9u3;UO6J1D3^OeDAPUBNZd7a0H=5rDdTmP!bf3YS7iN4!s+tgAM5F5E%%a zKqFF#bmF{i4dDtD{9k$m9x*B^3xgL3cff)%2)x4h7AgT{U{Wz@6JUD_*L&-{d#%Up$^YO}uFox)o>AlHAMF#kh9V3ioy=fkaK^aDfyvt$T;~T?-})D4ma-09 zz>)(#i1+BW0C)aYf`DZ2=c{6I1u{lLBtnYF8n71|JUkG$_fYvdcKiVH>OBI9)UPb< zMCgO4JB%kpGh;+@6;=^W25#o-9fQZa2 z!*i0c=Y(S(Y<3Hn&LkhkE0!)XT935L!VMW6);0h0hCSQaD)SR;DBAn%^g`v)|<$JoEZ2wz59!j31x z(MaWE=o^EDVG$WZ5aF55G$1qx$;-7(9PT1{FVHrlo_rD{VfGI|Z_p=BHRKI=_~w>* zA#OU5(TMFJlqYnrXc&xfMZ?t_kiH?72;39oJ(9e`A`m4bTv6U9gcHOCOA$>0ivgL% z(sKj(Tn3p3)a?y&JKDbU`hqjl3*Nx<1rMHqycJltzF}Z=Kabz8jkr~zf`{gJzUc1= zQxV6SQiRmJu82BQfOuey8Q#HXCWb%~I0ft>p)WEcF7#!kwS{&`^w)`zA~E_?9Dv|( zIAbvJGymWOpW=+@2{R-}-XIcAASeHM3y`Z9>~a86u;2-FgqhXEo*%P7Cu${(E7D$L zMHw+MX!L$CQt^AP5Tp|$}#5jPyJn zczZJOo%iE6UY}bqeT6q4^wa$2w_AV1GCYv{Zy8VCPCVzd0(J|o%o_~qc2>7Ji>n?q zan9kgdSSX46YlkYPtA+;TgyE^920}6SnlZB8L80qN?kWf4wUr78ZU_KAQhF)&_Tn& z$u6MtpRY2zA)__JoQYiw*{CsVuD${w? z!FLuezt#1*1=IKAf%FV+q{CDebUhv7|HtqjnD%W6RJY*!We>@*DMzN?>!$Y0-{FLLAJ_?E2{?5 z=+GB>6^=9rx%vA@2@6V03o&eD6IOfTq{3#e2;*Dy@G;~IeXZ2967NAy36dZhu_Jh$ zsWa`4iRUq=;hzYo-}?I8g6ZoB7a$TC2$evYKxMym1Y#aAas}xg91?lXKDSTUZUf#k z9#AuLNyc;f0oy?C61+IWiw(J7VO>BH)_`lkp5!g0!DdaR*l}Q&Cp3S^$bohg%5le0 z;nH`4en_`(6YLeut{gP3D0TKcS)fLtEU;`wm@ie=H+9fZ=uN05>|Ho&pq0qz{`qJE z7SEIt75Y52hq7p6v%(oztigy>g`kQ_XCy<{jOLD*5t90tX55oHWW5l09|_;4#YbW; zgjJDKK<)#Q66k`2g!SM@Txk$0gaY~o9{61ergst245AB(g-9Y12nb=0-SGk;6bK#S zF*~#QrHdC48bLSEJDhoAUTlQmSt(~^Q{+WLzH~;ud=KxO;5>t=D7*&yaMtEe25lng zhr**_Ie4NvL?5pR_L|PlLE}KFfQisV7{SLHBMWgcY=yBH&7gPB zM7;{N2aXmwTBNi@?+KMzza;f+&4(k|f)>u&>5NxtCWDt`PmWZvBfh7g& zud%YB1vC^0uTgo0h9@jmAAaWudGW45PtaG-=1CzI2)o~p+r6^4(A3|_U|I+)1Tu+8 zCbA$b9%3TiLEJ_JNgcILHhCYxE7)#+R&;yQ8tE)L1V76INM4W=gojLcZzFu|4*Bv1 zFB~4wL$MDBnmTIT83j6T+&;V##;{|csWAF0+Eb=J8OPBH{b1-$u$^I*N~d77XB15w zY86Ap5Sdt-L7gMRx(n?PYKhdIsAHw+N|yyq(Ss=r#5(i}wfnJL1SyfwB5{OCiy6cT ztcY|()QkA!Or?=_dU{QZnGyt%PT-)kQ|zEwXQ-f>G2Jk`2Yr~CKfO66LiCzbY<}K@ zH$0v1uFCn`gQ*basKjw3<&0;Q3C74kLdDVn%i%Mor8`I?v-@hn??h?%SyesAZuAt1 zNZ63L5_2SMnRLGr-&qM?S&(I1GFV!=j2k3hdb zFN4}UXd}?(xZ3Cw)(~?5PgIa@#HPfeXjJHm4x{aGT!pgx(bmRFw-Y^{8Os@6Gq#Ek zr9;7b!AgZ)v8$jlF-oGhNE`+-=v0PvbSOGHHq7~2(&?kn#+6+Yn9$7->kR8u?UW9! z2DD@}ZZPh^UO>C|ht)cSWp=&8zib8*c=_=g@TM;?m_8(OA@Rt<&Y8vogL_OIXc6fU zuimC^!zs5XgJt+@9KaF%2(%sa$=`(a2~^IUAb}v67*05M6cmt3ou!v|RbRkwH?t>IMqjw|@ ztPe(kb`|OosWs7cVZ<}KygYm963d2gf_?`*J31cKYf;^2%{nig z3=X@^G=2wy>DLJyB%WBfa-wpN@e<*VAOQ)L5F?rftQ9*KxCge`&bQl*jgIO~!kfRj z8&r@GTL~73npyh7)*@RS33CjgyL+Kt&W}ex*64D?9=Q2fu z6d*@qi4cY`Vn9b?kTF@oR%9fk1>_irG-MQ5ibDwlMpZ@;dI|KA9QaHpEZt-B63b`U z0^C(t7mS3~f^J67nG#QkYz&EHIJ>S1FV!KYcSIs4g&sPV zc32W@0XvIlr!c7))*QS!m;hwK!I}kZL2UE(xhGOlmVIz~G|sOo>-~k89vS_SgGZ!3 zVdIIx6_o=Y8tbHoI6pcB2vJA@k_b7#BK{pWgn-bSBkigu)b21&n6BQHl|i>1Yga}R znwTfP?w#D|f@(k|`V3N6%uc}Bzj;oc1QOV`cxM?;h}@%grfDRO3i^n4p|wV((6lm| zFw#7YI#?C7Bt~74z9J(dIAh|5_63WM>5hRd70m}Db}V)@6-W(>Hu=dDEVDcRVw)w$ zcR%VrHJEl17(wPVoH^K=jXk*0*B5g%7lG_n-( zxDfY8`okk_xf1(9d@O{gPB=QLc~%N$L~*z(WG89HJ0mYb5@)o8RR}o{awL?Dloi?; z>kst)62={7kufS&1EF_@6e@+o8r77p#xTc`%m&scxvv~P`MpN37vu2xVa^cl8R=bY zd5@9r30aW%X{mk_C^_v+rynXegBnFAjLLu1|EYax18BD4d z#@Nwb(blmt-ejVNA9dY)USVSfRP?&Un)g$404n8h@*7}Uoq2e@pbtB?9GF$u8@yJ+ zMUisj8bJ;|(hCc)f!w(+A%+0msjy>B=(f;~m0mVRInj?3y`CBMOzV+4k0D_}NpvEW zK(#azFFDa|+A^8*?}dQ#~_DKj!t2 zA4_Mj2h8pfi9~B;6EfCd&L6@Bpg^&eS#bjX8y@}2li_07uPN?7$b<%I6X8M(1 z`4#@AS4Q_MOF=n4w_yG7En6PcG`HN zw*x5+S=47`2@j%%b26T}GwPU17>GIHS7Z>hBK?XL=W;gfgjNVgp+6SN(rD5db-=EK zCqoM44x7RnJfI|bX`q{6CxH|g8mVhWw*@(`v~uDo8)aXq*G#=8#x*z|RA+yqeq)FE z5kD+Uo2Qi-+MW+IWpi{2n z*)LipSb`p$7o-F29ZIKPI%OFY88mgB9J@D)Yj=j^2`{vRz=%SNLd%hnBUKjqx)RPC zj1#?|DD}*7U8t`E^_5GS%s_MpEjT-IOc_-`&AF!f%sn+9QFOdoX7`qRM4u2BgYD$e zZXK%)u7VvEsZY?Ku-o@C|6Xs2{_s$q58r$MMFv5H7zSAeYfw&Muuf`j8x8A9t457R z5c4^D3zNtG1`6@-QQj)Rzu<1p-g!$>GjL*atx%%B$T906Yl!yww9VO%%olePr6eF3 zqDBZZGoI677gMq}vPOal+M$`!hRzR1soDA74bU#mJL*O8Ll2meN z^W+{S^DGK)^viu7)nQK3WilAeSk+j%kXt8ofq4=*OUlcxW|lZ3|BAW z8Qu`U8*FL_mL?y(s|WVprn15OHQT{&qy&9|!PMUKZxBqO8+}(obo%LWQ}qmE1Ve+m zBq`58T6tqHU_2jpGO6Avn6nvyQy;8V$cH~cqYG8V1IT>P+z8L=j)0-jl_5G5HARy# zqTjb6GP%-`L8i$@i+7CWC3=3r*fwXaEEVlLdac;&E4J-$+2_!@&7(SC%&MO4%^Q-* zWoP$dTqp6{`s7^&NNECJ9eA9~Imejb<5_PV-i$m8#(3thN@)IC(iO>pxDLXSkZf42 zs0_M-^az?ML#T(++n~#EuJ}-_%>knm*53q^{pMNea}B0)KL?-YStyGj9Yzsa@<~aX zXsiYt9J_(7w+uQ(0&>etIC9i2A+*Jcz{Xp z7jaxii2jD8`S^n2g}d1d9OhiC%r;Yi25}yJ>JUg@P! z4`nP!0anj*h9TL>#KWaiy2}ao)^oUVyA# z42SFNN{drzo;5oT9+<*2AzO~g+Y&DWXrKk@{n^naB4iZ@S%Daw(_~rElVW#@+(WvL zET=?Fz*x|-F`fopmC}_%2emcY6m6;=Z?HerHxP~e)`0|{XE3$5hX2j8wG5jR3~@@f z&C@`|#U{;L_D^J`Y14Fb9g)nGCwcc3;$-bML24lONRZ@zX*t1q2d5LjX@8J68Myf#B#DR!C6gVvmtXoENYDFH|87eqEhHt$5E zB5^^YD{D+9GOmhF7N|`nFL7oUw?^^lfEqMGN+MRT#+E+VY_JOCD1-%U6Z9@H9x}XK zi7(FN(@Ktk)D14`)Z=Q@!>HD2qcHk~p^umvcBAf^@`mpG*{;v80Dtr!)!xu!T)`+_ z>HA%tOrW=hmS91xnhRz_S}q&n!KQ}aNw3l*#@b}{5*g(t-horJq7V6ZnW zx;v&AG0v}Sj!|Ym-gWZH$52K$R9z%#NgWkGB{*fDVaAAzh?GQ#gCL1ymDL(270#lp zNql&Fj`aJ)c)7qkXTpmI@^&VLL?|$Zp}Uj%RE%Mrp;vUI8AyFX56|#A->umGwCnQ= zrVl-vJWoO@aP?Nxs`y-L2&Uc*_^8ehRnPvjb9+`e$^3B=FNZ)4Lg_-VK$nE2L|a!z zIzhNY@&j>sNlGuUaPlz)$DE!wx)n?hbU1Z37oQJCgPfE_l{hE0%{oVzv=T9NOql0M zW@BaEI*xfHFrMGI$UzJXDJVHBdELotA*>b4uKI67Y0X^U*KkGcjj?!rUz!u zVDK90?d!7z)aQpYefZY5RlGNFuGs1n;CCM5=+x8nCE2hQ=gW0(W07J(SBm>hG;caUC!+hX#3=Fz| zh-HG2yZQSsf;2-$_u;vXXGR4#;Q4>x`uu|FiHT=T)Yrlcpv7eyx+0dHj$%r*PPW22 z4%Q&E7xKH@F>1i+rXA<*WWjPkQ$q8~Sk5fV19^EMr3a$j!8oIRp*A7bSs?cV>VX_i z2M%=*)%{m;Fml0v6eF3DhGE@^Y6LaH1T}?#%n^tT7b7@L$TO3Jh~gcnbdnt85LhFe za_78M(pJz@B}ivPbt9Jot{5NVt2T7 zX+wxE4!1(q3#$odgA+=cPeucF@%m}bI2l>c9I(|pH(Ac4{6N^=VaX3Hb|Th{4xx;O z>15v?l86Cy zv963-(c=mdh%(4!CtVWq1nM6FIz23sc7`43bf0nA#2a<>G~RjZCjSoMOg`kgeY2@H zVZjv;FMyZa3?@VdQLZG3EWtR*;7n)19oninq|ntze0|)CcPHCfcE8IYccl2hnjYA6 znVyxd&)0IL?H77^Lfa$Kp9te1Yh$&+X6_}q3as2*q6Vg&Zk1vO)pnZQ49W^-liI?B z1RNlmC^XN(suE>@?X8j+l;i%o>2^UQ(Y3B>9e?f6>8ooA}=} zaL=BUXOI`bdms-;$Vlu+d?X|#g-QyA#R`kbvx>3Kv0boR!HyU2F}MtR6b>zz8`iNe59HVw1c`$H5D~m(QH3 zPpDHqJ)l0nVERNC@K!K!1|AT3iO5Tc+4&&ZVA&x(y4#xt&ptPm(_3B^!kk3hWOogg zVl2C{h>@FfbHY*4{fZqItUf}0jkF7496%w{*uvmWl>5mZ$P<(srZk&ceoIsfR2ek8 zGAy7p%;sproOUA5N88OU@0svGG{?`hLXCz_K)WA872)VU$L)upMX)+C5SbjODvrGg zo&g8Y1J)LF7+rmy)@evGKPw=S$vt0Sm$jKi=*kJa|{K?AW%WTcWR#76a7;x?l&?S(Nt1M2N-F zZ~V=#t+O4CO%9fEB_?lkN3FE$BW-`AmnTNQpt=JEkw%JxH7ciQ+y&!|$21Ccjt-ht z#KZ@})=(kwd@((q63W-*CozgSmtaUf}DB3!}|rG>Nbb znQ||ocXNgKJf_2ZN4I-t&w|=pKz-rCugj`{lYdNmj<%X70T%B{a|lb4n3$=>7hd&*s5?g!pVr6 z4~|!XlR|7p>+YtDb%q&Ll`;mU1ws|18d4Q$0~1K7ScQx#xRO0 zI+PZS65O}4S$v*Z>hvbZGVV)vHx0b*2%HIcK-uJiD~6g(5=y`{cpyy~zH~#zBp$}O ztznzPo&Yw2`=%WM3Bd)|)5uKhH`gAeo=bD-1`G0{xB$5b{Jg_7(DR)2VO|prC6CV1D$5X-ZUk!#OLP^ADu*@>EmRsUC?grN zf!+IOVlVI$cX1y0Hr$9yn>o&BZm?cdIC|$Ok+O>$w_Mc+uXR2;!w;ckj>xZUoY^>W zBH4lwC}wmUpo1Yw7k7mjB(KHn&2ZJjawBjiaVFWB#A&)XxR3fUuagFaK|!q|v%3uH z^pr3R8I!3xM_0@NL>umZ9zrC9NSEmRFJw-#!X)->Gk{(I@6I1@t~BvC-uzDsrl0wl zpWzpO@fWGJ5~Ke|zm;ar1L;j}y5*(Ym-hq%rJ>Ro>`tngB>_r%w)yfpVx4Gt7<(?;{gaVrKF zgQW)29tcl?@LFzFHDa#C7c)eKV2Q*+=ETB1EB7Q$lMC)9t61|HAAvze$rHbtG!yH9 z5m-sGkjX5hX$E2h>(2#urBJybJTh)vC#yv5L{p;LTsUOUiZM$<7k6a|BhfVzbfH`3 zn5vpSUk!8blW?8^^#HuYJxJHj9!!7c&-@vF_=kV^d8O!YrkMktb#u31x;a%Y7z@f` zjU*LeXIQ6crRjw_uJpd6tumw|IIv^}1P>@~6qGoFtD1wuX7p~<2^iJ^#|RV)G)W9` zXAemM*@Sp6#CJOK72qq1y$5>snV9s+M?1vnz$#=SNfV(h$R#7Mk|Uck`90d^*JKP# znTsD$4+t*_(@7YHcUV-fb0qpFGmkBvZL5#zUgG1MMU)er6NSYFQ52<7?cm@baHa4_ z;1RPckxH^o)wMjxA(+gJpZv+6tfArlhnc!by7Rppk>|~zix~fFc%06CU^$|UG zG*)-@5S>;!KIkO^i}RRdnM&c}nzW-vBqh2mm~AxOIBcO>cBi^92(c0JLD;T@`-Akd zklr)+ilJXJ#@(F%b?7^U!zJBaBSR9yR1%u{2PN7#laPviqyJX{q~KDp4wn?DN2DoQOnI)asP9?qF{C*tx*WEsJ&M z$PXhn!HN_AnhC80hD1tFv|~l~JE-?)zhl@+x4^It(=<9f-_6?wndZPP@R>8ykN^0O z^UZI5lVAMBU*z%e@i*1X!5iJ2%#X%1(_o!}bmo!O0)axaK+uNvM|yjt*DF0%pE==U z3hVP9o@(-Edad~6{zwUlD??5sIRp1}?x}W%TVq9n3{u#M>5;g;CY)Xq?u_)HgcpXs zFzOu~(YQ3w>|(@P0~d{4$Q*3cen(69sBTz^PK>)8kT28Z^uR0i8Yy%ZDjN$Yi)010LDV`k-40f0RG^UKo`UEw%m?S;?%bn^y)N@Pdo#x37lP#8Q-v)hJGH9uo8{jIhkbJW#PhxurUaBA(^1oskW2!$=6(s3Vm)Lx{V%g(K&$O<`08N8yV3R z(>t2IAlf@%@1piDgO{kC5eZY)^T-lEAg4^o=A^Sc2x+)4KrkL>QMw0^AQTZU7I{VH zi4(irQ{zi|xvpIB3KyiKtiJqdx_2`eroEZffPS$ zw8X7ac>!#J64=*@J@M=%-jD#rqw3TYx4p*s7a1_2UdBzn(SU9on8 zafj%cU_!Ekgic=q9~=2Q95l*YQ<7#E`5A)gfA;dX*m0eRc-p2SB>%Mx!X%K%RCZXH zP`9~Az(FIQZX+(ofe3ZMjx(dar1dWm#@E2U27DFlO9*E&fy9-K4>-xnf)cEtdeBEf zTcnl2by^7?4blNzp&r7u2QKD5Q&&4t>pi_bxJS5lFv0M7mT#}}54_^^(ZQ8UCZO}^ zQySG5`k=TE`rH#Ex1#L|=>qBM()MLwIWR)@<92M&Hl!`2mRNiA&qS4EZh@U{8)+o4 ziep9tagpj^+D`{JHtP``5z~lRB3gV4E0b9jb%aKzP&hcAF{baXQu!Q%=_Xj#IZ$s? z%Qk3I2@?Im3ch|+bMq&Co$hhW!KcErTke=s^yso@Z7cn_!}fO=^=pLj4Me{I_WcCD zOyca|Ca+l8*tD^Z*M#;l{aD}WI||~)9ha9vqfcXBax-%N%-2Pi_F=zrkmf3DB7cn^24lY$TIV#jjDwu;_$ z`n^%lfg_8P2Ad&!WW6q=>xr}{By>VBQm-`i!E=z{{O+8QC1b1mKUtnijnR?nMLwg# z*a;Ge(tK1(VdKh~g9XR0Dmu>^&%U$o{+$S>XF)C%BQQ|$WG~3%aT`Z8S~lm-t@yC$ z#u%8={Ac(4F}YEbFed$kw0oqzLl|Eo*q0Ie5<%XXK2R5snGtL-`hYU%SEW6gfByzG zPm*NutYl9xF3LER(L2K`g9Fto+QXDA}k9zWnv7(V{TLkC8{CqNJ`6oPko3{8M??34w&ubEE|>*On|l@20cD4BaaWVc zDK;B3OvY^W8@G{z&Hw)-b1`^8h2Xx2DPlRmn$VLl?gH&TAm@a{g)k!F5aKGNY@}#} z`4tZ%bV1A|%9<8*y+_XX*xicU32{A0lE~H2=pCUor}WA=h#HB)#K314XW&O$o~Yq% zj3UqP>I(}de+qq$KDbfEJi!FHCc53s9YPhd}mT?bYh zHp(2R8oqLT`U3Roz|^<}lM^<0<^TyceFr}I!OtIMa+YpVw>Kdx@OnRtpV8r zk~6^;;u1({kYXi0uEdF35T!obHhGiFAg%B^DWIItaR>Wnj5YhrR&e z))hJ+TrK#*LR)E^sdk>_{jFyD^y_mBrdyRC@QvKlqeyH@q}{Ya&1FN@1ontXhPn}S zg=4|?G}owftTpJ3FqEjilz9O*Gfo(*C{YFyDu^5?E+XgVAO_*(Td*I*&h10a<_Wb*4*XI^YGSi^o9+o#C&J`46Dl57> zpeRrS5)(mJB%BewL+WX6z&5m1tRL9e3By;EM)UQ2Cycurxh)u^!B~}gccI_ zKYz$rKW3yydfaKLQgf$erHiXVbi4^BnmH6v&9maoWXn3ww#9S1)3adWQ#JT+2i0%* zka+f;N=(|qlr)_mV0(!r_bAXFJS)o}WaAJcBw(rZt&;nW*b|zrj8VY6KTOgfG$gyb zaYD$Wq1FlAy#z%4DAv&F(SZsdr$owp`)EYRM%O#4Jy5x$vCIPe3l1hnC=x+3K^8(t zgn;vEV#QMNl$?CMDHaJe6UIj5jxb(8d%%trEyidEqd#Fo5Zejn{z?gB5uB(Yfh8tF zPV`fwJ}BkA!tw5v^6-S6F0`=ItN4Is=iF>ZFc94ldY3!cv}i&F8_Wb zQoQ};+fnx&N7~!hb5>%GiUcQBCAeSVd!>vIdH#2JxHloKgz#JCW%AgM#%(R^j3D}4J}TkwJA@xD^teK zTh1CYknAS0J>zXrnn82SS8~%k-OMXK!*wH~X7P=O=LmNh&eJPF?m+LblH8Ooc1Et$ zbYKz~qfwg?wL9p$6eLzs^l^(!uLN2rwq<%fG_+r#T?t1al};S$Ysw;oB#?w5)0(JF zEMA-m~NrX$^#}V(W#zW?B|ngdr6W;)vv)SS>r* zLRZGvXlfj)h%HoehuPcBy19p_%-wOp2+4?x5l@5oFo^GV!k0UC?hcw>)rT9!-mLw- zf?aUE-xU`sPN4J?fn+y=*gf1>0POQ!-$1d+yyT`mjYW-6-Kg`Jo2=D*lGSIQs0@bmIiA=lTu40mJT5e1yde1s z%Yki`zBQ^1j;2n$>*mXD9(fZ)x39uH+TISPNuzz{^sIOQMdm^~m_6l<7I=%9+x1*Q z$(vt4|E0j}sJtL{&PZ{8qCH?oMurcE3CRm{@eQ3449SXYgAfNJ43zGs8@*9VA@0Fn zVpXJfQg5yYL`gAVIZ)R`xl6QDqUA`}NIxPY8Y%^=PhbN}S3j0lVQ#1X4cF%uOea}b z*|-a5?$Zk%mKSWNOnm5!mzm?W&t#1RlH{Is#_Z^h82{%lSaJe`E)maMDfknB3{g>^a%{^?7Rl~Hp;@#LLX;(?X1@iNsr&eKKhU$PuO03 zID1j}*x(iL8V{rkZ-QtVWPEy)@^*Lj$%E_dkh?Lb-(FX93#40B<%TH&Be*Nf+p9hL zJb#;bAYFFYmkHI&%+FUOT(NjSuE^ppZ`Ksee)MaDz7OJ62)np~PPr9esP@U=_ez#2u#%rJY$za$kdS7WrXM1`gFCQTC z0quCAY=vf>zANJu=*Nal)7no25i=@o2I?~X$GPDKbIL#KHJShM&08|gvF%L|bVJ;3 zqa`0T&d{@fj^G*Ra5)nrPLtRU?SP}fo~C!>^b;}#b}T$Wcyg-4(x7#;D3aV$UAyC{ z;SjNT7)^!ti=n0veoV>te_Y zW`knItvkM!3U%2GM#v}%m2VVBvxm3gOyH|gYNuN0hX4#O>`DY)Bg z#%7Ui9Ner}Z&bGxdp-znOCj{FFvM_Wi3X%6BxGz~K~GRm?kT`T*<1IqO|_eEqXVZ@ ze97Pi`hkt+@}23ZBGaG2=bHBXz(|*Av+l}*p$Qu++U^+T1!I4YA$MpCUSUYJP z+_h_F@vf3nI&bCyDN-=^(n@{ZXO!YnLB^Ql|3qRhVCoRw6E@45Fb!^ ziS`F91;_M;lKRd{W))#G$TCPWyk7gFySZWkF~~Y^XuVK&XF2u&7j3*gl*?*d103o* zY~?D1$3WZzQX@Jt%~q;pYK-(O=&IP#u;d@4jfz(H6_D!uGV8`r-(W;(jM(W>eXwli z`ZRO!`gSJC&C^u;-;iSNnF;P0Zc4D*#9t>j=}omH_O?QBT1(HbTQJECI|tca0;cM^ zG_kogslPyE^@dLm2OD)@qaf2P!)*F3$&{OL`}HjL%!la#-q@>~XURNYBLvg}O%wFY zkPYO7>K#_@AqHO#sR!9M*0^vNGUpmti;xc^TxEg2~vw}$W z9c~1og*k4MXO^_^wy0WeA?OT=xg}k6L&j9@;tnd#R50I!U#0y1y5sVJOT7t46uoRJ$?tk-5JkuXmI`YX zmKcfC8fdW5M9;=lIoZSAf9DIoU`Krrv`V`wZSI~B#W43iY-JoOgDPp=(aI*b??j|Q`B%G z#4{ulG-^2%H~_%uFSF#WrK_wV@E|N39^{`>E9KA(Z#T!eYk%@zJPC!O1s zsyk^^V)Kf@%;>$+>Xj}Rnj8%6=)m3O+8O!_=cj#C~>*fq{q%FEWK}?r$Jg!0N0&`;;s6$uuAS z`osJt)7Dx}5H@-CO5Zx%7vUiY4;DCUBAI8Q-R{8dY1~uwo{Bd`Vyh3J#X!yw@4aJ8 z-2-=++c=syN$|4-(;xk#f0TdoZ~hIJ%Z0Cf?Q6fK2=jJT`TYbQ>3l%teIf<1E20(C zhK|ax!jJ>CiVnp_gpr_cUWnRD;9X2yifj_+6L>&uNL+{^lSXEZE9Vf|B(h8i`Y@<8 zqV5FZ9J5V>81>yLGRuow)6?yfjh)*yz*&U5084UMvPI#DfgC-H1P!znR>MWD7YPW#A&Gr_MlMwENcsA1N z0QI@S{S{skc^HZFSO|Ke={?9xioHv+_XxZs*aKqTHdw&-sH|snd%-vj+R3QXo59Qo z4ObL>mS6(-u^;;}{?6a|JD+KB;Y~L;7Wqcak zM0HFLM_cEkwqlm%XhV@PknUeWjNFBFWKNp7vxUtz0zT0`tdPhggC-rF2B2?b(~iZo z>efAVz67fKxu<4iHC9tjqHM7fQ)lGHvJ_GZ1S)Yfjyh;nvEs#fGYIMjT$qTrFGre= zy3rsik{TftqFh)^*tBsTnQa(x6cl%G3DTVP$eT4G;#-NZEXW#Q6JNe~4|z4;T2+lZ zg!`7bFAG~egFGO-%fWjHUqbjY@E*oX2Y$q0*~X^^E8?7I)O ze?u_+?9cuz|M5TmN1mRZ2qFBIBFx*R>B|6!yoRv&lFGSfg!Kx!di`QPjD*SSmK({; zrPoIHk-koaZ!K_KAU+~`A;>|jgS^95GUs;UJnl()Vi>GeW*Om9YL!B;K$?C@(uGls zvbzDTG%sGY<0F>?LUVha7#cYi@^ZyaiFoS790*bf?TQ=+a&6O?!u4$3PgUaW^Pck+ zi4d|9mI1jFQXvK*o3i%GsU}X%NOaT+QL3BY<7B?@Rv?iBvP4)TvIW8>gp-k0Wx?$g zPZhSDS>r|u_Y8R#m9H}R8iTLWd6ywCu(0_6cL6(f>;cBhNPn?XAI^l6iyfLA#IgIX z!0{QoxgY-FALfUD_=i9BFaKsDOj~?1Ws@&&N$0k}23z%|C!ssNvCI#7u1NJQkyp%~ zd|j6oM}^yfqd>SK`iPWAgbSfp;vvLiVX29AtSqzmSrHZvb~25~FrS+x!-J{$7%Yi5 zM&p=A&^(j$FkdCFzNAp+fHU1G;ch3a!8a?_E2Jw}rHsMo>Sht@wl>iw@m`#TdrQ80 zQrJ3S5w`E$~$z~7y&dj*X!L?>}A2;Yv{X5e-WwE*uLvRC<~PApz3D{rvKsV z?jwBgb-}^`xw&7>S^aBOB)U0;TyZCn-jL?O!pPnQSpjnsbd%u#ogG-)v2?|xVnbaN zSu$3Y(3Px4G-YKv>02yjf_Bum)M2-a>EAdYT( zjU=!o*j6}iz9VK?km%eh84dKnYDK%!MkeTr*_pyip7;_^e4VR&9|tdJ5*hJ8(`#%k zaC(CK173D`SK;Nrw&BNgbuTA>U7t}N{y%)Z`zv6tz}^QQ-M4ZzI0raKpE?Ic+5`iu zBB>*>AmK6@u?^)3)#M@2XPK8^rK?C-|LH9u+FxS2-`BsOl5Xt9%HBBV$RM$$n}ofM52TwxQ%jTbkJO=c*Fk%B>- zYsMy{Mb8GwSJwCFKw@bNqeiqTMH($c--nk3&5CV;oip5>;B*h`JuDk02^)hkc6!g$ za^)zE&?0ptR2L50c;t>J9=P&Ak$ajbW6nQG7i9SuIbDh8L3}RiyA$!$AWYLm`=W#C ztN*vF7{@EnCl7)YCs|%$nZA_fc!=MGK7FulQcMnpEv|ylNi}n-pcA}gxH8S4s=4x3 z1FIUF9bm#jC369;3<{0rD9yyb_=G}c*7(dXeWNiUDyZ|-eOs2XjZ7S_=GBgOqssJX zic=V*lc%oaDWXVo>*t8}?6cojH_7iKQO1HCXGS^^V@88fL*dY!8ibw%UAu3ox?kab z^I*y^A>3gkhIV>isZH36`c|qT1ieiYjSI$QPDO8z+F-}z3FKoWd_)XaQdE`@*+OK& zNH!1!?TOJcqkTr4>3{h8zF&44MXgh$C}5kTs$@i+`xA$dMQP^7ApDsRO3ehVghLm+r& z5DZ1Tpmi{gNJ~amH;t+_P!6G2W%OxHsEU}n|A`G&@%dh)`TIlP;$((t@v(4lLV2Ka zPhiEw=yFg(qC}%7L$kP*;5otF3a1Uy9cB+$&**6M0i{=D6oxep1|nYj)o|sI*BtUO z#TFVn&WQgA?0rN(B#e)V{UA43+sdh~Y((-%P*+CT=(RJ-7j0bewQo)?j=&{OF3##_ zM=eOtK60qWawuC>*4|0O?em5WDxHG}mM}7%xn+p#SwI{&>ei`rqIHsW7Rp;zo?c1c zoIGqwoE3ZKM&7oU)k(gz89%qAq~bzR8$KM?L+Bw;HIoe#R{9vUrt})<$p@}uWVD_b-9;5*o9F4U&`s&FF=9jWm}_Xmlt)c7Z%vpc zMlf1P5CfWomW`ALk}IQM>EmPCI8!yz>|n5=HZXmS^be5sA)!2x_Cda8wv%v9nbTNV z6ex*Sg;oQ-mQRa-&$DKFV29Rv>IRCBFnY&f54!Hpck&{Lp~ z$S8?^tk5?^PtHpUips)B2~E>dF$|}f;4q5Z_*McI8pBsMp`p{lAmEGx8^Y+3(SiCd55t=TpSzF z&Br9WKwaFIdK6N-dOLpv&PCXovjF>CdZUvyG$eXoCU=-Tpy31|V4>nZjz-W%(ZW?5 z+5Fhd_@v+MGCYy*V!rxrdNCcrNsc;exSJRF&AZ*LvBg(zG(j!WHB;M(5bhCq2h%%7 zU(haCQ_^r?bj+EWaS(VS6TJQJ1FeI0HHyMfp-g0RG|*zM-NPf>;^c3oxfnZmWhcEO zwIZ@3ahZ{yfsW-1Bk!mcdaj66ZmTw|8X7w~u9$v+;m#eRozNShDsfMw^+bBQCob<1 zh-mV|WLqo)?aJKn`FR;W3z>zRATOIw67J zpfthpRQm7*m(`%DFiaR0>6&P5LBcf2w-=bM*a!qwMCXDYY`*&4{f=yqtxj*KN@1`L zBNrnEmlT^;bht?dq8>bP$C8hch?DoNoluP+JD~-qf`kGR(KxWUW9fp*4rwwCn!?e2 zjaaB47fe>LH21qe7=hFlQh6XZZAje#Ht~MWkp1Idyq8ACsA}kWnFtzx3g=K{$I(;N$uMo#6n!n&+`nBfv`+IRR{Nv*O$9uqP#rhQ+1yz^ZG*hZNr6F-+t=Sa2+7#0o zk`37uP6j7sEDqwCD%y-*?0FZk)!?Ay5H5zSDDF}bJqGnyOCu9OyO#Wxy!vBC7sHo3>qi;~0OIDB_nYuLD= zcA@J*wZ?%GO);})621XZ9K{Zy6JsMT1v%~TAlQAtP6bPoeOLy1HF_4f1}sE|1WF4m zLkQ+mdJUbHI(;4Jc@XX(o`jf;JcRY|CJ8e{J&F4I2X6&nI8luQvkOjvGkpKA-4#tz}|?0Yige)t$l>Kq?+k}6w=tq(Q^ zSsiulK{&c^P_P@rpbzv=q|Eu_P=fj*My$K549Q8SmN|{VNsLn!au7oB%?G^%N`it? z`k=PX_%7$%=O-xmcJt$F?t-<}rYcpiF)%F*+v#>tX>X3+#p2GEqfkdd#}^C@#8=HRmZ6vqWYo!5a+!|lG_tcc z)4OuBExg(D1u^n7V(h^rOJJOpc7M=b3jJk6*6CPR;E)l*2Vd`F}n!48pLjd zYQ!DF?tY!A+^99qZNCwSVGtRt)+X?~v-T&Fy=Le~jJ~4BK@W|_YpOlb^vcjL7#hge zX}@f~mQ<~x8s>r-om+@p!rWb}F>j1L3q=g?u#zYt4kR`DTF4uuyGnQfdC`$mb$L|n z&YF)IlNIG2(|dHBk&(SU)BX8vPP1f$;ER4FkR%{M=n*U`Ep?7KsA0Gu@wVG!qEa?N z(OV#S>ku5#4v`)cEk()_*iX(i(o^D+PFyu{)J(njz7o|6+GduxRs8c?HHQd9yt$Au zMGk5Db8QgL<}AxR$5aKv#g|owsj*50*5*{+MIBFrH|B_fTj3|%DJAj>u@AsbST{7D z88Ya+qT4IFJ-Hj-7a2_FrYPduuyK@nn7lbDxOqC>R31*OeO9|p+(#$>hlm&7iz27t z8ltlycazp4-Fx^jJukG;0rur<8fMEQD=kzgy6|#3ct_>ELNE9G2-lPU)2r zh7#1AQE_W10?zviF2j=M-(I6j?Z!)^zXNjb%;(X3Bx+vDA?8x2WjI0MWK+{8{?Oi> zIMsZoFjg3Hq3g$pZFCD%DYSm3i_&>w@DXYsGq`@S!Ss&exUcz6lY;MUyB%nL-EbJs z5;y+X%|mO5@Aqgi5>lcsnYJ$UEi+ERkF+^7_H|f-bwRe^)Qok3b%S)m;^Lp?P^V?| z1*wbkm6FjULtoLb(1!07AXP@*>B~+XOiLfnOup9Q7EqnJS9~5&fhK`NB2SW8f+90J zt?PzuXZ9uYY8a1pu$P6F1Pe-#xwmept8asIplr;caF_vt|GY0y?jd{$##e!N-RGhd zs2?~^W>>ERV{qQq65tg4dCXHBTy@?foIlklU{?mO5F%X%hg>OaRDsTp@`&=&v$j!8#c*+ zLu3Y1a!$Tvn!uAJ7V&lRl@nKiCk~24bV38s8o|0#TIMoDQTIpc$<^U~0$AhZ6Dzrc z@g>OL2Yh)3*wr!Pe!zGQ@)~5S&y)E&w5o?KRZX1NHYY0uYEqU8l$H`;d3aeCl5 zZL}@Wb3;RR+@uX;x|@W<+=jn8jndKNF;}4$ck{3@*9b?Dg&}8*mkjO+%Yx(q@#4xA zX_!@18#?-IBHj)r|Ist58(_6yNDv4S6{VVR8Ijx;LOnC=jy3|NE{x{U630Q-Lej>f zgGHTP$uK8i1Ro?2p|3uzvIKU9_5$LUfUo-Rw-}I zwsseDYY4m7Q_^X2buny!Ap;G?$%%Zy(7@Nz%GHcKkA-eZHQxmk;S&h)BGb-MU<;5V z#3Wb>j2!7(rk+mJ(>>+5cjTBB0?tpujYYp?v zWpm4_Y+xs^LC2VMdc`R<1&$IZMJVR0)%xT;8P0B1QARBE5dHf{MPtLF`;va|GcDkeI$dCx55sroUlt{0H{94JcmGE@HqcAQ?KMq<;w4rouV28VVV9+Dh zq8Az?W41uZNES#5OU|RB`10vzkDLi=KJ!&>o~h9w2^#IWBP6;htc{Kyjiq<87-;>CJ`@yvV+laknqb0d+TE^m^)pjj_$Jr(=YtOFK|AeQB@uu9zL^f zZh80$O$T{2VlRYV2(1uCMaGR1;@dZJFfR%doEAvyBz|pBHyE4H?u7ObXlJ3V0Zm?Q zNezTuNUvf2*w{W+me+%H8E{ql-Y7>Vw1^FxvZ^r=>xd0;vBHRR$uhW7woh4acjZCU z*>B&1>8%mM%>(S_n>sb%7CkEsoY-$x3H@Mc2iqvDBa&+LMxdV2b_e4iSYgr1LL-{1 zb_YIc3g&@y6Z+}2nYNp8*(-k?_6YI_VfVKYaL9n{nXLsF>fDxNDPh7 zGJ1>R|K=v>3G7a>J8?(G&HwB(1=GL!SN{q#<3IhU|HL=G@r_UauimIOi<2XCquWvICc6%rJo9#6Dko#A=FB| zH1bntdo?(J+&MiC@)hhDv|5QRLU%_MwK=(P$&7#L?c6YDhxYE3tL!kQVf<**edKu( z{bt6y;C=R)H-40hI5paQ*>b3zbQIRIbJwms)Qxk?EG;5!Ve}KV-_dnP^@_~>aY{#R zIQM&M@EDGSJo|Azf_>=67|makTtL1J;UkDw$lW`)lpw62$r%tDBCTiIu~6ei3maV) zR1!MffGrzrhTVg`Fznp{?{;{n!b6+0r_T^fKl-CT3c&mCzfWuJGwbF){QCdS-YSpf z!u9f+dU-{ck5PLK7T(Yu@H8w@eXI4#+!tbudu5ylV{449FisbAi)iYY1*BDEe2l|uq*pI`D`P=Gp zu{Sid6!O@|vojl1cz`F%bzs&;tCx7xMf5Y$nX1cjw|Gxj7QVUO) zCmxrN*uw`j-LR375fK$n+ki|~z=%qn#bQ50c~Ht-r)`zKT+v0Gt}0j0B(=e@Bi9Sz zdJy+UIy$-d+=>P^#08gjTO4z<-bzdowSVb0^S?rwCR*|0ndvA{nooRH=YMV8B!b>REOzrGgXWrWoWdLA?KZ8z6X zMlP!rnpeo~AsaF|e@odxiLlli-d>UC49(61AJ@6a#)EfN29JUrmm zK8?o~$CVgRvRez)5|AT9U6Fo9^nqY637iQMowAh_T?{=Nx7{_XT(; z@E^XwVEXpMduXn-vC)rAUjvpQMkR;FBEBralewEns0qWAu2VH)#$j`>nShE@G|VD0 zoE~efqrEeRF9d8TMf7?5M8zq&U0u?oCYZUR#C;6h;EJ$%kilMYrV2r~m|GABuM z7I6!v?-TgF0>3xF*CV`Ky!bm-_;7XU?(2wc5lhKEIvI2-kunD|MwL266v~crOaw>m z8F_ei&5?-@b-cwd4bKj7J}&yEe@Hep2Njn2wk-ML=hy<_2!taMY9{uTXs+!^Bp0NL z`QD)jEh&8|^bPdX(Np(Y@l&qPEtp=zY%LdzGr|iZ2P;=jTsU)O;~-l_Oc^TFF(_?N zs!|T4>_*8(i^F}WIE_JT zV&vvyoJFS++Ps6rU&i>r zcwrl7+sIZ?GY*YB^-e1t-3KiVY8-Q@Xr7g1_-G{fi@#aEj7902QM?eK2d81s-Cwj$ zRJRS6d;2k{-u}nt+34nqbuN=?S-pc}lE9z)a-MAqpF6nwS#+NREZHpt`yCuF;c|~X zofwxuy&Bgt$jy+EsFuk@R02hmCkgBlDdIZ8W|5(Bu0sp>0#-O-EY3)w;xB6qcW>!) zt~rMe&p1?fPNLV5ahlEG<{Humu{Wd_mT~25nKQ;o5?OuT*=pdp3i~S~zh9~6;v>1YRGY++(sj5)%jx3AR-^ep zklzX>d+W^c?!jM67dQg)7?I0@Tu+4e@8Q3{gZ{_RM4HMl-?q$s0D4qJUz7qc&$ zjteaqLpo~iV`aS~!>gvCuBN?N+RH7F)H7Dp=jUt&l6W@6H5~-VV;=wFY{{+z(?-;T ztWTUZ@W9{&SI(MPhNJ%F2wWcn%VI3CQ^!twd_sTiifsk{k1sHo^qc<~=@B_Ha*4!? zvhIa*ed5KaysO4b%Vf?JTS1k*5Axh3j_sK3Zw*w(yl?Yydp3zQZ1bn@<^b1Lo;^cN zVRS>Q%>5q&lZMs9@ueWR%X;#=k&fvhHO;}}M0kUF;M_PgH%e0ho@9fSXGduQ__$x&_b-;IdOtgr5#v1fRCUa|v$ zWC{6CU=8w`J!@9T@`#ZoFZm-_d2cLPhy}0)dp8n-g=IY?(9?KkJl$QDAMtVD;zm}L zhkh`|-IhQdQmR`$JyF$_HzRJGd+s^yG_uS+tCCr3up|%9t~g7?84)f`+)J@wnz`4N zlg(0RhOE6-J;7vW{2rNuX)Y@y>qFD6R)0HF2Qm5Egdu<45~@5INGW|?nQHT#x%t_N zZ>BKKzAI1dy|}G4)%pnT9Q+f-eFpMakf#p)9MWZA92>}M;4bR#LVH431pH#4X`ayd zHBAnr3G!eFWg!`JfnxZOC?a*ca#ryv&TxHISG5^rpC6ivo*ECHa>QuOIAT3hC~J+J z(sL$(P1Ui=(QN6xIzoo!P{2FwTK%e{%QUaP=G&Pkw0g{T{AO)pCbwOrD zrlk}mGf-HVl3eXl6-`LVxLk}3zXlc+Z>)B2DjXcV92RVIttwR%SFHh&DcR9KTV<|7^rP>yhUUz39-bhq0lg zhEz^S{vIX#H`2XISaryH!Y!VV{mSUWLn8T+tRp2)q>>ob%uKbK@1$V%D$z?lB?_f@ zN_JLvOwKxZ6@wr4LimQ{JJj0hQ924MB5MZL>R4gcwP|LYsR^W&QIT!NH5nHY^0KwL ze|cOywj+44>UhbaqR@_B16%EIYDR@UcevCwJA|?9h=gkAJ|YZwop2?Ybn>LQfxXmn z?kcm;)rF(f7~1av*eR!A$DuI8?zVNJ(ayY+Sa66Q67GU^(+uu_7*-V)QZWPSog(MP zy9=Svltj0i=q`Gi=NI@3(91~fJJQ&a(+Vj|Z~>JyguA%zzwoP#DgLk4G(ice$AxpG zJV)bma=auakE^ue?2NLT4I&~nkb{wWG&ktfn1zOfN+VWssSl>ceY?D`IG<^E78WEH zI2H)4ntrDi>pAR-sSmrTv|edlA-ICSX}>v!=YWYp6%C}~8P&7ZhI1A;^$a+KB_pRi z=Ztgp^z165u}8VY=_O7*O%|FYgwen_T^Pk8E+MuD(%4=RuhKPe95r{&C}yp%X*~uV z`PY%H__?S5CSR+f>~lBU8f(jME|zjR#rr}VGmGuS>U>9c4&gGP{XiT$rnq1#kt{93 zP?1H@b{6bvk@bXxQ-~*+Y2gA5EyY_MM4XZT)Z2+w7m!zi&fFauBz49_*xB-(OVtr) z!eN>WGhNUENomP+2pxrX-=N)}6X~T;DdAjJRCX{L!L%1_GDb zA_$Qc=LkZ>#ELB&E@;?lU}P|HirnF7Xw}nY$07;cB!p4X4q^|(0*2a=TS=7UW9sJI zi8N8&E$dnNQ1Pu?j{2XIHx`-;wts1B`zK<#flJ29Zzfaxbs3Rj(uM`vTv*1)Iu5KR zB}@^GJyUAw%L37f2@kA6o|64tQdkokD6J#4Jt<6ND&?c*L^F9b*Q#l{+ zdV_)E66e>bZqU4hY~o*12GA|Wc2`=g0%6K}&~Y%C4?w8gilVhhtcIYh3Y-d@psdB( z(rz@ED26@QslxaH1Gy4hTLyAvF4(_~wDk4oQSHwXZDY1^NfC?gl9M($i8Qk$BaBB>Lq#RUg=SuQ>6@> zAh-^~LU1dGE&@e5%D6(y89JR&!p8VnE+MkT<7mjVL>45LWY(zRO8cOa>PV&a4p{># zSp9nA;LO1t2ODR7j#nW~*B4AWt4G?YAZ$tnsoDPm5=v=4F6p%RMxj#8+)oFRl>%K=IF*+!?r*LS)j<8S?)! z<;SG_oRWK_Sf#*>erm6r!}KwO5E_THhNfy+K{y5f3fZ5cuEE4J8XC#8Isl!&&VFj%3}3GaV7(`@yaH^Nq`e?rJAQV8ch4(nGd3M$;6`AGk2vF}!y=lL;#~N05kf(N!->Obz!$qhaRgPGl55gyI!$I`nkkZ)+(`BpB>xx< z8`L#qDWtL^md}al2`PP!<}JQVNJ)4t1W`IQAC$Hp{hQ)(IT$XImryo{o{{v7q75*a zf*o5Hw5;%KaBOg<5st&k-NF*72115qWhvyU!`U{EMpFtl4AQAOoI>hQvY(Ym`%zGlf^W3EUTtivc%vxUs{JjkRB=LMS6n zN3)G5QN#kNxJdB>$zKw~XC%Kg>oYNJyPS9Ad`ZcdNZz4kMAVFie9@YpSwkg*;DbFu zDoypCS>fWk67HJ1(zi4k=x{7*5FH275~0J_MBId3ebwR2E~BB~yFyrkTPwPTys(er zYY3)!%y_L?`%z^>nfgl8|$u*#24W{|vR_ZpV zb3JV)1yju0b0 z_xN&Iw{+XkizKo{qD)MFXT5s9A&IZzVwqAV>*3EMQbyct6E;%_lA0ghAUJOV6gt{x zE3w7(9d77wu_I#^)=?vV) zI%TdJHxNvRLe^~SSH+C6U3uPb#!3M@M|lohAci~{a9Q&!&_nTFs}wfYkT7yM#hX%Y z1vDEDMcJ7fH}ALemGuBrNj5iy)dZW-l7l*y6?pnLLuLlp>ldZ!^!CDYM9iBh>rS)dq1#6gvis8%#K?0bM9%qr|&{ zyQh@1dPceVV490X$NNf$?mZ%x*m5CQX1qf_RWc0c}f$9cE5QMY4i{$5qn_w zx?%MQEv@~6o?yQ`>-Egcv-#@`MwV!&wr*$GO_MA!Xd`r6kR61{N<0zV;PKlIf3c$3 zZ3z9EFgko3@wwp2fa(s9Cn#nws|Bwaa)Kz7vFdy8n*>a`L{dRhvFbpvKpKHQ&o{7w ziId6RbmH;Zt$opG|+Xh_U z;K#PE2#V&0qJg5WYI8zXrIcbT&}3yoZN#`%-avWXDGZqQWR<>2b!&I}0y`iBPUxIi zgN6cG$vG2CWGabRYFWmphI_}lgs&;b3gHap4v`H)Ln0Emq`}ilVc`VR)dGzQJ@N#`~e;-iCt5eItTYG-Vgmt_3Xyq9ujCt5o zxU9{a^McYH&9>mMYt3@ZBUT}Cb#D zdL4TJzt(LHdlIV+iY9_a8Xc{Aq!G~_x_4@b7)SvHkE129qUDUhJ%lqNO=W_1D1}BR z8XZ{a$XXL^5uBkOXgc)n{lj27Jw2uCI>vGQQhjq*kK@P+ zGn0F+jx3MjL)mx!>@%?$VwzuhKK;)WaKh0=ee6Ts=8_tPv7ZiLX0JvF%E`1jKNq+V zbeqsVQpS-S2U2oO*&{mQb);w{VN#0&N?(*+HZtGj)K;rp3)OS^XGZKt_EehL*yWg_{TryTi^QDKMj|;Z|?sea`gSY&&-L7L6_#4l8lZ^ zVMODq)lAQjjdGpCdKYyE*{p z&H`yjlxZTTu`2$W;kD|7R8L*A*~p9(237iE1cf-+U#|H8uXMkej2$v!Z3ql7$nvXvEh;Mx38(dsmyj0)Zc6+t>(;GbP2U5A>W%e)tM`ahS-G4nJA2Eh+5Rzpsn#KK{bXmY2O-pe>IgKnm=NK@=w7rq0*!jQ;#g(JItq)r_{G*fmaFc$sP3gb z{^jw@U;dI$KmC+9-gx7yPIcv&bF%M#_q#VT|JU&v|4H+2w_Cpb?QiqiYp;Fbcdq3$ zcmMu!%5*z!2V4zg-j3VxFO6&2=Kk=9KU|4HeU*oj<6U2Ny#N0DJbLtqX`1lfU;F$YfBZ4O{q1l0?6c2Utyb4Q zKYTeItuH&iEZzv$v@(70!3R8i_>fm$eU-#etV`t<3I%>TXbeGl(F+wJy7=6_K~>$>NE5pR5%`Ts4pIT88i zzw{-Cs`Bo;@A9J`{pd#K2l&ZPe!{(b_ikkV-~8q`{OM1B;?IBnGvj!*W58D#Z@&2^ zZ@>LEzx&FMdU&p(d`4<68U z-L=m@g%Ye`Pa&LoX-) zAUj`2khYz_x4W0WyQ@8Yprfy!tCy!Bk2nt>7rl$Wzqg+RFR#b{KfvSV>&z>Jt_BAh z!S+@&_XB}2S^v9Gw2jurK_Kd3btQSD;J0TbmhMKYUE@f%rq_k7()&HF&oOPV$Zgn6 zMGru%5(fIq@#u%u1lrYVxfXQ>mWw7HuKF(6NjXBiu3H?&uWTUg28@ep`Y(*-!gWv? zc0eM5Eq39em<0WNtFJ9jq{ITf!>KmH=1VCq%4}?V`&Nqj2KsJ025y>=sxUArEp0K@ z&s&g`)R$D}YuaM`d4#kGP3dCI;NT!A0{L&d=Hc&B@6{A%TFNB!QO%Qp-muH>82@GD zW5DTFA611OR(pGU@HK)k;xgH0@Z!(@#`cp4KDAiXbUyEX$G{^u*Yq%zS8sI3b_CVq z<6~`Oo-gxN)Lj!AF6IL*B@tvZrg-hr4mvQsn7HBvLZseyMD4KJrO$E9^oV zBB4gv@mm_pf;CKSMz%Vsh zkzIguqUmzO;W9Ie-`ARC*Tsd{V`xU#^d*{#2!0uMsy$3ZZc;bJSGIDDibc(vId~A z=!n1B>ROlVpT;+QxI>=jwNRg^GS@r*;>^fpm?lr(joe4#CaQ^*22peTeYL^_v;ELb z7M3IZi_m6$j0rYaAdabkRd2SObs3I>lc<)F@K|F_CehvQc~?g_KQSMbHPgK-nAj*3 zih~o&mXUz|=JVKLj$MJ=%C;b_vYUStU0@E}%OJCZHket{g8 zey_*6Dj#j{_ODCHgtD*f-iPag91F_<(j!d`YU5FEH+qhW8q4M7t@HOvo>D1!9hL4< ze(Y5Ihv@y~L^#Z;odV<2U}i*HCeS^IlQNG=O=&nHKo-<}pFRg|-yN3CXxb>Ykio}#g?~Z1sTlcr^N2%pMY$eO|c|EmZ zrSI&Ij(0X}YrCQ(;--nE$KakH=q+g|FnJ66yf7s)%>d7FL7L*2xkf|6Kf5hBU+>yK zKrm;AgwkfrsYUNOOkzE*L{2k;W+F*<31LlqV@*EZjh+)|V<4W-ycP2Z%bPWKc?c1X zIRlQl-8!ST-!($IlVku!=;!ba%4Ik72E(wxmAbG^ECL?Gm7sXK@MfOF3KotQ& z7|{O7Vd!0gm-cybwiu05@-M7+Qb`nUm}dKXEnrlF#%J{OwGO%Ff2t--+(qw+Og_@n zmYoQ(J!fP=C2TKVY)3g zZh)_~g>{Ut34SaaTh)17l3)pQYJWomuJhd6Axp2m*r4s#WFYUoQEa}kUqjF6kVeAnJC>=E-4mnc)@y4Hr%*x<$O-`k1|()XLfU@)H@u!Dwl9BbUfr3>e> zOslJ_F<5m=N56h0L5&8ze0^ga4zw4HeEj2cdupDOmuGKp57~sy;y6GmO`6-!RDf0Y zbnCcW3+?^c%Av2XhO0lCU&q+?QoV-mJFh5jww%4w8bUX4J0pHkzI^LXrdDV+L- z4gi$oiB91yK5JIw;SIXbV2K2KIZnJFYY4zSh754?|YXYSfRYP zS85Zwnd6=dhp=a@v3O|y2xiS*3%#ZVaASd6y7C&fGL@xbV`wONF&$X~hl{8_J3#lL zRR!R|z@+==I)OaB!t_EJqhZh_q(PCFLz++R{y@xsPL|`C$j~?gh zuxdE@NWuTMGz_plKBN1rbpg~bxxROYcpM8A5`4>Phk7yAVS4Zj_e{L8tmaXox{*Zk zfQHSExq&6YRbfQ>krv050j^=0FZ2^#w_Nc=4uxH=VZ9r{GF9`pwlV-$hVz}WSS>9~ z!$`}Fm@l@lE#%rSm_d{5(opvHpMfTS!6%l*nm-J3XaHyowyZrSB#Y5nrW3OIn6(&i z$c|!-h=hoIs7cAS1D5UBFw5`X#vz#b6FNXx?xjxxV^&(pK5Ieb9Qv!qjc>rHuEYgO zDv^rXE_Vw5P~}OLnv8ntmOw`J?~JWOD)E`WuAHhG76jYM$VFiS6Lxo7^{`4``Y;N> zR754+pX;!e9V=m(^lA#xaUaGIw2t--apvsOMZwqDobr`l#2P$%VFV{a3Rtidu+tr0 zQ4Y(U%YFR|CscrNDL56oj7kk9TLsn;j*v=y5C5h}T+;j02%uOktDB1F5$0gzA(9tV zSUugQ4?MJ~0wF3|c2WxY3Samc#UB-aq$pIU?P}<KCjsfHU1kj>7<8`7o?Y zZc>YwGA@OIq-`&x@K+7eX^nggdUDxqYAh|HqeA&e{Dl6VPhTrpFz`K4FdOS0W-WRO znn$=~DaHGpm7DU|oJ0D^MGRae7+9!b7I^7>Ns&AsClI_=ol(urJcw@mF((_^{H0zW zR7LX6yiNd%_iB+YV{uas5NLca+Izm;8k@ty{Er4kn#*et4|9Pg+qMr1&Zx=FQ^~yj z2}?D)P+wmG=GkKLR@Q%*N`&PoBnxj+I5QrE{BYd2f4zndDd?ea%M}lJ$lHmh&)8=! z;xc$YJ%sc=TunE3lBV!cU(pN$Lli!a!ea2{C_UZ9Oa5N$cRaT{6!Q%6^QQt6k-72j zG3?fwo2qt;;bHDG(kw_K+RtchvX|FtCCVz1x2NA63FaJTsfbxv8nF=yss>EEF@Q$b z%+6daz9}bq5g#6Ji-jg%WXVh~eah^dLR*mf)`632z^xn%!~zrI1OAgsIE8{#IMn0VE4hEP&u7 zB=hNrhe|e60h>LJ@8c`zM=)wfR>wlA%Q)Kap*P#od~`d%4~`LEB@EM;2*GMfZM!(UF4~RSilTS zaEMs|y<^l`Z@h`7!G`}Ae~5_G=-xP2Uehdi*jaSeYxGVfRX?Q=s+ZW!CT!0x&<=WQ z-6^v)u7u}9J6}J}#(>WK3ipZiB2sUNgYGC3u(cyT+D=c|Ejl_B)`X)``y$ zE8U_(TcOposT%u5IJ1NwuaB;*25$=c^<9JS2vtN!cfnf= zW1wr~DlpYwjmU2~*;DA~M3Cqzwc#j!WC+lcZ(X?IsiiinxSvSbFV@bY&LmWeD6jke z9{`#$;)N~uB_tt=JRMZX@&Dliij&I>7ztear^e0*1(`{4+r6vr1MAnGTyk1Og z7qN5M{`09*O+o>Mghca_F9it|+Xx#UI1(^^aR7Pm^!?cyiAH;-XJL^99Z%b(L*6}G zUtjk=o}+nqIF=p%&kR01goj^*ZI7I9_8cRhvUO#KioN<1$Hr88e+}bK-)0L6N;9CN zP?}|Iluo7;RRZ;uc`%7lW~z}Zb|C4{aM^^8rpV8!LdgA5n#G29m9cE*RnAeidDBcT z;`DS>H{y!;E?mZ?DDLOaS9o}M<>lp|ZRFD*26ZBtdjumoDvFV@am7Gnh41fD>9Ai7 zHjg-GXJ?#1w<=}Y-a`O-@AgN~5$E;X1XI`nP(SMTbnITn$5DmTtM?1F^hP?F4CkfU za*da(>yYF{FPimgCScZFgp%^|^8U_hfxI8?d`9i?;!FQbCwtgjLG999ybjIYP;eEKXnUY&+9SQE&yZ+p99`*;?3wHe{f^!r*6^i1ppx? ztm2|TD24O9t$uYSL`f{w*oNUamQ2U>mhnw|_t!FCr`nU1Dwx#%+mOeXFF50atI9r6 zzGmnAtdG>MRVvXX7Tb{o;pUz3@A=kIA9}+Z=M{f5&=wI?fE_azi52v?jV0HSzN5i- zBS(B3=GAB2@`SWQDZ<(C0CZL}hig&K?}OEv`C$C`%zRPKq1LRmrMKtWZYNXgKL$_n z9s(*DRo1eoSOHM76O^*7mQYj-4e+HQc_(%+WKN(s20s58M3~7aHiN&h(#)=9rD-(+ zKylhJO-IwFf7?WjZ1697a%jS5d-g9`2J|_votW6?PBB|C2GK5r7ymOAT zfG1;UwNc`$I!fWJzG7djeyf|gtEpYh{)tkF031h0M|Ka0nD2p_M%#l0J_HQDJPWv5A z2Z>0}ss0RMj0;W)CX3`wD7`erkV!vfb*>T9c0nwMEjnK_Fc3ATLcdf0hjO}|9dfzs z2d7#Wj50zTsbFxk3KeDw3P8SW|(7D zr0hs?g@BnvGPk6G4#5Z9fWNG9i0A9OVxg6}l`%39iD{%ewfEXds5LCkAC7978v<4) zjq-7o9~RA)4v|B1O_j-M71;Pnn->#Zp;51pGJX`qA)d%S)0nnXwOmsJvrqr{EbOba ztIg-$hHd4rcQg!)IeAYoa6zjR_)vfRpMc#k7!y@!rxVX6W!zgZWNA_6Ez242#}@_| zft)=RL;*K39!ks{By6rHq7ca4z$kK8V-Lr_LUMi8s(H3WyjM4XTgGDLo-}9b&f$K;|SW~@_Eg1}qZWg#Y zFl$k{&EkPZ7XEE-OlJt=88b=D7(L){7VF`F@wL8MTTYP%!kyn@8s>FmPFjrSZr6VI zKYmngBgR|!fuCsoU9GVw7|>n4%`0sZ#TQrD-Ecgf;Ged8iwhft?`sM_zAwk6Ic)sC zC{={*hx(r+{wrj$Y=Z7$yZOY8BG1nmhfN&f;8Cq8EuykQ1t5u)?n_lb56nc)EH;H` zzgl{OL$h&scH)@~enVP~LH*s&j`7f?wSP&d{YVPXcJ<50(s zTL=m}z%9o_{W5@g0A+}!rc*ZkDMHY&fLtKXJ^zrNse|{+5)gPavwo0-co6Lpeqvqs z(d^0f&Fpj^717*%105@h9?Q3B!nL#mEIY8Ae$hMGsy|aGD}mR`XTSErrtCmnYFC$U zw#It-E{JmDd$GD35ddUB5lNkd1x_j#gPz+af5K4yG;wJ>WmT=y9HZsno_)HlVLqjC z`SwEmu9mmvL%N3wFCjK;CP|b#UgK-XO_WuWl{49N(=DbsSlj*P1>!~yTe-^!h<5e? zI0YO75N-ILTX4HxBLQOb`FnNv_xTRzkF6QNgnH9ofSA;#ORF*n71CE?B&z2!A5yW> z2T53xTrwHSG{_{_5Dyv+U{neJi42dhnw>ZMK)x}mgcp|E_dCGx7%osV0s{OuV4Z+3 z0slRxTtI4FH#etqaCjK)SDJM``*w{u3S-JqbC*%(u6922GdaM}fQ(2drehptgm(au z7LfW~x?UL;45Nd9)S#2RCn0Q1dgH#;)(&vp1^jWf1=^Y-v5A1V_IIJ-gSDCGfb{03 z{;M@**U3p?n6ls_Ks3%Tg2KY|4WS>*^iqvtfILFrbUQ;^USAHE9`{=&Is)PQSI`rp=P5V`-X!!D-}x?j$5cK`YjE1jdP#KC|VHPx%gv zy1ydEWmjQQA1EzV%)T1d$K^(IJYQyWZC;;C`+Yb~wH4eJn=-&dsDLG(!S#KOn~yWsBh#r~R~R`c7%3xzefw<8uBKdmLxmsdE% zKmw!?>autN1B1ak6ZkP5zr#I9sptC&Tb+KX?-3s?g-abrdKyRm>-V7ph{WP+-t%1J zvT4(@78YRgmJ22@RUqb*XXb!US(2UBazGt}WiCP{M*I1yvahKg0ZR?!W;L%P<>7mG z#?16>UHt5U!miqT*pVN1TNsDj%ScGjaiJS1l!Kljq5y7@ov|~+=mB_b=CCCLbSNS1 zCCS&7sC~iL85oz_Pe>i}Gq|#ce5zHOT@$$i*6qJrMY4NHtYD&=ixJi8*HSPlDHmT$ zibciu%2WAVqr*=49JQ6VQWG0a!E5%FA1vJSUa~heQhrB^Qa}U*hk=*{8q`{78B6$B zy}5lnWQs3KeF(_*$@o1}y@IU*_i7?G=1e(&jPK6spA`3Napfb)CLVtKG}{^K;tw!j zBaenq`_Kl3IGjU>zpcCr_RQF<49dXlLW$TbP8_>(h3n&3}qlvQI|-x@9x@GkKz&x1zuf!*wBz@j{1L$#Vo66HA=HOWu63Z(|^gs zq3dD@E390nS+jV{@qKoHAvLrANSuPa{!Eg-Qt}rbiaCt&*#qtHqC1j(C#TO3qPXV5 zOB)!Bc>k6k!C=>I3MDzcDGj?}EjMP;o&rtp=_49rb5Y<4X8&beQC4#S5ABT=UIz43 zjY>|>TpGYJRWg3tz2y@begF1{q2mPka0Wj*ki^+5!`q8zeIpFRIZjNg1*Z_lSKDYp z8HDNc5?!7n3Qb)qZC_HJ{Y2}CZ)FgeLR`26@j2jE7<eH8n zb8KZV8E_7GBR{4Bs&T65$|R{R3fTTO$a_kM>2Z6P~os!{(5B)eGv7mG6 zt&n3HEiJ8?hO%$i%tB{Aev;*|J(n&0j`-b7NOHWRrHR+nnk4dN-bK@BZdQjfqGv$6ORD+*;Gx%I_Z(#3;L6_&KkO zY6_ohs!7K8hve7Q)yn@nh)KKN7fN=zCuY^z8BHYCr0spXW*86Yg zd0S&+S-fg>nBL%hAoApv8?1Tlk%m!a4L|28kY0UAmn>Re}|_y>`cl1c$qaJBZs zhj^~(@|#~z zx;iGH=iwCMid5e5S7{_g}%2WLwl zFN%B};v$%i^13dsUZFAnH(o1i&h@`~C9tj_Eura|b=2zRNIC0w@*vbYEFSEKj)70x z5d2U+T`i+fSF2N9Hkasd$TI6MP_WKGc2QuXtwNc_KVJiwp3T*uGQ= zbV8iWvHB;0-uKAH$h0?#I7{3w0}f~M_V&Kr`1owvVhHqD_F~TDw_EQeG1>$^uLB(^ z!`q)0+qW`T9+@mhv+|mnGCsasd~9E+=XPL7Qg}lt$YtJcuG^r`>3+O~*`c8-CLzdrx<1;gTvu)RW z>t}-9-QDs4OtfIpyUir9%|2>xATqk(<+Ah}A6OiD^BOUNb6!S6utaHTzf8(z<%Q>$({WZ~kwsQsq= zyw0f{<@67I+-aGE-}mKZLstWum4NIt`Hxv+h#TS+5SW^smDA85#1p@z4LD~Gkfg<5 zZF@LI{?Xc+J?3}pa_cL;e{isOaIkUm@_e#)c-Yb;@cEOVylqzW(SSTmT^`W55i5cb zD!K8{>sG*cd3AL*;^3cDUd}0gH=tVWb}&;S1}g~?^SQ=9ZcVk*TiSTIZcSj&IWt4E z%y=^^MEKdza1`8jo?IdpV}ryW^%uq7EJ>sSfWd9h}L&wE|-2) zM|Z0Z8Tk7npdovb)oqCieYjh&F5*RfDKIX3;E9iI6 z8z~)hFB{62kaj*Ege9ru-=14nSC_T9ub&59zg2{5jY4fYm_-LODY*Aw{>tstrAmu?9LTkP#4{%D~hB|exkv_a> z4We&_LItY)?KXTJs9;11I#~_L`7JCFb)RTteSAdtgNVS#157tN*Rd1nwu@^Yjl>4z z1p^^gOKuO|A^6byAYbeG<|;e%0)_qqC;E`p6!W1ub)&ekqZ}T(*E;Tr-N7IGw7F-( z#G{K*@0JxG$P{#Z5BqR6lZnhJPp;E6!QrIj(P~<6Gc5i-##gdGZ-qWd6DUwGwWSc_ zPNR|uP}nlc0}`sT*_OF`3J=I_v^8B0rAgSnRER3kp@^rVGbs=nOkv#lUbK~)_eTu* ze|lSk?zjf`VF=TVT2A;wWa6C8MAKx69qa0tZgDMf5u%Q@==Av$NSFn-OSxy+5mjAX z_%790MaX|xCFLX>%EZKBDgtgXjZm24wV$iXqSDhFck?q`VxqbIC&TOGzzj2idQ(jB z1^e)qe`QGcQ49ZYAQYBo`dQGg@4)>cQjhbUFV%218n30*1i2}`^+?1 zUH|ZKV{2BW#@>yi1U$2IRjR(LL%5wYxFL@CaVQW@Gt%!O=E_sWK)w8{um#{w{Bz&RF%SWuY8_);6$V5Fotri*Wgw`5Q?vJAooPHLAv zO;%wx6%kB@)p^TyVa}8lk2rcLWr0N@##ZETh}qKjn3T?#L(!#rY>*PknCrzlp+gY- zEz6aL9(Xieu~qb?PhlaWWWLRPFB-@QfkiM75;(O}bi#F6em*tgbi8yex0I)6m+EUTOZl8z3F}(F3m&^1{mAH)NK=iQ7=%edv z_ity@_k^mVoeA`Ef-c)QFZUHM!v2RrjV?b?-~yMqxxh=U3f{mSZ*&7MSN#Fa6^t zP0mhZ|Bbodo2NC>Na)BXUWpIMgSVMwuI%Q~FSdMF^zP07F0ZduxYdCymaH7W^}Lr0 zcmLai;7c4iuLlzvfife|q_fVVekeC`eFBI(O9_c9Q9kwr3^HEi4>LA<+E0wMU$0tm z6xJTxAM{vrA|AC$14>HJ&^(-8J8`4$-7%X~{5UjP*Wugk{3tz(U9(RiWyX^8>>-$v ztONe5<(!(F|L1^icZevdZj#xF$uA1`>~FD#@6&J}H*^0qRuiyVMdc{#;+$kRPS4`)FMnOib!11L_irG+8k#Uk&uJ7@;eq!gJ@7EIini{zXsGf3o4 z!AC0Q=ugPD`a;O~NNz4yxd~Ca$jc)3A5wW+&mKg5ayt}1NtZ4dhflqfvo8ds4#D^I z6-cU1$*0(|`mz=pJiB8&-3N0K7J%`VVrD_|M<6(OVwaz1T(G^F&kAK$zNYjM2qXKG zo<;@Ot6!HV@;loNv7z{e*42~n&H-0kqG!(p_F_MXqzI=1QpQ{a@4Wee*?BxN)Zhz= z=nMO>Z_yAQD#7mFfLF>nKRTViDCU6^pk^iV6Fl}POaK|OF6Q< zg%hUgba(q|h`O5FU>>ueEG?w!{hK*p6|rSetN+ z&VB5@C}hukpl*x|OUY;%jni}TpT1&6%=wD8*KtDAm0*lp!l8!SL#{18x*i;98*;S8 zSwpq3+49ovqYo2TY&+74YW~wNtXQGoi^>+)U+1FR~N&iPCLW4o4(=@0E znO9hP<{%5<{&!}2>LsS}e9W&;L?L0K5@Z67mp6%T>sf^Qe8Kdbq&)Ww5Qxhbj|T8Q zd;W_FX_=CfCc5)`>^xdj`}FnFThi(lsX5`C-__5Wa)z{eKHa&OZ)f)k-B0v3{pKs; zWoV`j`*$&DdJgHV=|z3|yba)dL?OdHi+{jO#tEMLzOfW#Je-Z(ghO&OyAlX}12=YV zcyW@=!p4D(R*zo2%7Ab1;kX)@bd;`HOJX3k6|>5EJAJ5mJFAX{39+!uCc^BIy%}oTZx5=dY#ekYy@XTUZ1VO z7T>WCy9%@SB`$6l1i1$Amwg}J@RmAzYk24v;|#_WVH3!=5$5O>*l1e%qftT42yv^j zh42xZONbg1E~cb~#cmYJ%eEu`B3$&RdG1g~ZFG=>kmZjINGbu_Zr`KhQEv_~=Wd3Q z7?Q5!Ca2Am8w~l~{Oey`rS1G0FGu4C1^TXGzw*5sGTL5TdFF3zZPnD&R1yLc27x@{ zzn{^u>NhTWEJN@mi9#r{V)Qo?rxoj`$HEe|P?r2%m5Z8Ma()epvk&D}i-xS6tZmTD zM;-|hk)Xy;qLp6K4rW|j(MPga(vw6oO!4#SSr=}G$SVr*S+CCMIhc=R@JSa;L|0k( zk=MNK-ef?cL-F`q@1$YQ=AXLr3qMIc%4||nn4ki9BP0)q{DQd@Cd5eOo>$5Hd93l4 zkyRh&Il?yB*nhMh*j)e0TJ({Lv?brwCL?g@!KmtpTK~9DNX!yO1iWb--IHcEcp~{( z#7i|+Xc1@WeOT2GZ_uhHmm>}fml6S#Gag4WT0RcH?OY6_7%OG=Qh{_HXO+(KZCspAL;riOPLE3lxY28g+ zCyysPU8)tHGJ;V2h9#MXr8_cn=KPP8@L_S8DuH2Y`m!Qn_5KvXo~$f#zE^t$Qp(E8 zWDnKiWnp1qP6m(62K08KVESpd)DYo_5XmHoJ?F~HHrTD z%Fr3r$BNV4LEY*) z&~IuQ1AicOA+)9P1?A=P0xaZzd=8D*>vy`{!n-~l#=AK`p9=r*oO9$=Zlo>cl@J$G z2CmlRRqoz8SJf-}8qE{)R8DR=NbJPc>CklZ)e7z8*HWy_+e<<^6x=`#QZ)fOFjXI~ z)2yUXta|^aQt61UQdyTsdL&26fhw2QkPQWutP@iiQ+SddMj&&%zqa8E@s09_Cvntu zx&X~W3RWKt%Mr#g{=7c=4^M!XU~P3Vo895(54g!TS9Ih1Y(D+8%H-d3zQ~*YP|=`` z-+7Sv)+f_}Vh-G7c0F~zt!1b>gZeuH8C`RyzLg+7>!T5WC_6kFcKST*{}<$*5Az& zIm#I#mo%pTfUXyt_IqA93sBhny8X1c&|>sAanMA?{Y_zK{wc{iKP29)}<1cja;YpH+d=~dX|HQ-w)$Q&*hd( z76`x@zGA}QmeAMn1|f5R6-6%J>hTO-{jihcnZBAgkZYtp!q53P23|9em^`z-9u0%qU0lR=+p>zQ);CAhtv3xzYcqa`1Yld>bm~jdwbn1%IMSK^yI(0dd0XF_Z(D=@JA*Q@9^+&uhl^G z-MziMiV6qc7e10Y?A0?nhV;_{2#nO_kK+JpQ+Z9(B`4+H*IlM!=t!LZ6B9l6cRc}i zb4$I87DHQ3CL9>izXXX=B!^FhwfND*Hb3Vf#}lPxW(jN)meNGp(hCj)}GF!1LBH9H_-rhbFPe z_u)5c3$fj*oQIZN#Fn(g8>3-vs5?qYBxD)zDlY}8zP8AC5Un|UE;G7Juik>AB|w>A z%7nI+Rvy2K0Qw-A7S<`%Ebuml8Gvan8xcm&t&(#LhbRa!_CjRDa^%lv<~gkgr8JWb zsvZ#kT(=cKY|C`-0FR&VW5HCvOpFXHc20ryUuS9mF#12S;6uL+oD%Tg%F;WKUr0KJ z_Mn3BRIUB%F*FzveI%)0Bo}UHHENuzk`Vuk`@PnE8Ox1keN~k;z%Wiks*@&j_Tt4W zx?c50P2*p$-0?~*9Me$5IRzZ2PrK`OR{Y>@VYb6Rlj^*$L=%I{$TLGT>z0j9ULRjs z<|PBOAS?y}uG^f#2#dV;j=op^u6obhs;<1QwEu+{cNh&J<(}tTFZI-xwU!5!YuQLN z$+SSm1HYXe0G6^;MWsDsfU?|bf)`I#Bgmg{rXPZl7;m>O`B;(-4KzP=w!9KZCA!4EQ6)YGl6d)+=Jrb20> zw{P_e(7nOmUy`XL+ON?uiZUeu8gMyc-lL8rs8{#Bge#Nfm)631b?MR zp>kjQp<~<2J||l)xet`c$jB(c5^blUz*9sEeDPO9^{f}@pHN->!Q00N!G9IqGJrVU zt7|ZkTeFK8+mb<;(dJN^0_`7@cJRyv2(rZ zGG-f`5z#=rJ(71qpFGkuF2{1dm4ijC@Wzqhy!}=0mm?3BtJx9bK5TRVmF&7$S`$fb zGYgvp23>=b7LUjmiXvS2&_h#fnCtK8BYp1Bw6@1Bpv&F4xv;P>MtFCP1)SZi9Xh`% zPa((nI+P=+hWg{%Ca#H4>|)AjGRjIYcIyNWWt?~HIG9SSx0+6tKt1+FVN4F)6t zwzs>}+irV-Dn6%5=>RaYzX!KXFioEtzXTX?Hs_^OLtx7L!v0U1IaHmX`_EzCx=Z++ z=Hr>`*4x?j5HgY2?!7kqE``jxg#}%rlcnA|b*y5d1@?EBuJVpZ1}+84Me(YW>F0oE zuaVO2Z6xI>maC#>RB1VC^hO%y=`i4}8lpp~BV;IeL6LS^${l8LWDMyz$jOd?stZmk zh}beuYewP7o!u@8$q+M?qJ%#L$iNO)>n*`+JO-tu>EbT%u9axT2Nw8oXEqGGv*)hW zb&Y{4SZt})j8PFZ2p8j|{ufuBi)6>qwu)Y)3&HxcDtm%I)#1v`=hAJ}avgGwGZ`7U zkIFE)hw5BmClt0Ooo3|TF8%4oS?)@?L5Tzeh@}?n{;t}hpe2P#2x(x!(;%JeWiksz@dNjQ!7=Uqt>vGbRmfvr^(|bIZpmI? zzzm+#w5o;-0)*wOp}kLt-fGS!mj`!g^~=LXM;?UXd%ZhGVGxz!{@GJ>{?)UiA3*o8 zwm*$6JSHh9*1AUy#x_4=R)6XcNmlOlZT)rl#kA*fR%1lTl#}N5di8_;l6qTq3`_w3 zwD2FF!5ZE3NjLCf_WYaU77@bdg{H`M;9_=ZXt>a-u^0?)I82FyJOs{_Pq-3X_9d-t$2 zVRJJ~pPiIkp8vhKlngi&Wj1+vZDq?Z9vC3`mJ6`e@N?w)n;0ja0Jw;R+B6_leKLJ+ zthx+69({0^qEvCVH0yPOd=|~iBW?(j=P9zDWAm$n*+W*U&E#jaL=~)$^D?<`h-hVL zOmiv2>#^hpjw#;T)uQFlWTiE{?RS)di`@R+<+_8O=E#6zwpA0FP!N@#!b(-7ymc9l#3it=s&T`8o;1R*E( zWk9`mqk)h}$OrwfOzT;1*(??ED0^B#DZc@3O)#;`#i;qEJ6_FAYIt$6n8Knv)xVzw zRi1l6L=|Zo>eBB}MZxTX+>i`eDahU-dhEIx7^q6)Poa$QTq!6GvAW;#QUAh1lmc2C zUX{)m#lb)InX77y49sI0O4+s-o8`?mXRQALm5x}7=-|Ssv$uv{*b00MxMTRY6J=%f zM0&-+f$GqM3@P{fv4az=9ZgY^y5Y8Em8ON1(3!sc8XmaUatZSX8cucW@Smw@==tC% z=XoC!39JdU-0I)S0rs0l%a6)FJp+^><-#>W9Zr;q#YEJpHu}tLC7kd)KIH~2Tt$+Q z_*Y#6r0s{|QPOH_nALgJfO2Ok*_Vd_MTEhVVv zI*+U~ak&P=dT}o@hV-9w6Qr{)=!1EQlP#{8L4dosf_}lj)7{7BBeo78u*ZeB*wRI_ z05Y~->VQoHitrUc|i(&51T!=N)|Z(A;-oFKJ($ z9YF;Fx*{>?Bjsm=wvqG7bvi?p94xz41MW4%FVe@%4)u6PiDF4TabpU%x=w=V*ksi~ zXzJ9y$0%oMsk3-F`M(9FZ7zoU25ni+y)eE1k!H{eewPzBV$gWk-ja5mwW76>q9&>} zI~4yO!*jg0OUpPNGif?#?G_`>$p77E(a$dEFt5^OIBodPK!nreE571cu*1zx&W##0 zG#_yKgn}$nT5?{~wYBv$*;fT#bHRE4%IQhc+_}IuYfqs44uQk3j~3;?V3O(Sy@pbNs^81=Q_E==Ho{9; zyMX;e3H_sT_Tabln@xFM62sM%Uz7LjHb<((t&t|KV9mv*)W{1f**jDE(a*D`5xh(y zxti1aM+mHWA7yoxOpmYW$@4u4ldg|m1i@V1Cn|A4Bh(1lsSUJ7Xb{R#`eq`@+Hzc6 zs`%!H3n%G=1EvSKjm_Y2D7Cr%_2fX|M=z6)4y|rcZfuVmCiXaz+>F+%HMvjhul)S{ zHbPKzUEw2{0tE-c6IQ=h=K)E6P9T8n`ri2WnLcOn=LD)!Be7^mR=0wfs^Ru58G`w6 zk*`^eF(wv2K16#O*o!Yz2gm)k9UI^<7U6$!wxxNwoH;mb_O)O39+rkSq?lmQFj@>X zR3j__w>;~8sk;*NCit3DW<*s4D1A!+x9xtq!ot=v89h$?p2g^XGWDHrWBGwF&BAiu z-LpiB7NSBmM;2#(9xbuD=Re0J%f_tRICWrtGoBjlq^L#ToO!N#I$RP#$|d`;`dbJ? z$Pd-*S;(t~yv_N{zxB{3eG6lhlULeNCu_NNbQ^<420YBa32MeQ`g_U14vkCBIMXTDF zu``6=6Wu)$-NB9i#iWzJD`3(PXJN7QX6|&Nd~@xZ=xST1EpQ6j(amZ0jgmEb8&^7? zl6UK^&na3SsxYq&npotEaHV9m8L}Y$h?sLZLPH5ThR0tfK>upZ7z2tLHLyb4bV72o zMWT=fm&6HFnk+0N7l=S!W2nJT^rFvk?Gg&T=r4<9H8UN*Z8Kj)J#r$<&+bE2tiOM^ z+*Cnm-p_)F=h42}rlcISZDzj{CBqT#y=Fp_+7zx+52v#sIr4q$*a&d(siIV9oK@9u z0w{M9SU=`PQ`U!u!6RL-jC!N(_WD;p;MiOJEGj5U9dYAnj)$s0TCeC)ica?~z>=5z zOc%SNdOeUVL0iY5nv+!9g^#vLC{-vxysnW5 z$@)d!Oqe0bFM2wAy^+BTwV~~V_dOoTH-7yp*p<|N^vd~W#JoG-iK%=1uH!+#1a>Hd z-VTwJ>W-leQ@j}KvU{&zJSW0aQI`n#+e9>A-TNK+;ns1LGn!Hx3IZ$(y5(3ds;W_6 zo((M4jr^PvjO~6AJ!5i+!*LAQu%6K}k8{aL!}zKT7^QCb#%3|$_bYW~D>N=<|Coj@ zQU0~JO0mk6u^O>sH`h65mFU|u&8XNm?;s0YkwHC6(V3%vA;}m&_RLB9hVOOV{805q z_$VW&3YN68)G+P&qg!wF8`9|+NlcB^L6c-u7E~HTxAFj?4~{{j(s&|6{P?cuI<2rs z7v?4_18_U-y}FAZstq|?2j-3$dD^b?LuUJ3=GoJCNhL^n)~8}3xRAhGbjQN>%p`gO z1FNO^_^zzZR}fOp=#GV~zG=6Zk!c&<>@Q)9j^+^rFd~2m!lkkYGYKHNX}5U-Yeh*Y zny>BwuF~cQxbbPiHvecePBn~G#8{`U*QPtX3aWf2AAgig(Xcjpb=cC1YZrz#ei7Jz zWODIK79lOGo#6v&$e98T+lZ`G#$=Sgnk&Pr1>3y`!@A{;sEVQ!*|gzDyWZqm ztfc~#Dc(nWnlSr9;&adF<5R9RPa4*Hj;vdOh@NJ6mL(}#jMf&)lLn71NspBXuMc{MjU@6D8yEEN1IT6dm?8G2%Rqd)GtpNt{#`7R91u`D7|fLf|XHW$V8C5h?g|2)>f5{7Ah(KU(jr- zap_%X8p6=`yK9X|{|Zy5+&`u>UZgxYg0=mw9M!BC8Y$AYH|Ki0aSB{nKOK748J{{M zT5^uCiR=U_WI^c~O5Opk@!Fu(+$`y|DMHJ&w6Mw*>ew$~1={6Q> zJ*2uBmac&uY^(c?_dOj7aVP#(Anmq$c!HwwmzkNou|OWPgYj&{*Eu#+)oEfM4&>|E zUw`L1LMuL_&5=I%u#nUrPWkIc3n67|EdAaQMc)QD0BgjSJ_8rQ+iIOK^(ZWV@@T%o zptJg7IMn~#KwyyDc=wgAyjFBveo#m$@@=O(K|?e?-m$)%y7~qz zr_)%KKzvYi{XCKup0t-q$RjX4_^C4|c0*8_Ab8IG8n}BoZ|{G4D^Bg>yYf4<%2q_- z!$cgkW~I|4)V{D}{vG`8!j|Rg9pMM6_tqg76Rod89>ujeiaR%QGV{QI-b`Tj>SHAJ z_47VPtk^{j)l_<4fR`sSR$n@&93FkRVZpHez^}w+1dN=l7Mz?4K2qOHwn(`XidBAp z{t_r0m)%5S^!0&{c&COGU$Outi3yyR#;V_3sil{(dqQp83r7p#cX{}*rUb$p3+5q? zPNO8w>y1uuRzZB{IXUdm6V=aq5(-@BE$5pG!iml6*RNEb2BC>R!d^$JH>Tq)L()4; zKZ9a2`I!7&)md}I=D)#_G*N^9l6*~I4f#=3n{q#<>EyIYi;nIuU2tPD_&Bc>Vk2pJ z9;vNzW4Qg-k*t2tI0b3|pkq(ojO|O}4MnT;dI!iXlGfO7id=gVqiiNj^|A1YR*$+c zcsl&(l(Wx7Bcb|DYl{zc2+4^ueAz1P+fYov@2-JyN0e{t`FC>GIr~Rjk){q+zgmKt zwaY$D@v7Iw*om+M?N%KT^h>K&V>Hvo6^@UD;<|{g@Q1F6Ae`j_IN)D7%!^CA1K9Sz zpUOR7{|=1p_a1}p^QL4^w)hG?_^?eWr?rB0T$=<=$O$Gc1g(J}Glx(F^w}szu>dyl z&7`!+T|0AA=M|MisFHGi8FU8D@e#hcay@S!rGgD_ZB5xVM@0ABJ*yfbzeePX+y_BY zHFzJom*TE3KS1Rtm;ZJ9qyU4VBdR%`!Hx`_=*#VadSCS5q^2-vi|;;8|HFgV%iz}E z=oz#T;wHGp&az6}v_zT5)TK#zd=`XVPXnaDMY=g5N3(@(bc+;kA9{A=5RkdNyc}JK z^eT6#m&X1d+E2i!e#K(UkS7hmM#R3rZ29h1663$1m=C|?pgUKn1numg)f8>F#z`S@ zC^WG^W{Yujmt(el_c`y??qh2yI+>()GM$WNgYDqpcBkHt{f{R;2GMQmzkqdZ~oK)MG%`l8W_0u8Mee02RIXkg+KYr9Un3}|ZT z?H^>lS1iCx5iq%wF4rIQ*$RsG)UO>w+4W^63|T*5LLwei-$zB4XUF;$GK3-*^m|PC zM|tNbM5UOdyW3)Ny>oNwbP+)^Gd$Oh@xSk)g6tEr99Q(rHd_4J6!Shs?<^QUH}0bj$3aG zj}d`cr!QpPza31R>zIZ)rU3>$23rgZSXII0((lS!PH`x|~Dbi5VT@^sR*0Pehg0O#cRVh*2P)BNS zNL5G`QI2u+%{K$eJ&?`fj4)%{ULy)^q0Ucort zRoi)o$VyrmIwZFbo<_nxCy21j@gB#zVO-{egRa!rr;w|%JC94K`Ue6<4^@_3>+<+k4a7l6Lr=Yh>V^ z?$Vu&TDHSW5rV(8FgZM;49PT+?ddR^9b|7{3PRABjkRA>kqtt8$pdfQwiqAl7d`SJ zd$1dXAyaDD$OVy^N0_-Woo_}Ng!gK~e&T8GKCvf}G%a@=)myFTZSG(VrQ55^r;gmxJoqa6N8#OR z3FEejf{n^~+hB9X&nc-tKf99|AL{;0e`6!-DFS~bzoQ_~0AyO`I9Lj1t>%c#7!)V+ zzXk!8R0yWBC5L&Iy=_}G2m5$k=wp9ZXTKm7=E0)2C{pn;yfBG(UJ{hk!!^m4SZK0~ zznUEXD@i>c$EDadqd)hF($RCrQ;gy-fGjku*hX_nEo0R|{eo&?H(!-BsC7XK>?wct zr{kp8vjpd*z&98qKwJy4mA|;4dkT3S$-*dspwz`v;}a8||7D~Gl9+j}Nuqi9w3A2XtQ-;2LKML0v!Hr0b7y!w zSC@7F%N1;I8djkV=v|&wtv(rXvDMcb8xb&j_X!ux&vtpdLtao}D)GG3lr(G<-67q` zI41H!XnI0{J#kahhnd2?OnR1tj>04M&u-hQ(g$2YYtY}0g<6}&lePzVv<8<+6aKR5 z;PQu)_3XN%(|f&aX~K8ms)q1#=j&30ubZ73t^i{|>=V?@TSU9Z8@Ix=Ixg|V$-lci z9&wWIgLaga5NSfs6XAuM9Fvhc3*MszTn=rJ)3zrw@M)bE7Ss20lsAyV0GE3REzsS&!&pXQky#uLVnU_j@n$uU-e*6Qr*^)E_w_C7%pG%kF~o7 z73By;Z*0Dm$wigK_;oPkHHN(S(;!ig>0~;{m7zj>O5!FoaDvKt(82j{$Nz?{ZjU|M zGtWw~NRZ2uHw1RHcc&5(pB+it7h}ebEQWqD@mJUuw(bzbQ3(&nMtR#?Tcw}DpeZ49 zwPH10Ao055>d99bOe}0yBlH^3{HK=LBBLk%V-Bsa z#8FCNH27SxhRBDncpGS#j)&sgep%@`H`nrPJ(d|r7(#`ZmTQb>>?YBLe4UD6i=61M zGVT8)#5FNfn-RBPb$nnFFFu=fQtKodff z5|-zD`%8;Iu{rdO|6U?lUK-kB-lHH@?->7nL_F4^y{V%4N@HXysy*Rw7rP!pW!ZTV zU;NPgF0H{|)#kc!i#kSfk+V{fvjU^ZulBzy-dU-Aaov8PEeVp)PAQVID?##qgCniR z^1hfw%gktE+aM8hjOV|X2qe3 z5@D^mH!CnD37Ho~!-Q*VYn92CY`lNEfeGpgG`)%f0^!EvU$m0199;aq8P;dU*)DOl z`4zbnax<%#a=m!wMn**QpGOYo<*AC}OIG&yAg2;M6kpkL3eqV6gv92;nNXB_!OoQZ zM@oeCEaKz-~KHU@BO7!;!MtvAg{DIU#!^#NZ?#<^1aOXEN@igV7LQRW;&JS z_?6gwNVTN*L)TpH%O=<*NF}5}bG*u@i~jL&Ne$&R4`#C~>X#z@{&lRha;biToOSBw zq?vzXVmFiyEoW$Fo|LmzJ`sx!?T9>t=XSc~XOJ@HVhdv75@NMZ?+MVB`42bJNZ(a? zuM>&$VJLNdV`I7wd*3YeJxNpyo4U!LW;`A;G#hqO6_U;M|=xg-+MbxUDk1t+r&sP<3dOFw>G@d5H7diwIbbXsYA~`>SWU(=k-~;|PBq@I+ z?DyMDki!rO##s>CXs_gv_wn)hl1tv}>e2ugoRTJ*jOI!ieLU^PvlN9&8q;}`@92GE z01unp&;N!HdWd;rpLmV=u%wTFo~U-NDOTpj)mPSUa18wTQ1h{%vC*p~cGK1(QvVbN z^otxy1=hG7Xmc8X7HzP&InCbb;Dt(}Qc@O(z*h<~D$?bWUxWS! D$m3|G literal 0 HcmV?d00001 diff --git a/tests/_images/Points_datashader_matplotlib_stack.png b/tests/_images/Points_datashader_matplotlib_stack.png new file mode 100644 index 0000000000000000000000000000000000000000..58713a708607aad9c68e1a70497101ffe843b20d GIT binary patch literal 14365 zcmV+&IO4~NP)V>bE;%kREFfrfbZ~PzFE4FjbZ~5MbZlv2E^l&Y zFEj*5wg3PC32;bRa{vGaCjbBjCjm`V*v9|>H;PF_K~#90?VWdc6xG+qzcX9YD-c2w zAfXFL7X(41>rVk2px99q6a_&+DN;m5K#KG#D2NCll+YoxP(w>dfW$x&(tEOfcjml* zOn{JVwj{Gj^yTwBd6e0id+*Ggx#iqb7$PD(<}o6|)^r0N^WWpKP(4O;9t+iDMCY+k z*}&bqck$hK-`RY;H1EFqu0z^Xn~~OZ18OTbZrs4!xpUFHd2=KtCZbWJM###_!k1ru z2`48f3>-KR9XfPCL_`F>`Q{sN9EYEN`U#5{FUFlaci`vehsl#C*Ys`F2#;ug1=uI{^S~+O)ykxpSe_YEf8N2%hKR=H`a5u&|nH zdX4i~sI0&+%!B_YqDTH+mV$x;+_`fHvuDqSmzP&{G_O|jc=)S;4?g$+pMU;2f`fyR zlaqr+jT&Ltuwj@uaUxt@T`^?H5CA}QbTp<+nF42LX9NTUASWjWKmYtQQd3iF_6BMY z=CShAJUKZzSiE>KGBPso(MKPlbLY-A-K2jFkA>*w|ROxVVT>++%7v zYPA}zTD5||zrWevMG>ml*jOYcCZbK7HewW4c}|@|-uh4+xYY=ZiwpWZN2ojLUEtOkz+TZEL&ZdDg5B0#Sf+Hk)tg8^U|KqfPP{+u~DdQygTy&NNZ zNdPgQ7taG}X@Ek3%z`>7C;+^@#Vg(({+%j?&`@0e{W|9SC;@-)@WbQ1d&56K0tN(t zP{&^eNToo#c4*pK4hc&LYNixNp#J<5`1Dg?-0G@{R zNGbQT(J1_xJ1z?O`J{dAEo$TP554Yxg6QsDqGiiya^o!&?379=KmBg@GJSFexhM=I zmzpT%ie99~2yoNCqabAlxl4;FY2`lQU~Vo^LIRP&AZlIK-uJ(=ry4pGb-g)AC4?SB z>y-OhDi!|B9Tx>>XDrJZhG=zT>`3npldC7tv?+eM`w9}Z{usS;fZ5CBf%Tx$$)MFs z@I*VgNR1INbbl79dS4WoTyQAA3y_&9^qdnW0BhF@y=!TRfVemyEKC?S#iTkU_rFMH zd;!(>)l;W-B|VPGM1%VjkdLR4u3R$JNSmfj($Fh|4a9Lobm0Px4GJR#qoY=?<;FR+E(lELb2sm&{D@%lGiW$&1by(8qvdr_HD?8a9Nneh>hFcXq8< zpyJ|yZQFo^1e>dNtWc>`FdaOKW;=dBtJ%ZBA3h3bG{DA^q0})e?BK~BC`2zCtfFFLq z{(0AtH|tkn#IXU#jz>;TIHBl=1g$>+!sR*h(CJ1r_>@%m9{2+qPKgVBW6<{9h834l zTnzm56C$pAqh9^1Xu08Q`xG!HlTe1cJ0Ou%`Ibb((-ZGRkARMsW8}sl?EIwUxWNFN zJP8yQ0?$1sC`uhcdOG;IKO!T(EezR*pjg}7yqaYq=EL7K z8lr*%BAz#Y?Y(PL_njVb|Ho2BsWsn_EMANlP4*+-!R(PW5dJ9ns;&ou%#@>Av$snzAr;dvVKzq`bOqx#;u2hZ8sL|uHLG8l+*a)?Z&RzNBh z4Q&}q08cJT0};>Du0ORjca9CKB%Y^JOD>ayF_DwpM4YKafKpS5#*QU=^G%{dhsnmwE7_1B4(E`3C>j5-rdo2;j6*D6;>c6OPMDLkBJHQ7Ke5;c9( z^@#c2YUfUpLQleqfurK*UZP)`ub{#i3o7p2Dnh3d8euB)A>`%~oj*^Mm1XxmnMI|d z=6?CaG8~O~$x1IS0x2mZs37d1dWGulUAozOGp~(fYg$;m`-E4;t9E@F|K(k;Y zc7&PPzy$z>0*9k~!SgH*obg0RQx+12(5Sgm*!Z=Tv^4BL0J!elh+Z>#JzNE8G}wEg z4!n>21@Cq9fYelA;X+|$8x>`}+~j27yYGORGll)a%a?(e7_h{_SLCh1YDH|x)$Z6w zt(|U>D_cxkhW=^)FymkyaONCFf6KIXKU2ILCmiTRNpIz%*V7}0z05yi(p zxF1@yh|+tHpv?chF9=o{s((YdZClLwRs(kzBgTDh#}HGN^mJipm^>M{d9(86q@@AJ zj#XKE?;YrZ$S5T~7+q$lWf%;8RVFAHgM%@?`91(VR3;@@aj2Mj>pr`7QBFk$Ah zIxsOZ@L-T%ztyAXUAqR18wY&$84ww1xlDth{N`~5@TVi;7rGp7n|?*|ky|APdBk5+ zBQRnU0`{#$MqHjvWj%rz0`HIw<*i$3Ve=K_E6bo}P3~AO)kN&2N;LwdnPaPVjl2XxA>H zk3SX`6J$SN+4ESe!aV?ot&2L;>Dr z6mzzanMs>FFCrEua$oxBTk)T9lDuW6bLs{=GtFGSPxB zQ_07}Kzl+fY0$=T^!wOIVqqd@xzT=wjJ@dfl>PP>fNytvi))c9-Hf-xCR1f5fVS5Houxuw@ItviR(|8}ML@P*;}sa2}eR4D8zn+`U`% zw%CZ-4ds}a7)+fy74zrM$CfQy008~__lHCxL3nsLi0DDIle1^fB04(yQBRuf+p%G_ z0qfJA#U7{A@EkG7;WuK5C)!q%9T&r((>DX~GpQkSR-kovkE)f&^MWn>(W5~7_Q2@T zfJ}yd3kIS6u|gQq^3XG47W&_g#$V4Dfp>L>ig!gJ?}jZW8(`$qD&YI?vExDjIyF0i zhO57Q=p^DuD(b>EY}h~-E?l663l|a*(X?sPXxg-CL`1Y?$r8#cnV&>PM$)QPs~+tR zSFIu?E6f1~b*{zxY_*;Ie7gJMr{pWkr;Z-Ci0<5Z_*-Eib@#hX0FjeKPo&ckmz6^e zywv3GY@n2XbVNBhG_mPz7N}F+p5JD&vS6C|J(f;RMed4i z7Ed*6PMngGg0H{+8o|NAkjZ3VSr+~J^#jYYIB?(qhzK)g%y@AABj%;|UxHF;-v6y8 z*t-|_`)?p92$(ogu&uI{&``lSpmAg1%P#>}SI4!}44RB0bnhiWx2G8_-(Z`$Z&~tk zc({AARkP+kb`0JA=Z$mPC*dJIg8e65D>NrITTNhK;DfnTS-N!TV!862JB61J5>g{C zp&UjdOnVQY;nts_Cpi*&-iBx2PPSdE*FzO{6d_v{f?&+b`EN0t>eTA=oT(Vp2-dV zBm3aFqu=7Hu?0M3YW%U;@}=Smm0XUDz)r|Cc!I#-OmYwgi&{Y?U>Jltn~Y|0ary1uKbw06Z{{h>Jtr{Trd-l-OABD4y;q z!TB&Afq?+51Y!|as8|*OQ{KkscMFiIa>3wtJVhrmg3Y*|bxcqVN5J#8>E(<@VE%mEibzMj#?J6u^c~>o3H0s_6NJ6Uoce7s zZ=)aP5dn-VN~X{iBEAAgp}@k8|B?s4c{6@rufv(b)>yXbJlI~(ls|_^P%iGYOv4mEOO^52?y>rvybjASH^mIX} zTDHW^I33~>jOa7SPdL{`8?>5rd)H`$hg1`B*N6a}4)4cI#LDCsFmF>m^XKVw$l7-i zZrj(v^XG4%iM><3?JSK(1I(I*Gvoe5&gA*Rzsr)4fZZ!n;QPxr`1#9p!N>AGMkA1& z{iu4x%?*Y&9RV=lS*b7tmJrlbi9Unu?;vVUsIsyIvzD=Ah1Ftp!!W=bZ=kNb7D^e1 zzV8N_KTn~+&AL6I=9DNlIU$t|61N>3XZuH0T#Olqp2johen!ZKsphx9+}toGW-M}z z9+-3We-GW#^m<_4JmGmSUHXV%jhY*R^9S>BH(6qS3ScKS=Z78}i}7K-uq@>zbUvAa zGcR)0@0S-8gnI?||53rCU|ID1`9*xk7Qn?}$#~~$H%$B51V~GQo@Hla$)9yF@x~tr*`JBHw+jIuAK=Il zV8aGnw=OcZ-XOT|TIEI=zV4<;hgXy}C|C2n#uz2Hj6 z^Q50LgWM!)QZPE&8}t!%c8{ZrZ!EAFeyH|^5;7PBYfK*>oWAOYvG1yI?xLNV_oyhu z%nL`m3B7T8#-WG9M@|lK;slVAB7VEYhGlVue-eq>`bg5&!?~0Omd`(bUI-JgU;$c& z7~$$>X3~7{Aoww#g8yOyaHq`N&&dfpFEzsWQAId;l7!t3hCxy5ZqRTF7)XY;>MYMXg5Qaw#3JdAE_PNxnQ}IJ5S#y|7WSsC71u4_1z9NUz z2hBJ}@H`p6nNBXOhGdL^t}QxhTN^lz{$86$%YG`P$*;$glgvOL4=ZJ{nwv|se!bv` zTO?S4)-EHxtSpM^Gm&JBft)2;GTyys{TjW#yjO}&N8b&@T}2@I1v#*-R`<#DJp% zYs`Bt7w8Q-Q!5RWhW8UHlZ>YXwN&_i2+Aa)O*CW(Nnj*^r`aQ7h^V#+ z6%moaKs(pxQDlS_f#Nu#$&;z6b27O}R1~|`jN6hKj9R1Rp`;>ACMx{wTauvUOaCFd zq(GxEv}(2V+4J$_s?gJ~Kbn~c|JN`Y$@t~h)Y&bT+?Yc0c2-kTtY8ms?R{6jN%&A6Y`@Y3(`_hTs;uNKNhEN7b10lvnqc6i2<4zT(xmO%@Q0ZC?ahMg_x28 zEM5#~wSp;e+iJAXrqz3NpKtRtBn=h^1LZwt<%F+Pq=ZQ&a1#Dm-_T|)p zo)p-Z+@|URrBVq(l9&h#86xb>S>p6*bY1%at`-I3m$c*fB-BAdMI>zbD3y=gziZ1t zq*TI}^xqls+EL`?309&PE(kqJvu23cr$N82ML2${%=ed{4=i5}{PN4A-p|*s1K)oy z#K|mcz8oKquhun#V#5mb|2X*3eK~W4gdj%M2zruZ%PkuX3+?6phWXgBg1}pAlWk>I z9Z~)M^m<`gzjrTDVKuI@c%JCTA8B@z^)xDQKUoc=_19nYSzsu6%L?e!7k>-K3RHOA5WxuaF9Js~A(Fq8JQ?g>