content: from __future__ import annotations
  from numpy import nan
  
  import numpy as np
  import pandas as pd
  import xarray as xr
  
  import datashader as ds
  
  import pytest
  
  from datashader.datatypes import RaggedDtype
  
  try:
      import spatialpandas as sp
      from spatialpandas.geometry import LineDtype
  except ImportError:
      LineDtype = None
      sp = None
  
  try:
      import cudf
      import cupy
  except ImportError:
      cupy, cudf = None, None
  
  def _pandas():
      """
      x          0  0   0  0  0   0   0  0  0  0   1   1  1   1  1    1  1   1  1   1
      y          0  0   0  0  0   1   1  1  1  1   0   0  0   0  0    1  1   1  1   1
      i32        0  1   2  3  4   5   6  7  8  9  10  11 12  13 14   15 16  17 18  19
      f32        0  1 nan  3  4   5   6  7  8  9  10  11 12  13 14   15 16  17 18  19
      reverse   20 19  18 17 16  15 nan 13 12 11  10   9  8   7  6    5  4   3  2   1
      plusminus  0 -1 nan -3  4  -5   6 -7  8 -9  10 -11 12 -13 14  -15 16 -17 18 -19
      cat2       a  b   c  d  a   b   c  d  a  b   c   d  a   b  c    d  a   b  c   d
      """
      df_pd = pd.DataFrame({'x': np.array([0.] * 10 + [1] * 10),
                        'y': np.array([0.] * 5 + [1] * 5 + [0] * 5 + [1] * 5),
                        'log_x': np.array([1.] * 10 + [10] * 10),
                        'log_y': np.array([1.] * 5 + [10] * 5 + [1] * 5 + [10] * 5),
                        'i32': np.arange(20, dtype='i4'),
                        'i64': np.arange(20, dtype='i8'),
                        'f32': np.arange(20, dtype='f4'),
                        'f64': np.arange(20, dtype='f8'),
                        'reverse': np.arange(20, 0, -1),
                        'plusminus': np.arange(20, dtype='f8')*([1, -1]*10),
                        'empty_bin': np.array([0.] * 15 + [np.nan] * 5),
                        'cat': ['a']*5 + ['b']*5 + ['c']*5 + ['d']*5,
                        'cat2': ['a', 'b', 'c', 'd']*5,
                        'onecat': ['one']*20,
                        'cat_int': np.array([10]*5 + [11]*5 + [12]*5 + [13]*5)})
      df_pd.cat = df_pd.cat.astype('category')
      df_pd.cat2 = df_pd.cat2.astype('category')
      df_pd.onecat = df_pd.onecat.astype('category')
      df_pd.at[2, 'f32'] = nan
      df_pd.at[2, 'f64'] = nan
      df_pd.at[6, 'reverse'] = nan
      df_pd.at[2, 'plusminus'] = nan
      return df_pd
  
  
  def _cudf():
      return cudf.DataFrame(_pandas())
  
  
  _backends = [
      pytest.param(_pandas, id="pandas"),
      pytest.param(_cudf, marks=pytest.mark.gpu, id="cudf"),
  ]
  
  @pytest.fixture(params=_backends)
  def df(request):
      return request.param()
  
  
  def _pandas_DataFrame(*args, **kwargs):
      if kwargs.pop("geo", False):
          return sp.GeoDataFrame(*args, **kwargs)
      else:
          return pd.DataFrame(*args, **kwargs)
  
  
  def _cudf_DataFrame(*args, **kwargs):
      import cudf
      if kwargs.pop("geo", False):
          pytest.skip("cudf currently does not work with spatialpandas")
      return cudf.DataFrame(
          pd.DataFrame(*args, **kwargs), nan_as_null=False
      )
  
  _backends = [
      pytest.param(_pandas_DataFrame, id="pandas"),
      pytest.param(_cudf_DataFrame, marks=pytest.mark.gpu, id="cudf"),
  ]
  
  @pytest.fixture(params=_backends)
  def DataFrame(request):
      return request.param
  
  c = ds.Canvas(plot_width=2, plot_height=2, x_range=(0, 1), y_range=(0, 1))
  c_logx = ds.Canvas(plot_width=2, plot_height=2, x_range=(1, 10),
                     y_range=(0, 1), x_axis_type='log')
  c_logy = ds.Canvas(plot_width=2, plot_height=2, x_range=(0, 1),
                     y_range=(1, 10), y_axis_type='log')
  c_logxy = ds.Canvas(plot_width=2, plot_height=2, x_range=(1, 10),
                      y_range=(1, 10), x_axis_type='log', y_axis_type='log')
  
  axis = ds.core.LinearAxis()
  lincoords = axis.compute_index(axis.compute_scale_and_translate((0, 1), 2), 2)
  coords = [lincoords, lincoords]
  dims = ['y', 'x']
  
  
  def assert_eq_xr(agg, b, close=False):
      """Assert that two xarray DataArrays are equal, handling the possibility
      that the two DataArrays may be backed by ndarrays of different types"""
      if cupy:
          if isinstance(agg.data, cupy.ndarray):
              agg = xr.DataArray(
                  cupy.asnumpy(agg.data), coords=agg.coords, dims=agg.dims
              )
          if isinstance(b.data, cupy.ndarray):
              b = xr.DataArray(
                  cupy.asnumpy(b.data), coords=b.coords, dims=b.dims
              )
      if close:
          xr.testing.assert_allclose(agg, b)
      else:
          xr.testing.assert_equal(agg, b)
  
  def assert_eq_ndarray(data, b, close=False):
      """Assert that two ndarrays are equal, handling the possibility that the
      ndarrays are of different types"""
      if cupy:
          if isinstance(data, cupy.ndarray):
              data = cupy.asnumpy(data)
          if isinstance(b, cupy.ndarray):
              b = cupy.asnumpy(b)
  
      if close:
          np.testing.assert_array_almost_equal(data, b, decimal=5)
      else:
          np.testing.assert_equal(data, b)
  
  
  def assert_image_close(image0, image1, tolerance):
      def to_rgba_xr(image):
          data = image.data
          if cupy and isinstance(data, cupy.ndarray):
              data = cupy.asnumpy(data)
          shape = data.shape
          data = data.view(np.uint8).reshape(shape + (4,))
          return xr.DataArray(data, dims=image.dims + ("rgba",), coords=image.coords)
  
      da0 = to_rgba_xr(image0)
      da1 = to_rgba_xr(image1)
      xr.testing.assert_allclose(da0, da1, atol=tolerance, rtol=0)
  
  
  def floats(n):
      """Returns contiguous list of floats from initial point"""
      while True:
          yield n
          n = n + np.spacing(n)
  
  
  def values(s):
      """Get numpy array of values from pandas-like Series, handling Series
      of different types"""
      if cudf and isinstance(s, cudf.Series):
          try:
              return s.to_numpy(na_value=np.nan)
          except AttributeError:
              # to_array is deprecated from cudf 22.02
              return s.to_array(fillna=np.nan)
      else:
          return s.values
  
  
  @pytest.mark.gpu
  def test_cudf_concat():
      # Testing if a newer version of cuDF implements the possibility to
      # concatenate multiple columns with the same name.
      # Currently, a workaround for this is
      # implemented in `datashader.glyphs.Glyph.to_cupy_array`.
      # For details, see: https://github.com/holoviz/datashader/pull/1050
  
      with pytest.raises(NotImplementedError):
          dfp = pd.DataFrame({'y': [0, 1]})
          dfc = cudf.DataFrame(dfp)
          cudf.concat((dfc["y"], dfc["y"]), axis=1)
  
  
  def test_count(df):
      out = xr.DataArray(np.array([[5, 5], [5, 5]], dtype='i4'),
                         coords=coords, dims=dims)
      assert_eq_xr(c.points(df, 'x', 'y', ds.count('i32')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.count('i64')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.count()), out)
      out = xr.DataArray(np.array([[4, 5], [5, 5]], dtype='i4'),
                         coords=coords, dims=dims)
      assert_eq_xr(c.points(df, 'x', 'y', ds.count('f32')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.count('f64')), out)
  
  
  def test_any(df):
      out = xr.DataArray(np.array([[True, True], [True, True]]),
                         coords=coords, dims=dims)
      assert_eq_xr(c.points(df, 'x', 'y', ds.any('i64')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.any('f64')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.any()), out)
      out = xr.DataArray(np.array([[True, True], [True, False]]),
                         coords=coords, dims=dims)
      assert_eq_xr(c.points(df, 'x', 'y', ds.any('empty_bin')), out)
  
  
  def test_sum(df):
      out = xr.DataArray(values(df.i32).reshape((2, 2, 5)).sum(axis=2, dtype='f8').T,
                         coords=coords, dims=dims)
      assert_eq_xr(c.points(df, 'x', 'y', ds.sum('i32')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.sum('i64')), out)
      out = xr.DataArray(np.nansum(values(df.f64).reshape((2, 2, 5)), axis=2).T,
                         coords=coords, dims=dims)
      assert_eq_xr(c.points(df, 'x', 'y', ds.sum('f32')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.sum('f64')), out)
  
  
  def test_min(df):
      out = xr.DataArray(values(df.i64).reshape((2, 2, 5)).min(axis=2).astype('f8').T,
                         coords=coords, dims=dims)
      assert_eq_xr(c.points(df, 'x', 'y', ds.min('i32')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.min('i64')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.min('f32')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.min('f64')), out)
  
  
  def test_max(df):
      out = xr.DataArray(values(df.i64).reshape((2, 2, 5)).max(axis=2).astype('f8').T,
                         coords=coords, dims=dims)
      assert_eq_xr(c.points(df, 'x', 'y', ds.max('i32')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.max('i64')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.max('f32')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.max('f64')), out)
  
  
  def test_min_n(df):
      solution = np.array([[[-3, -1, 0, 4, nan, nan], [-13, -11, 10, 12, 14, nan]],
                           [[-9, -7, -5, 6, 8, nan], [-19, -17, -15, 16, 18, nan]]])
      for n in range(1, 7):
          agg = c.points(df, 'x', 'y', ds.min_n('plusminus', n=n))
          out = solution[:, :, :n]
          assert_eq_ndarray(agg.data, out)
          if n == 1:
              assert_eq_ndarray(agg[:, :, 0].data, c.points(df, 'x', 'y', ds.min('plusminus')).data)
  
  
  def test_max_n(df):
      solution = np.array([[[4, 0, -1, -3, nan, nan], [14, 12, 10, -11, -13, nan]],
                           [[8, 6, -5, -7, -9, nan], [18, 16, -15, -17, -19, nan]]])
      for n in range(1, 7):
          agg = c.points(df, 'x', 'y', ds.max_n('plusminus', n=n))
          out = solution[:, :, :n]
          assert_eq_ndarray(agg.data, out)
          if n == 1:
              assert_eq_ndarray(agg[:, :, 0].data, c.points(df, 'x', 'y', ds.max('plusminus')).data)
  
  
  def test_categorical_min(df):
      sol_int = np.array([[[0, 1, 2, 3], [12, 13, 10, 11]], [[8, 5, 6, 7], [16, 17, 18, 15]]],
                         dtype=np.float64)
      sol_float = np.array([[[0, 1, nan, 3], [12, 13, 10, 11]], [[8, 5, 6, 7], [16, 17, 18, 15]]])
      assert_eq_ndarray(c.points(df, 'x', 'y', ds.by('cat2', ds.min('i32'))).data, sol_int)
      assert_eq_ndarray(c.points(df, 'x', 'y', ds.by('cat2', ds.min('i64'))).data, sol_int)
      assert_eq_ndarray(c.points(df, 'x', 'y', ds.by('cat2', ds.min('f32'))).data, sol_float)
      assert_eq_ndarray(c.points(df, 'x', 'y', ds.by('cat2', ds.min('f64'))).data, sol_float)
  
  
  def test_categorical_max(df):
      sol_int = np.array([[[4, 1, 2, 3], [12, 13, 14, 11]], [[8, 9, 6, 7], [16, 17, 18, 19]]],
                         dtype=np.float64)
      sol_float = np.array([[[4, 1, nan, 3], [12, 13, 14, 11]], [[8, 9, 6, 7], [16, 17, 18, 19]]])
      assert_eq_ndarray(c.points(df, 'x', 'y', ds.by('cat2', ds.max('i32'))).data, sol_int)
      assert_eq_ndarray(c.points(df, 'x', 'y', ds.by('cat2', ds.max('i64'))).data, sol_int)
      assert_eq_ndarray(c.points(df, 'x', 'y', ds.by('cat2', ds.max('f32'))).data, sol_float)
      assert_eq_ndarray(c.points(df, 'x', 'y', ds.by('cat2', ds.max('f64'))).data, sol_float)
  
  
  def test_categorical_min_n(df):
      solution = np.array([[[[0, 4, nan], [1, nan, nan], [nan, nan, nan], [3, nan, nan]],
                            [[12, nan, nan], [13, nan, nan], [10, 14, nan], [11, nan, nan]]],
                           [[[8, nan, nan], [5, 9, nan], [6, nan, nan], [7, nan, nan]],
                            [[16, nan, nan], [17, nan, nan], [18, nan, nan], [15, 19, nan]]]])
      for n in range(1, 3):
          agg = c.points(df, 'x', 'y', ds.by('cat2', ds.min_n('f32', n=n)))
          out = solution[:, :, :, :n]
          assert_eq_ndarray(agg.data, out)
          if n == 1:
              assert_eq_ndarray(agg[..., 0].data,
                                c.points(df, 'x', 'y', ds.by('cat2', ds.min('f32'))).data)
  
  
  def test_categorical_max_n(df):
      solution = np.array([[[[4, 0, nan], [1, nan, nan], [nan, nan, nan], [3, nan, nan]],
                            [[12, nan, nan], [13, nan, nan], [14, 10, nan], [11, nan, nan]]],
                           [[[8, nan, nan], [9, 5, nan], [6, nan, nan], [7, nan, nan]],
                            [[16, nan, nan], [17, nan, nan], [18, nan, nan], [19, 15, nan]]]])
      for n in range(1, 3):
          agg = c.points(df, 'x', 'y', ds.by('cat2', ds.max_n('f32', n=n)))
          out = solution[:, :, :, :n]
          assert_eq_ndarray(agg.data, out)
          if n == 1:
              assert_eq_ndarray(agg[..., 0].data,
                                c.points(df, 'x', 'y', ds.by('cat2', ds.max('f32'))).data)
  
  
  def test_categorical_min_row_index(df):
      solution = np.array([[[0, 1, 2, 3], [12, 13, 10, 11]], [[8, 5, 6, 7], [16, 17, 18, 15]]])
      agg = c.points(df, 'x', 'y', ds.by('cat2', ds._min_row_index()))
      assert_eq_ndarray(agg.data, solution)
  
  
  def test_categorical_max_row_index(df):
      solution = np.array([[[4, 1, 2, 3], [12, 13, 14, 11]], [[8, 9, 6, 7], [16, 17, 18, 19]]])
      agg = c.points(df, 'x', 'y', ds.by('cat2', ds._max_row_index()))
      assert_eq_ndarray(agg.data, solution)
  
  
  def test_categorical_min_n_row_index(df):
      solution = np.array([[[[0, 4, -1], [1, -1, -1], [2, -1, -1], [3, -1, -1]],
                            [[12, -1, -1], [13, -1, -1], [10, 14, -1], [11, -1, -1]]],
                           [[[8, -1, -1], [5, 9, -1], [6, -1, -1], [7, -1, -1]],
                            [[16, -1, -1], [17, -1, -1], [18, -1, -1], [15, 19, -1]]]])
      for n in range(1, 3):
          agg = c.points(df, 'x', 'y', ds.by('cat2', ds._min_n_row_index(n=n)))
          out = solution[:, :, :, :n]
          assert_eq_ndarray(agg.data, out)
          if n == 1:
              assert_eq_ndarray(agg[..., 0].data,
                                c.points(df, 'x', 'y', ds.by('cat2', ds._min_row_index())).data)
  
  
  def test_categorical_max_n_row_index(df):
      solution = np.array([[[[4, 0, -1], [1, -1, -1], [2, -1, -1], [3, -1, -1]],
                            [[12, -1, -1], [13, -1, -1], [14, 10, -1], [11, -1, -1]]],
                           [[[8, -1, -1], [9, 5, -1], [6, -1, -1], [7, -1, -1]],
                            [[16, -1, -1], [17, -1, -1], [18, -1, -1], [19, 15, -1]]]])
      for n in range(1, 3):
          agg = c.points(df, 'x', 'y', ds.by('cat2', ds._max_n_row_index(n=n)))
          out = solution[:, :, :, :n]
          assert_eq_ndarray(agg.data, out)
          if n == 1:
              assert_eq_ndarray(agg[..., 0].data,
                                c.points(df, 'x', 'y', ds.by('cat2', ds._max_row_index())).data)
  
  
  def test_categorical_first(df):
      solution = np.array([[[0, -1, nan, -3],
                            [12, -13, 10, -11]],
                           [[8, -5, 6, -7],
                            [16, -17, 18, -15]]])
      for n in range(1, 3):
          agg = c.points(df, 'x', 'y', ds.by('cat2', ds.first("plusminus")))
          assert_eq_ndarray(agg.data, solution)
  
  
  def test_categorical_last(df):
      solution = np.array([[[4, -1, nan, -3],
                            [12, -13, 14, -11]],
                           [[8, -9, 6, -7],
                            [16, -17, 18, -19]]])
      for n in range(1, 3):
          agg = c.points(df, 'x', 'y', ds.by('cat2', ds.last("plusminus")))
          assert_eq_ndarray(agg.data, solution)
  
  
  def test_categorical_first_n(df):
      solution = np.array([[[[0, 4, nan], [-1, nan, nan], [nan, nan, nan], [-3, nan, nan]],
                            [[12, nan, nan], [-13, nan, nan], [10, 14, nan], [-11, nan, nan]]],
                           [[[8, nan, nan], [-5, -9, nan], [6, nan, nan], [-7, nan, nan]],
                            [[16, nan, nan], [-17, nan, nan], [18, nan, nan], [-15, -19, nan]]]])
      for n in range(1, 3):
          agg = c.points(df, 'x', 'y', ds.by('cat2', ds.first_n("plusminus", n=n)))
          out = solution[:, :, :, :n]
          assert_eq_ndarray(agg.data, out)
          if n == 1:
              assert_eq_ndarray(agg[..., 0].data,
                                c.points(df, 'x', 'y', ds.by('cat2', ds.first("plusminus"))).data)
  
  
  def test_categorical_last_n(df):
      solution = np.array([[[[4, 0, nan], [-1, nan, nan], [nan, nan, nan], [-3, nan, nan]],
                            [[12, nan, nan], [-13, nan, nan], [14, 10, nan], [-11, nan, nan]]],
                           [[[8, nan, nan], [-9, -5, nan], [6, nan, nan], [-7, nan, nan]],
                            [[16, nan, nan], [-17, nan, nan], [18, nan, nan], [-19, -15, nan]]]])
      for n in range(1, 3):
          agg = c.points(df, 'x', 'y', ds.by('cat2', ds.last_n("plusminus", n=n)))
          out = solution[:, :, :, :n]
          assert_eq_ndarray(agg.data, out)
          if n == 1:
              assert_eq_ndarray(agg[..., 0].data,
                                c.points(df, 'x', 'y', ds.by('cat2', ds.last("plusminus"))).data)
  
  
  def test_where_min_row_index(df):
      out = xr.DataArray([[0, 10], [-5, -15]], coords=coords, dims=dims)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds._min_row_index(), 'plusminus')), out)
  
  
  def test_where_max_row_index(df):
      out = xr.DataArray([[4, 14], [-9, -19]], coords=coords, dims=dims)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds._max_row_index(), 'plusminus')), out)
  
  
  def test_where_min_n_row_index(df):
      sol = np.array([[[  0,  -1, nan,  -3,   4, nan],
                       [ 10, -11,  12, -13,  14, nan]],
                      [[ -5,   6,  -7,   8,  -9, nan],
                       [-15,  16, -17,  18, -19, nan]]])
      for n in range(1, 7):
          agg = c.points(df, 'x', 'y', ds.where(ds._min_n_row_index(n=n), 'plusminus'))
          out = sol[:, :, :n]
          print(n, agg.data.tolist())
          print(' ', out.tolist())
          assert_eq_ndarray(agg.data, out)
          if n == 1:
              assert_eq_ndarray(agg[:, :, 0].data,
                                c.points(df, 'x', 'y',
                                         ds.where(ds._min_row_index(), 'plusminus')).data)
  
  
  def test_where_max_n_row_index(df):
      sol = np.array([[[  4,  -3, nan,  -1,   0, nan],
                       [ 14, -13,  12, -11,  10, nan]],
                      [[ -9,   8,  -7,   6,  -5, nan],
                       [-19,  18, -17,  16, -15, nan]]])
      for n in range(1, 7):
          agg = c.points(df, 'x', 'y', ds.where(ds._max_n_row_index(n=n), 'plusminus'))
          out = sol[:, :, :n]
          print(n, agg.data.tolist())
          print(' ', out.tolist())
          assert_eq_ndarray(agg.data, out)
          if n == 1:
              assert_eq_ndarray(agg[:, :, 0].data,
                                c.points(df, 'x', 'y',
                                         ds.where(ds._max_row_index(), 'plusminus')).data)
  
  
  def test_where_first(df):
      # Note reductions like ds.where(ds.first('i32'), 'reverse') are supported,
      # but the same results can be achieved using the simpler ds.first('reverse')
      out = xr.DataArray([[20, 10], [15, 5]], coords=coords, dims=dims)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.first('i32'), 'reverse')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.first('i64'), 'reverse')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.first('f32'), 'reverse')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.first('f64'), 'reverse')), out)
  
      # Using row index.
      out = xr.DataArray([[0, 10], [5, 15]], coords=coords, dims=dims)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.first('i32'))), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.first('i64'))), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.first('f64'))), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.first('f32'))), out)
  
  
  def test_where_last(df):
      # Note reductions like ds.where(ds.last('i32'), 'reverse') are supported,
      # but the same results can be achieved using the simpler ds.last('reverse')
      out = xr.DataArray([[16, 6], [11, 1]], coords=coords, dims=dims)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.last('i32'), 'reverse')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.last('i64'), 'reverse')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.last('f32'), 'reverse')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.last('f64'), 'reverse')), out)
  
      # Using row index.
      out = xr.DataArray([[4, 14], [9, 19]], coords=coords, dims=dims)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.last('i32'))), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.last('i64'))), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.last('f64'))), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.last('f32'))), out)
  
  
  def test_where_max(df):
      out = xr.DataArray([[16, 6], [11, 1]], coords=coords, dims=dims)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.max('i32'), 'reverse')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.max('i64'), 'reverse')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.max('f32'), 'reverse')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.max('f64'), 'reverse')), out)
  
      # Using row index.
      out = xr.DataArray([[4, 14], [9, 19]], coords=coords, dims=dims)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.max('i32'))), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.max('i64'))), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.max('f64'))), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.max('f32'))), out)
  
  
  def test_where_min(df):
      out = xr.DataArray([[20, 10], [15, 5]], coords=coords, dims=dims)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.min('i32'), 'reverse')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.min('i64'), 'reverse')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.min('f32'), 'reverse')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.min('f64'), 'reverse')), out)
  
      # Using row index.
      out = xr.DataArray([[0, 10], [5, 15]], coords=coords, dims=dims)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.min('i32'))), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.min('i64'))), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.min('f64'))), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.where(ds.min('f32'))), out)
  
  
  def test_where_first_n(df):
      sol_rowindex = np.array([[[ 0,  1,  3,  4, -1, -1],
                                [10, 11, 12, 13, 14, -1]],
                               [[ 5,  6,  7,  8,  9, -1],
                                [15, 16, 17, 18, 19, -1]]])
      sol_reverse = xr.where(np.logical_or(sol_rowindex < 0, sol_rowindex == 6), np.nan,
                             20 - sol_rowindex)
  
      for n in range(1, 7):
          # Using row index.
          agg = c.points(df, 'x', 'y', ds.where(ds.first_n('plusminus', n=n)))
          out = sol_rowindex[:, :, :n]
          assert_eq_ndarray(agg.data, out)
          if n == 1:
              assert_eq_ndarray(agg[:, :, 0].data,
                                c.points(df, 'x', 'y', ds.where(ds.first('plusminus'))).data)
  
          # Using another column
          agg = c.points(df, 'x', 'y', ds.where(ds.first_n('plusminus', n=n), 'reverse'))
          out = sol_reverse[:, :, :n]
          assert_eq_ndarray(agg.data, out)
          if n == 1:
              assert_eq_ndarray(agg[:, :, 0].data,
                                c.points(df, 'x', 'y',
                                         ds.where(ds.first('plusminus'), 'reverse')).data)
  
  
  def test_where_last_n(df):
      sol_rowindex = np.array([[[ 4,  3,  1,  0, -1, -1],
                                [14, 13, 12, 11, 10, -1]],
                               [[ 9,  8,  7,  6,  5, -1],
                                [19, 18, 17, 16, 15, -1]]])
      sol_reverse = xr.where(np.logical_or(sol_rowindex < 0, sol_rowindex == 6), np.nan,
                             20 - sol_rowindex)
  
      for n in range(1, 7):
          # Using row index.
          agg = c.points(df, 'x', 'y', ds.where(ds.last_n('plusminus', n=n)))
          out = sol_rowindex[:, :, :n]
          assert_eq_ndarray(agg.data, out)
          if n == 1:
              assert_eq_ndarray(agg[:, :, 0].data,
                                c.points(df, 'x', 'y', ds.where(ds.last('plusminus'))).data)
  
          # Using another column
          agg = c.points(df, 'x', 'y', ds.where(ds.last_n('plusminus', n=n), 'reverse'))
          out = sol_reverse[:, :, :n]
          assert_eq_ndarray(agg.data, out)
          if n == 1:
              assert_eq_ndarray(agg[:, :, 0].data,
                                c.points(df, 'x', 'y',
                                         ds.where(ds.last('plusminus'), 'reverse')).data)
  
  
  def test_where_max_n(df):
      sol_rowindex = np.array([[[ 4,  0,  1,  3, -1, -1],
                                [14, 12, 10, 11, 13, -1]],
                               [[ 8,  6,  5,  7,  9, -1],
                                [18, 16, 15, 17, 19, -1]]])
      sol_reverse = xr.where(np.logical_or(sol_rowindex < 0, sol_rowindex == 6), np.nan,
                             20 - sol_rowindex)
  
      for n in range(1, 7):
          # Using row index.
          agg = c.points(df, 'x', 'y', ds.where(ds.max_n('plusminus', n=n)))
          out = sol_rowindex[:, :, :n]
          assert_eq_ndarray(agg.data, out)
          if n == 1:
              assert_eq_ndarray(agg[:, :, 0].data, c.points(df, 'x', 'y',
                                                            ds.where(ds.max('plusminus'))).data)
  
          # Using another column
          agg = c.points(df, 'x', 'y', ds.where(ds.max_n('plusminus', n=n), 'reverse'))
          out = sol_reverse[:, :, :n]
          assert_eq_ndarray(agg.data, out)
          if n == 1:
              assert_eq_ndarray(agg[:, :, 0].data,
                                c.points(df, 'x', 'y', ds.where(ds.max('plusminus'), 'reverse')).data)
  
  
  def test_where_min_n(df):
      sol_rowindex = np.array([[[3,  1,  0,  4, -1, -1],
                                [13, 11, 10, 12, 14, -1]],
                               [[ 9,  7,  5,  6,  8, -1],
                                [19, 17, 15, 16, 18, -1]]])
      sol_reverse = xr.where(np.logical_or(sol_rowindex < 0, sol_rowindex == 6), np.nan,
                             20 - sol_rowindex)
  
      for n in range(1, 7):
          # Using row index.
          agg = c.points(df, 'x', 'y', ds.where(ds.min_n('plusminus', n=n)))
          out = sol_rowindex[:, :, :n]
          assert_eq_ndarray(agg.data, out)
          if n == 1:
              assert_eq_ndarray(agg[:, :, 0].data,
                                c.points(df, 'x', 'y', ds.where(ds.min('plusminus'))).data)
  
          # Using another column
          agg = c.points(df, 'x', 'y', ds.where(ds.min_n('plusminus', n=n), 'reverse'))
          out = sol_reverse[:, :, :n]
          assert_eq_ndarray(agg.data, out)
          if n == 1:
              assert_eq_ndarray(agg[:, :, 0].data,
                                c.points(df, 'x', 'y', ds.where(ds.min('plusminus'), 'reverse')).data)
  
  
  def test_summary_by(df):
      # summary(by)
      agg_summary = c.points(df, 'x', 'y', ds.summary(by=ds.by("cat")))
      agg_by = c.points(df, 'x', 'y', ds.by("cat"))
      assert_eq_xr(agg_summary["by"], agg_by)
  
      # summary(by, other_reduction)
      agg_summary = c.points(df, 'x', 'y', ds.summary(by=ds.by("cat"), max=ds.max("plusminus")))
      agg_max = c.points(df, 'x', 'y', ds.max("plusminus"))
      assert_eq_xr(agg_summary["by"], agg_by)
      assert_eq_xr(agg_summary["max"], agg_max)
  
      # summary(other_reduction, by)
      agg_summary = c.points(df, 'x', 'y', ds.summary(max=ds.max("plusminus"), by=ds.by("cat")))
      assert_eq_xr(agg_summary["by"], agg_by)
      assert_eq_xr(agg_summary["max"], agg_max)
  
      # summary(by, by)
      agg_summary = c.points(df, 'x', 'y', ds.summary(by=ds.by("cat"), by_any=ds.by("cat", ds.any())))
      agg_by_any = c.points(df, 'x', 'y', ds.by("cat", ds.any()))
      assert_eq_xr(agg_summary["by"], agg_by)
      assert_eq_xr(agg_summary["by_any"], agg_by_any)
  
      # summary(by("cat1"), by("cat2"))
      agg_summary = c.points(df, 'x', 'y', ds.summary(by=ds.by("cat"), by2=ds.by("cat2")))
      agg_by2 = c.points(df, 'x', 'y', ds.by("cat2"))
      assert_eq_xr(agg_summary["by"], agg_by)
      assert_eq_xr(agg_summary["by2"], agg_by2)
  
  
  def test_summary_where_n(df):
      sol_min_n_rowindex = np.array([[[ 3,  1,  0,  4, -1],
                                      [13, 11, 10, 12, 14]],
                                     [[ 9,  7,  5,  6,  8],
                                      [19, 17, 15, 16, 18]]])
      sol_max_n_rowindex = np.array([[[ 4,  0,  1,  3, -1],
                                      [14, 12, 10, 11, 13]],
                                     [[ 8,  6,  5,  7,  9],
                                      [18, 16, 15, 17, 19]]])
      sol_max_n_reverse = np.where(np.logical_or(sol_max_n_rowindex < 0, sol_max_n_rowindex == 6),
                                   np.nan, 20 - sol_max_n_rowindex)
  
      agg = c.points(df, 'x', 'y', ds.summary(
          count=ds.count(),
          min_n=ds.where(ds.min_n('plusminus', 5)),
          max_n=ds.where(ds.max_n('plusminus', 5), 'reverse'),
      ))
      assert_eq_ndarray(agg.coords['n'], np.arange(5))
  
      assert agg['count'].dims == ('y', 'x')
      assert agg['min_n'].dims == ('y', 'x', 'n')
      assert agg['max_n'].dims == ('y', 'x', 'n')
  
      assert agg['count'].dtype == np.dtype('uint32')
      assert agg['min_n'].dtype == np.dtype('int64')
      assert agg['max_n'].dtype == np.dtype('float64')
  
      assert_eq_ndarray(agg['count'].data, [[5, 5], [5, 5]])
      assert_eq_ndarray(agg['min_n'].data, sol_min_n_rowindex)
      assert_eq_ndarray(agg['max_n'].data, sol_max_n_reverse)
  
      # Issue #1270: Support summary reduction containing multiple where
      # reductions that use the same selector.
      agg = c.points(df, 'x', 'y', ds.summary(
          max1=ds.where(ds.max_n('plusminus', 5)),
          max2=ds.where(ds.max_n('plusminus', 5), 'reverse'),
      ))
      assert_eq_ndarray(agg['max1'].data, sol_max_n_rowindex)
      assert_eq_ndarray(agg['max2'].data, sol_max_n_reverse)
  
  
  def test_summary_different_n(df):
      msg = 'Using multiple FloatingNReductions with different n values is not supported'
      with pytest.raises(ValueError, match=msg):
          c.points(df, 'x', 'y', ds.summary(
              min_n=ds.where(ds.min_n('plusminus', 2)),
              max_n=ds.where(ds.max_n('plusminus', 3)),
          ))
  
  
  def test_mean(df):
      out = xr.DataArray(values(df.i32).reshape((2, 2, 5)).mean(axis=2, dtype='f8').T,
                         coords=coords, dims=dims)
      assert_eq_xr(c.points(df, 'x', 'y', ds.mean('i32')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.mean('i64')), out)
      out = xr.DataArray(np.nanmean(values(df.f64).reshape((2, 2, 5)), axis=2).T,
                         coords=coords, dims=dims)
      assert_eq_xr(c.points(df, 'x', 'y', ds.mean('f32')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.mean('f64')), out)
  
  
  def test_var(df):
      out = xr.DataArray(values(df.i32).reshape((2, 2, 5)).var(axis=2, dtype='f8').T,
                         coords=coords, dims=dims)
      assert_eq_xr(c.points(df, 'x', 'y', ds.var('i32')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.var('i64')), out)
      out = xr.DataArray(np.nanvar(values(df.f64).reshape((2, 2, 5)), axis=2).T,
                         coords=coords, dims=dims)
      assert_eq_xr(c.points(df, 'x', 'y', ds.var('f32')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.var('f64')), out)
  
  
  def test_std(df):
      out = xr.DataArray(values(df.i32).reshape((2, 2, 5)).std(axis=2, dtype='f8').T,
                         coords=coords, dims=dims)
      assert_eq_xr(c.points(df, 'x', 'y', ds.std('i32')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.std('i64')), out)
      out = xr.DataArray(np.nanstd(values(df.f64).reshape((2, 2, 5)), axis=2).T,
                         coords=coords, dims=dims)
      assert_eq_xr(c.points(df, 'x', 'y', ds.std('f32')), out)
      assert_eq_xr(c.points(df, 'x', 'y', ds.std('f64')), out)
  
  
  def test_count_cat(df):
      sol = np.array([[[5, 0, 0, 0],
                       [0, 0, 5, 0]],
                      [[0, 5, 0, 0],
                       [0, 0, 0, 5]]])
      out = xr.DataArray(sol, coords=coords + [['a', 'b', 'c', 'd']], dims=(dims + ['cat']))
      agg = c.points(df, 'x', 'y', ds.count_cat('cat'))
      assert_eq_xr(agg, out)
      assert_eq_ndarray(agg.x_range, (0, 1), close=True)
      assert_eq_ndarray(agg.y_range, (0, 1), close=True)
  
  
  def test_categorical_count(df):
      sol = np.array([[[5, 0, 0, 0],
                       [0, 0, 5, 0]],
                      [[0, 5, 0, 0],
                       [0, 0, 0, 5]]])
      out = xr.DataArray(sol, coords=coords + [['a', 'b', 'c', 'd']], dims=(dims + ['cat']))
      agg = c.points(df, 'x', 'y', ds.by('cat', ds.count('i32')))
      assert_eq_xr(agg, out)
  
      # ds.summary(name=ds.by("cat")) should give same result as ds.by("cat"). Issue 1252
      dataset = c.points(df, 'x', 'y', ds.summary(name=ds.by('cat', ds.count('i32'))))
      assert_eq_xr(dataset["name"], out)
  
      # categorizing by (cat_int-10)%4 ought to give the same result
      out = xr.DataArray(sol, coords=coords + [range(4)], dims=(dims + ['cat_int']))
      agg = c.points(df, 'x', 'y',
  