Wrong units when using `da.integrate()`

[astoeriko]

Xarray has the nice xarray.DataArray.integrate method to integrate a data array over one of its coordinates. I automatically accounts for the values of the coordinates. The integration results should have units of the original array times the units of the coordinates. However, if I integrate a unit-aware Dataarray, it does not change its units. Here is a small example:

import numpy as np
import xarray as xr
import pint_xarray

t = xr.IndexVariable("time", np.arange(20), attrs={"units": "s"})
velocity = xr.DataArray(
    np.random.rand(20),
    dims="time",
    coords={"time": t},
    attrs={"units": "m/s"}
).pint.quantify()

velocity.integrate("time")

xarray.DataArray

    9.440035502696144 m/s
    Coordinates: (0)
    Indexes: (0)
    Attributes: (0)

I guess the problem is that only the data array is unit-aware, but not its coordinate. If this is not an easy thing to fix,it would be good to strip the unit completely or at least warn the user that the unit may be incorrect.