Multi-dimensional thresholds

[aulemahal]

Pull Request Checklist:

• This PR addresses an already opened issue (for bug fixes / features)
  • This PR fixes per-gridcell heat wave indice calculations? #1093
• Tests for the changes have been added (for bug fixes / features)
  • (If applicable) Documentation has been added / updated (for bug fixes / features)
• HISTORY.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?

Thresholds and other quantities passed as parameters of indicators can now be multi-dimensional DataArrays. xArray broadcasting mechanisms will apply. Those parameters are now annotated as "Quantity" in the signatures (xclim.core.utils.Quantity), instead of "str" as before. Attributes where such thresholds where included will now read "<an array>" (french: "<une matrice>") for these new cases.

Does this PR introduce a breaking change?

No, I don't think so!

As expected, it was quite easy to implement this! Only a few indice needed modifications, and only a few are not compatible :
fire_weather_indices (because of the use of map_blocks) and [snow|rain]fall_approximation (but only with the "brown" method). Both cases shouldn't be where this feature is used the most.

Other information:

I did not add any checks about the dimensions of these thresholds. I thought we didn't need to be stricter than xarray. And that xclim users should expect xarray behaviour.

For example: if tas is ('lon', 'lat', 'time') and thresh is ('lon',) then it gets broadcasted along the lat and time dimensions automatically. But if the two lon coordinates differ, then a combined (union) lon coordinate is found in the output. And I made the choice that it's not our fault, nor is it for us to detect.

@JeremyFyke if you have time to test this branch for your heat wave need, that would be wonderful. Writing tests for this will be harder than implementing it, simply because of the large number of indicators involved...

[bzah]

Removing PercentileDataArray would break a few things in icclim but nothing too serious I believe.
ping @pagecp

[aulemahal]

Indeed, sorry I should have tagged you guys when I made the change.

On our end, there were two motivations:

• Having a subclass of xr.DataArray breaks some xarray methods that either strictly expect DataArray instances or return DataArray instances no matter the input's class. Maybe @tlogan2000 remembers the bug better. I guess that this is a bug upstream but
• the class seemed not to bring much functionality, especially with the new "Quantity" annotation. To me the combo xc.core.utils.is_percentile_dataarray and xc.core.formatting.get_percentile_metadata seemed like a smaller impact version of the same idea.

But if this breaks to many things in icclim, we can figure out a smoother way to fix the upstream bug while retaining the use you make of it.

[Zeitsperre]

Neat addition!

[bzah]

Neat indeed!
I would only suggest to use another name for Quantity as it's shadowing pint's Quantity and could be confusing (especially since it is already confusing between pint.Quantity and units.Quantity)
Maybe "Threshold" would make sense or "QuantifiedThreshold" ?

[aulemahal]

Indeed you are right on this.
I'm not sure about "Threshold" because it has become more than that here (ex: wb_cal in water_budget ). So this type is for anything that isn't a "timeseries" variable.

• QuantityType ? (we aren't using the "type" suffix anywhere else though)
• Parameter ? (now confusing with the Indicator class itself where all arguments are called "parameters", but that's hidden for most users. It also doesn't reflect the need for units)
  • ParamVar?
• Quantified ?

[tlogan2000]

I haven't dived into the code yet but basically the 'threshold' in all cases will now be an array right? ... How about ArrayValue or something? or ArrayParam ?

[aulemahal]

That's the other thing : the goal is to allow both arrays and scalar values. So the object can be either one of

• a string (value + units) (the most common case)
• a pint.Quantity (scalar with units)
• a DataArray with units (the novelty of this PR)

In most cases, it wouldn't be an array.

[tlogan2000]

I think the small bug with the PercentileDataArray was when I tried to plot an annual cycle directly ... something like doy as x-axis with percentiles values over the year in different curves .. e.g. doy_perc.plot.line(x='dayofyear') . I don't think it is critical to change as can be worked around fairly easily and is not something that is done too often.

[JeremyFyke]

Sorry for late arrival here - thanks @aulemahal for this work! I will definitely try this out during my next foray into the heat wave work. In general, this capability is really relevant to on-the-ground work, even outside of heatwaves.

[review-notebook-app]

Check out this pull request on 

See visual diffs & provide feedback on Jupyter Notebooks.

---

Powered by ReviewNB

[bzah]

I like "Quantified" best.

As for PercentileDataArray, it would make sense to simply migrate it into icclim code base as it is used in icclim own Threshold class.