The CROPS Plugin ================ **authors:** | Florian Elsäßer, Etor E. Lucio-Eceiza [1]_\ | Zentrum für internationale Entwicklungs- und Umweltforschung (ZEU), Justus-Liebig-Universität Gießen, Deutsches Klimarechenzentrum (DKRZ)\ | Version v2506.1.0\ | `Link to GitLab repository `_ **CROPS** is a Freva plugin to quantify the impact of **extreme** and **compound** climate stressors on agricultural productivity, focusing on oat, rye, sugar beet, summer barley, winter barley, and winter wheat. The plugin takes in account the phenologycal cycles of crops and calculates the impact to them in regards to extremes events such as heatwaves, water surplus/scarcity or compounds among them. These impacts can vary depending on the time (e.g., heatwaves affect the flowering, while water when harvesting). Three climate stress indices can be computed: **HMD** (Heat Magnitude Day), **SPEI** (Standardized Precipitation-Evapotranspiration Index), **CSI** (Compound Stress Index). .. admonition:: Info - The core of the plugin is a modified version of the ``crop_indices`` python package (installable via PyPi) developed by Florian Elsäßer that can be found in GitHub [2]_. The actual SPEI code is taken from the ``climate_indices`` package [3]_. - The index calculations are based on Zampieri et al. (2017) [4]_ but where the CSI has been updated to standard deviations. Theoretical Background ----------------------- HMD (Heat Magnitude Day) ~~~~~~~~~~~~~~~~~~~~~~~~ Is calculated using the methodology of Zampieri et al. (2017), adapted to the crop-specific **harvest month** using phenological data. It is computed using daily maximum temperature (``tasmax``) values. Steps: 1. Identify harvest date per region/year based on phenology. 2. Extract a pre-harvest time window (user-defined, e.g., 3 months). 3. Calculate rolling **T90**, **T75**, and **T25** quantile thresholds from a reference period. 4. Mark days exceeding the **T90** threshold. 5. Detect heatwave events that meet the minimum number of consecutive exceedance days. 6. For each heatwave event, calculate: - Duration - Absolute exceedance magnitude - Top 3-day exceedances - Interquartile range (``T75 - T25``) - Normalized and standardized exceedances The **final HMD value** is the sum of all normalized exceedance intensities for a given year. Mathematically, the daily magnitude of exceedance is defined as: .. math:: M_d(T) = \begin{cases} \frac{T - T_{25}}{T_{75} - T_{25}} & \text{if } T > T_{90} \\ 0 & \text{otherwise} \end{cases} SPEI (Standardized Precipitation Evapotranspiration Index) ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ It follows Vicente-Serrano et al. (2010) [5]_ , and is adapted for phenology-aware crop phases. It evaluates water balance stress through the difference between precipitation and potential evapotranspiration, tailored to crop growth stages. Steps: 1. Compute monthly mean for `pr`, `tas`, `tasmin`, `tasmax`. 2. Estimate **evapotranspiration** via Hargreaves-Samani (1985) [6]_ method. 3. Compute **D-value**: ``D = precipitation - PET`` 4. Accumulate D-values over a user-defined pre-harvest window (e.g. 6 months). 5. Apply gamma distribution to fit and calculate the SPEI time series. CSI (Compound Stress Index) ~~~~~~~~~~~~~~~~~~~~~~~~~~~~ CSI combines HMD and SPEI into a single stress indicator using Ridge regression, optionally rescaled to match yield anomalies. Is based on Zampieri et al. (2017). Steps to take: 1. Use **standardized** and **detrended** HMD and SPEI indices. 2. Fit **Ridge regression** to yield data. 3. Compute: .. math:: \text{CSI}_t = -a \cdot \text{HMD}_t + -b \cdot \text{SPEI}_t 4. Optionally **rescale CSI** back to yield units using the historical yield mean and std. Agricultural Data ----------------- This plugin uses agricultural data produced at https://cropdata.de/request_data/: - **Yield productivity data** [7]_: Contains productivity data (in decitons/hectare) for crops including oat, rye, sugar beet, summer barley, winter barley, and winter wheat from 1989 to 2020. - **Phenology and management data** [8]_: Provides BBCH-phase data (sowing, flowering, harvesting) and day-of-year values for all relevant crops from 1989 to 2020. - **Natural areas** [8]_: A classification of 86 distinct environmental zones based on soil, water climatology, and water regimes across Germany (1989–2020). Usage ----- The indices are calculated by processing climate, yield production and phenology datasets across selected regions (shapefiles / ``regionmask`` library) and timeframes. The tool expects input files in NetCDF format and the Output files are NetCDF (`.nc`) and animated MP4s (`.mp4`) for visual inspection. The tool will run a whole index workflow per climate-data ``ensemble`` member **sequentially** .. important:: - The crop data is currently limited to that available e.g. via `freva/XCES `_ **from 1989 to 2020**: .. code:: bash freva databrowser model=crop-data project=observations or .. code:: python freva.databrowser(model="crop-data", project="observations") - The **region** of analysis is limited to **Germany** (whole, states, municipalitites). Input ~~~~~ - Data selection: - Climatological data via ``project``, ``product``, ``institute``, ``model``, ``experiment``, ``ensemble``: chooses automatically ``daily`` data of mean temperature (``tas``), maximum temperature (``tasmax``), minimum temperature (``tasmin``) and precipitation (``pr``) depending on the index. - Crop data via ``crop_type``: e.g., ``winter-wheat``, ``maize``. This automatically selects **natural areas** file, **yield** and **phenological cycle**. .. warning:: - Natural areas and crop yield data are regridded to the climatological grid using the `XESMF regridder `_ with the method ``nearest_s2d``. - Phenological phase data is averaged to a single value per natural area. - The phenological phase used corresponds to ``bbch99``, `representing the post-harvest or storage treatment stage. `_ - ``index``: Choose ``HMD``, ``SPEI``, or ``CSI`` (``CSI`` triggers the calculation of ``HMD`` and ``SPEI``) - ``HMD`` relevant parameters: - ``hmd_time_range``: Relevant period for HMD (in month before harvest month). E.g.: ``0,3`` where 0 is harvest month, e.g. July, and 3 is the starting Month, e.g. May. This would therefore consider May, June and July. - ``exceedence_period``: Number of days in exceedence analysis window (before and after the actual day). - ``exceedence_days``: Minimum days of threshold exceedence to qualify as a heatwave. - ``SPEI`` relevant parameters: - ``spei_time_range``: Relevant period for SPEI (in month before harvest month). E.g.: ``0,6`` where 0 is harvest month, e.g. July, and 6 is the starting Month, e.g. February. This would therefore consider February to July. - ``CSI`` relevant parameters: - ``ridge_mode``: CSI regression method, either ``cv`` (Cross-validation to find optimal α, relies on `RidgeCV `_ ) or ``min`` (Select α at first local minimum in MSE, relies on `Ridge `_). .. note:: ``cv`` option is considerable faster but the results may differ. - ``rescale_to_yield``: ``True`` or ``False`` to output CSI in yield units. - ``file_hmd`` / ``file_spei``: Optional NetCDF files containing precomputed HMD/SPEI data for use with ``CSI``. Required **only** when selecting ``CSI`` as the index. Must match the spatial grid and time range of the other input data. - The spatial region for analysis: - ``region_type``: Select the type of region mask. Options include: - ``countries``: Predefined countries using the `regionmask` Python package. - ``ger_states``: German federal states. - ``ger_munip``: German municipalities. - ``lonlat``: A custom user-defined box specified by geographic coordinates. - ``region_name``: Specify one or more region labels depending on the ``region_type``: - For ``countries``: **limited to** ``Germany`` - For ``ger_states`` or ``ger_munip``: e.g., ``Hamburg, Berlin`` - For ``lonlat``: Provide bounding box as ``lonmin,lonmax,latmin,latmax`` .. note:: All selected regions are processed **jointly**, not individually. - ``mask_method``: Defines how the spatial mask is applied: - ``centres``: Only includes grid cells whose **centres** lie inside the region. - ``corners``: Includes any grid cell if **any of its corners** touch the region polygon. This is less restrictive and often includes more boundary cells. - ``seldate``: Time range of analysis. Provide a comma-separated string in the format ``YYYY,YYYY``. For example, ``1989,2020`` will restrict the calculation to those years. By default, it uses the full available crop period: ``1989,2020``. Output ~~~~~~ The plugin generates the following output files, depending on the selected index: - **HMD**: ``hmd.nc`` (yearly HMD time series) and ``hmd.mp4`` ( animated visualization). - **SPEI**: ``spei.nc`` (yearly SPEI time series) and ``spei.mp4`` ( animated visualization). - **CSI**: ``csi.nc`` (yearly CSI time series) and ``csi.mp4`` ( animated visualization). - Also includes ``hmd.nc``, ``spei.nc``, ``hmd.mp4``, and ``spei.mp4`` if HMD and SPEI were internally computed. Each NetCDF file contains a **yearly time series** per region and crop. The MP4 files provide animated spatial visualizations for inspection and communication. Installation ============ To install the conda environment required by the plugin: 1. Clone and build the environment: .. code:: bash $ git clone https://gitlab.dkrz.de/freva/plugins4freva/crops.git $ cd crops $ make all 2. Activate the environment: .. code:: bash $ source $(pwd)/plugin_env/bin/activate 4. To deactivate: .. code:: bash $ conda deactivate Updating the Plugin with Git ---------------------------- 1. ``cd `` 2. ``git checkout -b `` 3. ``git push --set-upstream origin `` 4. Open a merge request on GitLab and wait for review Contact & References ==================== .. [1] `lucio-eceiza@dkrz.de `_ .. [2] `GitHub code of crop_indices `_ .. [3] `Original Python adaptation and SPEI from climate_indices `_ .. [4] Zampieri, M., Ceglar, A., Dentener, F., Toreti, A., 2017. *Wheat yield loss attributable to heat waves, drought and water excess at the global, national and subnational scales.* Environ. Res. Lett. 12, 064008. https://doi.org/10.1088/1748-9326/aa723b .. [5] Vicente-Serrano, S.M., Beguería, S., López-Moreno, J.I., 2010. *A Multiscalar Drought Index Sensitive to Global Warming: The Standardized Precipitation Evapotranspiration Index.* J. Clim. 23, 1696–1718. https://doi.org/10.1175/2009JCLI2909.1 .. [6] Hargreaves, G.H., Samani, Z.A., 1985. *Reference Crop Evapotranspiration from Temperature.* Appl. Eng. Agric. 1(2): 96–99. https://doi.org/10.13031/2013.26773 .. [7] cropdata – spatial yield productivity data base for the ten most cultivated crops in Germany from 1989 to 2020 - version 1.0. http://dx.doi.org/10.22029/jlupub-7177 .. [8] cropdata – supplementary data (for spatial yield productivity data base for the ten most cultivated crops in Germany from 1989 to 2020) - version 1.0. http://dx.doi.org/10.22029/jlupub-7203