Analysis/event_detect

[tommylees112]

Event detection using thresholds and climatology values. An analysis piece to determine drought events / extremes in the data we're looking at

[tommylees112]

In [1]: from pathlib import Path
   ...: from src.analysis import EventDetector
   ...:
   ...: data_dir = Path('data')
   ...: interim_dir = data_dir / 'interim'
   ...: p_dir = interim_dir / 'chirps_preprocessed' / 'chirps_19812019_kenya.nc'
   ...: v_dir = interim_dir / 'vhi_preprocessed' / 'vhi_preprocess_kenya.nc'
   ...:
   ...: e = EventDetector(p_dir)
   ...: e.detect(
   ...:     variable='precip', time_period='dayofyear', hilo='low', method='std'
   ...: )
   ...:
WARNING:py.warnings:/soge-home/users/chri4118/.conda/envs/crp/lib/python3.7/site-packages/xarray/core/options.py:50: FutureWarning: The enable_cftimeindex option is now a no-op and will be removed in a future version of xarray.
  FutureWarning)

chirps_19812019_kenya.nc read!
Detecting precip exceedences (low) of threshold:         std. The threshold is unique for each dayofyear
WARNING:py.warnings:/soge-home/users/chri4118/.conda/envs/crp/lib/python3.7/site-packages/xarray/core/nanops.py:159: RuntimeWarning: Mean of empty slice
  return np.nanmean(a, axis=axis, dtype=dtype)

Calculated climatology (mean for each dayofyear) - `clim`
Calculated threshold - `thresh`
WARNING:py.warnings:/soge-home/projects/crop_yield/ml_drought/scripts/eng_utils.py:391: UserWarning: assumes that all of the null values from the HOLAPS file are valid null values (e.g. water bodies). Could also be invalid nulls due to poor data processing / lack of satellite input data for a pixel!
  warnings.warn('assumes that all of the null values from the HOLAPS file are valid null values (e.g. water bodies). Could also be invalid nulls due to poor data processing / lack of satellite input data for a pixel!')

WARNING:py.warnings:/soge-home/projects/crop_yield/ml_drought/scripts/eng_utils.py:392: UserWarning: How to collapse the time dimension in the holaps mask? Here we just select the first time because all of the valid pixels are constant for first, last second last. Need to check this is true for all timesteps
  warnings.warn('How to collapse the time dimension in the holaps mask? Here we just select the first time because all of the valid pixels are constant for first, last second last. Need to check this is true for all timesteps')

** exceedences calculated **

In [2]: e
Out[2]: <src.analysis.event_detector.EventDetector at 0x7f73bc17f898>

In [3]: e.ds
Out[3]:
<xarray.Dataset>
Dimensions:    (latitude: 224, longitude: 176, time: 2760)
Coordinates:
  * latitude   (latitude) float32 -5.175003 -5.125 ... 5.924999 5.9749985
  * longitude  (longitude) float32 33.524994 33.574997 ... 42.22499 42.274994
  * time       (time) datetime64[ns] 1981-01-01 1981-01-06 ... 2019-04-26
Data variables:
    precip     (time, latitude, longitude) float32 19.300377 ... 7.446422
    dayofyear  (time) int64 1 6 11 16 21 26 32 37 ... 85 91 96 101 106 111 116
Attributes:
    Conventions:       CF-1.6
    title:             CHIRPS Version 2.0
    history:           created by Climate Hazards Group
    version:           Version 2.0
    date_created:      2015-02-27
    creator_name:      Pete Peterson
    creator_email:     pete@geog.ucsb.edu
    institution:       Climate Hazards Group.  University of California at Sa...
    documentation:     http://pubs.usgs.gov/ds/832/
    reference:         Funk, C.C., Peterson, P.J., Landsfeld, M.F., Pedreros,...
    comments:           time variable denotes the first day of the given pentad.
    acknowledgements:  The Climate Hazards Group InfraRed Precipitation with ...
    ftp_url:           ftp://chg-ftpout.geog.ucsb.edu/pub/org/chg/products/CH...
    website:           http://chg.geog.ucsb.edu/data/chirps/index.html
    faq:               http://chg-wiki.geog.ucsb.edu/wiki/CHIRPS_FAQ

In [4]: e.clim
Out[4]:
<xarray.Dataset>
Dimensions:    (dayofyear: 132, latitude: 224, longitude: 176)
Coordinates:
  * latitude   (latitude) float32 -5.175003 -5.125 ... 5.924999 5.9749985
  * longitude  (longitude) float32 33.524994 33.574997 ... 42.22499 42.274994
  * dayofyear  (dayofyear) int64 1 6 11 16 21 26 32 ... 350 351 355 356 360 361
Data variables:
    precip     (dayofyear, latitude, longitude) float32 24.191504 ... 0.66290134

In [5]: e.thresh
Out[5]:
<xarray.Dataset>
Dimensions:    (dayofyear: 132, latitude: 224, longitude: 176)
Coordinates:
  * latitude   (latitude) float32 -5.175003 -5.125 ... 5.924999 5.9749985
  * longitude  (longitude) float32 33.524994 33.574997 ... 42.22499 42.274994
  * dayofyear  (dayofyear) int64 1 6 11 16 21 26 32 ... 350 351 355 356 360 361
Data variables:
    precip     (dayofyear, latitude, longitude) float32 10.857923 ... 0.5121623

In [6]: e.exceedences
Out[6]:
<xarray.DataArray 'precip' (time: 2760, latitude: 224, longitude: 176)>
array([[[False, False, ...,  True,  True],
        [False, False, ...,  True,  True],
        ...,
        [False, False, ...,  True,  True],
        [False, False, ...,  True,  True]],

       [[False, False, ...,  True,  True],
        [False, False, ...,  True,  True],
        ...,
        [False, False, ...,  True,  True],
        [False, False, ...,  True,  True]],

       ...,

       [[False, False, ...,  True,  True],
        [False, False, ...,  True,  True],
        ...,
        [False, False, ..., False, False],
        [False, False, ..., False, False]],

       [[False, False, ...,  True,  True],
        [False, False, ...,  True,  True],
        ...,
        [False, False, ...,  True,  True],
        [False, False, ...,  True,  True]]])
Coordinates:
  * latitude   (latitude) float32 -5.175003 -5.125 ... 5.924999 5.9749985
  * longitude  (longitude) float32 33.524994 33.574997 ... 42.22499 42.274994
  * time       (time) datetime64[ns] 1981-01-01 1981-01-06 ... 2019-04-26

In [7]: runs = e.calculate_runs()
WARNING:py.warnings:/soge-home/users/chri4118/.conda/envs/crp/lib/python3.7/site-packages/dask/array/blockwise.py:204: UserWarning: The da.atop function has moved to da.blockwise
  warnings.warn("The da.atop function has moved to da.blockwise")

In [8]: runs
Out[8]:
<xarray.DataArray (time: 2761, latitude: 224, longitude: 176)>
array([[[0.00e+00, 0.00e+00, ..., 2.76e+03, 2.76e+03],
        [0.00e+00, 0.00e+00, ..., 2.76e+03, 2.76e+03],
        ...,
        [0.00e+00, 0.00e+00, ..., 2.00e+00, 2.00e+00],
        [0.00e+00, 0.00e+00, ..., 2.00e+00, 2.00e+00]],

       [[0.00e+00, 0.00e+00, ...,      nan,      nan],
        [0.00e+00, 0.00e+00, ...,      nan,      nan],
        ...,
        [0.00e+00, 0.00e+00, ...,      nan,      nan],
        [0.00e+00, 0.00e+00, ...,      nan,      nan]],

       ...,

       [[0.00e+00, 0.00e+00, ...,      nan,      nan],
        [0.00e+00, 0.00e+00, ...,      nan,      nan],
        ...,
        [0.00e+00, 0.00e+00, ..., 1.00e+00, 1.00e+00],
        [0.00e+00, 0.00e+00, ..., 1.00e+00, 1.00e+00]],

       [[0.00e+00, 0.00e+00, ...,      nan,      nan],
        [0.00e+00, 0.00e+00, ...,      nan,      nan],
        ...,
        [0.00e+00, 0.00e+00, ...,      nan,      nan],
        [0.00e+00, 0.00e+00, ...,      nan,      nan]]])
Coordinates:
  * latitude   (latitude) float32 -5.175003 -5.125 ... 5.924999 5.9749985
  * longitude  (longitude) float32 33.524994 33.574997 ... 42.22499 42.274994
  * time       (time) datetime64[ns] 1981-01-01 1981-01-06 ... 2019-04-26

In [9]: runs.max()
Out[9]:
<xarray.DataArray ()>
array(2760.)

[tommylees112]

Something's not working properly because the runs are just calculating the number of timesteps ... that's going to be a really annoying bug i can tell already ...

In [9]: runs.max()
Out[9]:
<xarray.DataArray ()>
array(2760.)

[tommylees112]

Problems with working on environments:

Opened an xclim issue

[tommylees112]

WRITE TESTS!