|
| 1 | +import numpy as np |
| 2 | +import xarray as xr |
| 3 | + |
| 4 | +from xrspatial.dasymetric import disaggregate |
| 5 | + |
| 6 | +from .common import get_xr_dataarray |
| 7 | + |
| 8 | + |
| 9 | +class Dasymetric: |
| 10 | + params = ([100, 300, 1000], ["numpy", "cupy", "dask"]) |
| 11 | + param_names = ("nx", "type") |
| 12 | + |
| 13 | + def setup(self, nx, type): |
| 14 | + ny = nx // 2 |
| 15 | + |
| 16 | + # Zones: 4 rectangular blocks. |
| 17 | + zones_np = np.zeros((ny, nx), dtype=np.int32) |
| 18 | + zones_np[: ny // 2, : nx // 2] = 1 |
| 19 | + zones_np[: ny // 2, nx // 2 :] = 2 |
| 20 | + zones_np[ny // 2 :, : nx // 2] = 3 |
| 21 | + zones_np[ny // 2 :, nx // 2 :] = 4 |
| 22 | + |
| 23 | + x = np.linspace(-180, 180, nx) |
| 24 | + y = np.linspace(-90, 90, ny) |
| 25 | + |
| 26 | + if type == "dask": |
| 27 | + import dask.array as da |
| 28 | + zdata = da.from_array(zones_np, chunks=(max(1, ny // 2), max(1, nx // 2))) |
| 29 | + elif type == "cupy": |
| 30 | + from xrspatial.utils import has_cuda_and_cupy |
| 31 | + if not has_cuda_and_cupy(): |
| 32 | + raise NotImplementedError() |
| 33 | + import cupy |
| 34 | + zdata = cupy.asarray(zones_np) |
| 35 | + else: |
| 36 | + zdata = zones_np |
| 37 | + |
| 38 | + self.zones = xr.DataArray(zdata, coords=dict(y=y, x=x), dims=["y", "x"]) |
| 39 | + |
| 40 | + # Values: one total per zone. |
| 41 | + self.values = {1: 1000.0, 2: 2000.0, 3: 1500.0, 4: 2500.0} |
| 42 | + |
| 43 | + # Weight surface: use the standard Gaussian bump. |
| 44 | + weight = get_xr_dataarray((ny, nx), type) |
| 45 | + # Make weights non-negative. |
| 46 | + if type == "dask": |
| 47 | + import dask.array as da |
| 48 | + weight.data = da.fabs(weight.data) + 0.01 |
| 49 | + elif type == "cupy": |
| 50 | + import cupy |
| 51 | + weight.data = cupy.abs(weight.data) + 0.01 |
| 52 | + else: |
| 53 | + weight.data = np.abs(weight.data) + 0.01 |
| 54 | + self.weight = weight |
| 55 | + |
| 56 | + def time_disaggregate_weighted(self, nx, type): |
| 57 | + disaggregate(self.zones, self.values, self.weight, method="weighted") |
| 58 | + |
| 59 | + def time_disaggregate_binary(self, nx, type): |
| 60 | + disaggregate(self.zones, self.values, self.weight, method="binary") |
0 commit comments