Fix conversion of units of multivariate DataArray in sdba (#1972)

<!--Please ensure the PR fulfills the following requirements! -->
<!-- If this is your first PR, make sure to add your details to the
AUTHORS.rst! -->
### Pull Request Checklist:
- [ ] This PR addresses an already opened issue (for bug fixes /
features)
    - This PR fixes #xyz
- [ ] Tests for the changes have been added (for bug fixes / features)
- [ ] (If applicable) Documentation has been added / updated (for bug
fixes / features)
- [ ] CHANGELOG.rst has been updated (with summary of main changes)
- [ ] Link to issue (:issue:`number`) and pull request (:pull:`number`)
has been added

### What kind of change does this PR introduce?

Conversion of units of multivariate DataArray is now properly handled
`sdba.TrainAdjust` and `sdba.Adjust`. There was a bug where the units
could be changed before a conversion of the magntitudes could occur.

e.g. 

```python
out = sdba.dOTC.adjust(sim=sim,ref=ref,hist=hist)
```
`hist` would be converted to `ref` units properly, and then, the units
of `sim` would be changed too, without any actual conversion taking
place. This may have something to do with the weirdness of Generators.



### Does this PR introduce a breaking change?

No 

### Other information:

Just to show where I'm getting at, you can test: 

```python
from xclim import sdba
from xclim.testing import open_dataset
ds = open_dataset('sdba/ahccd_1950-2013.nc').isel(location=0)
ds2 = open_dataset('sdba/CanESM2_1950-2100.nc').isel(location=0)
da = sdba.stack_variables(ds)
da2 = sdba.stack_variables(ds2)
ref = da.sel(time=slice("1950", "1979"))
hist = da2.sel(time=slice("1950", "1979"))
sim = da2.sel(time=slice("1980", "2009"))
out = sdba.dOTC.adjust(sim=sim.copy(), ref=ref.copy(), hist=hist.copy()) 
print(out["multivar"].attrs["_units"], out.sel(multivar="tasmax").mean().values)

>>> ['mm day-1', 'degC'] 14.25703882896165  # new xclim
>>>['mm day-1', 'degC'] 282.0954876255738  # old xclim
```

You can see that `scen_tasmax` ~ 15 with the new xclim, as it should,
it's reasonable temperature for degC, but with old xclim it's around
`scen_tasmax`~282, clearly still magnitudes of K. This happens because
internally, the units of `sim` were changed, but the values were not
converted.

CHANGELOG.rst

@@ -2,6 +2,14 @@
 Changelog
 =========
 
+v0.54 (unpublished)
+--------------------
+Contributors to this version: Éric (:user:`coxipi`).
+
+Bug fixes
+^^^^^^^^^
+* Conversion of units of multivariate DataArray is now properly handled in `sdba.TrainAdjust` and `sdba.Adjust`. There was a bug where the units could be changed before a conversion of the magntitudes could occur. (:pull:`1972`).
+
 v0.53.1 (2024-10-21)
 --------------------
 Contributors to this version: Trevor James Smith (:user:`Zeitsperre`).

xclim/sdba/adjustment.py

@@ -180,12 +180,21 @@ def __convert_units_to(_input_da, _internal_dim, _internal_target):
                     v: _inputs[0][_dim].attrs["_units"][iv]
                     for iv, v in enumerate(_inputs[0][_dim].values)
                 }
-            return (
-                __convert_units_to(
-                    _input_da, _internal_dim=_dim, _internal_target=_target
+
+            # `__convert_units_to`` was changing the units of the 3rd dataset during the 2nd loop
+            # This explicit loop is designed to avoid this
+            _outputs = []
+            original_units = list(
+                [_inp[_dim].attrs["_units"].copy() for _inp in _inputs]
+            )
+            for _inp, units in zip(_inputs, original_units, strict=False):
+                _inp[_dim].attrs["_units"] = units
+                _outputs.append(
+                    __convert_units_to(
+                        _inp, _internal_dim=_dim, _internal_target=_target
+                    )
                 )
-                for _input_da in _inputs
-            ), _target
+            return _outputs, _target
 
         for dim, crd in inputs[0].coords.items():
             if crd.attrs.get("is_variables"):