Information content example: compute knees

%load_ext autoreload
%autoreload
import fcpy

ds = fcpy.open_dataset("../data/cams_q_20191201_v3.nc")

ds = ds[["q"]].isel(lev=range(0, 10))
ds

<xarray.Dataset>
Dimensions:  (time: 1, lev: 10, lat: 451, lon: 900)
Coordinates:
  * time     (time) datetime64[ns] 2019-12-01T12:00:00
  * lon      (lon) float64 0.0 0.4 0.8 1.2 1.6 ... 358.0 358.4 358.8 359.2 359.6
  * lat      (lat) float64 -90.0 -89.6 -89.2 -88.8 -88.4 ... 88.8 89.2 89.6 90.0
  * lev      (lev) float64 1.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0 9.0 10.0
Data variables:
    q        (time, lev, lat, lon) float32 dask.array<chunksize=(1, 10, 451, 900), meta=np.ndarray>
Attributes:
    CDI:          Climate Data Interface version 1.9.8 (https://mpimet.mpg.de...
    Conventions:  CF-1.6
    history:      Mon Feb 15 18:37:28 2021: cdo -f nc4 copy tmp.grib data/mil...
    institution:  European Centre for Medium-Range Weather Forecasts
    CDO:          Climate Data Operators version 1.9.8 (https://mpimet.mpg.de...

xarray.Dataset

• Dimensions:
  • time: 1
  • lev: 10
  • lat: 451
  • lon: 900
• Coordinates: (4)
  • time
    (time)
    datetime64[ns]
    2019-12-01T12:00:00

    standard_name :

    time

    axis :

    T

    array(['2019-12-01T12:00:00.000000000'], dtype='datetime64[ns]')

  • lon
    (lon)
    float64
    0.0 0.4 0.8 ... 358.8 359.2 359.6

    standard_name :

    longitude

    long_name :

    longitude

    units :

    degrees_east

    axis :

    X

    array([  0. ,   0.4,   0.8, ..., 358.8, 359.2, 359.6])

  • lat
    (lat)
    float64
    -90.0 -89.6 -89.2 ... 89.6 90.0

    standard_name :

    latitude

    long_name :

    latitude

    units :

    degrees_north

    axis :

    Y

    array([-90. , -89.6, -89.2, ...,  89.2,  89.6,  90. ])

  • lev
    (lev)
    float64
    1.0 2.0 3.0 4.0 ... 8.0 9.0 10.0

    standard_name :

    hybrid_sigma_pressure

    long_name :

    hybrid level at layer midpoints

    formula :

    hyam hybm (mlev=hyam+hybm*aps)

    formula_terms :

    ap: hyam b: hybm ps: aps

    units :

    level

    positive :

    down

    array([ 1.,  2.,  3.,  4.,  5.,  6.,  7.,  8.,  9., 10.])

• Data variables: (1)
  • q
    (time, lev, lat, lon)
    float32
    dask.array<chunksize=(1, 10, 451, 900), meta=np.ndarray>

    standard_name :

    specific_humidity

    long_name :

    Specific humidity

    units :

    kg kg**-1

    param :

    0.1.0

    institution :

    ECMWF

    Array Chunk Bytes 15.48 MiB 15.48 MiB Shape (1, 10, 451, 900) (1, 10, 451, 900) Dask graph 1 chunks in 3 graph layers Data type float32 numpy.ndarray

• Indexes: (4)
  • time
    PandasIndex

    PandasIndex(DatetimeIndex(['2019-12-01 12:00:00'], dtype='datetime64[ns]', name='time', freq=None))

  • lon
    PandasIndex

    PandasIndex(Index([               0.0,                0.4,                0.8,
           1.2000000000000002,                1.6,                2.0,
           2.4000000000000004, 2.8000000000000003,                3.2,
                          3.6,
           ...
                        356.0, 356.40000000000003,              356.8,
           357.20000000000005,              357.6,              358.0,
           358.40000000000003,              358.8, 359.20000000000005,
                        359.6],
          dtype='float64', name='lon', length=900))

  • lat
    PandasIndex

    PandasIndex(Index([            -90.0,             -89.6,             -89.2,
                       -88.8,             -88.4,             -88.0,
                       -87.6,             -87.2,             -86.8,
                       -86.4,
           ...
                        86.4, 86.80000000000001, 87.20000000000002,
           87.60000000000002,              88.0,              88.4,
           88.80000000000001, 89.20000000000002, 89.60000000000002,
                        90.0],
          dtype='float64', name='lat', length=451))

  • lev
    PandasIndex

    PandasIndex(Index([1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0, 10.0], dtype='float64', name='lev'))

• Attributes: (5)

  CDI :

  Climate Data Interface version 1.9.8 (https://mpimet.mpg.de/cdi)

  Conventions :

  CF-1.6

  history :

  Mon Feb 15 18:37:28 2021: cdo -f nc4 copy tmp.grib data/milan/gg/cams_q_20191201.grib_v3.nc

  institution :

  European Centre for Medium-Range Weather Forecasts

  CDO :

  Climate Data Operators version 1.9.8 (https://mpimet.mpg.de/cdo)

da_reference = ds["q"].isel(lev=0)
da_reference.plot()

<matplotlib.collections.QuadMesh at 0x7f43d8bcfe50>

bits_params = fcpy.get_bits_params(da_reference)
bits_min = fcpy.compute_min_bits(da_reference, bits_params)
da_sigmas = fcpy.compute_sigmas(
    da_reference, compressors=[fcpy.Round(), fcpy.LinQuantization()]
)
fcpy.compute_knees_field(da_reference, da_sigmas, plot=True);

Compressor: Round

  0%|          | 0/29 [00:00<?, ?it/s]  Activating project at `~/work/field-compression/field-compression`
WARNING: method definition for == at /usr/share/miniconda/envs/fcpy/share/julia/packages/ChainRulesCore/ctmSK/src/tangent_types/tangent.jl:68 declares type variable T but does not use it.
WARNING: method definition for getindex at /usr/share/miniconda/envs/fcpy/share/julia/packages/ChainRulesCore/ctmSK/src/tangent_types/tangent.jl:120 declares type variable T but does not use it.
WARNING: method definition for getindex at /usr/share/miniconda/envs/fcpy/share/julia/packages/ChainRulesCore/ctmSK/src/tangent_types/tangent.jl:120 declares type variable P but does not use it.
WARNING: method definition for canonicalize at /usr/share/miniconda/envs/fcpy/share/julia/packages/ChainRulesCore/ctmSK/src/tangent_types/tangent.jl:240 declares type variable L but does not use it.
WARNING: method definition for canonicalize at /usr/share/miniconda/envs/fcpy/share/julia/packages/ChainRulesCore/ctmSK/src/tangent_types/tangent.jl:241 declares type variable L but does not use it.
 76%|███████▌  | 22/29 [01:04<00:04,  1.73it/s]/usr/share/miniconda/envs/fcpy/lib/python3.10/site-packages/dask/core.py:119: RuntimeWarning: invalid value encountered in divide
  return func(*(_execute_task(a, cache) for a in args))
 79%|███████▉  | 23/29 [01:05<00:03,  1.73it/s]/usr/share/miniconda/envs/fcpy/lib/python3.10/site-packages/dask/core.py:119: RuntimeWarning: invalid value encountered in divide
  return func(*(_execute_task(a, cache) for a in args))
 83%|████████▎ | 24/29 [01:05<00:02,  1.76it/s]/usr/share/miniconda/envs/fcpy/lib/python3.10/site-packages/dask/core.py:119: RuntimeWarning: invalid value encountered in divide
  return func(*(_execute_task(a, cache) for a in args))
 86%|████████▌ | 25/29 [01:06<00:02,  1.77it/s]/usr/share/miniconda/envs/fcpy/lib/python3.10/site-packages/dask/core.py:119: RuntimeWarning: invalid value encountered in divide
  return func(*(_execute_task(a, cache) for a in args))
 90%|████████▉ | 26/29 [01:06<00:01,  1.78it/s]/usr/share/miniconda/envs/fcpy/lib/python3.10/site-packages/dask/core.py:119: RuntimeWarning: invalid value encountered in divide
  return func(*(_execute_task(a, cache) for a in args))
 93%|█████████▎| 27/29 [01:07<00:01,  1.78it/s]/usr/share/miniconda/envs/fcpy/lib/python3.10/site-packages/dask/core.py:119: RuntimeWarning: invalid value encountered in divide
  return func(*(_execute_task(a, cache) for a in args))
 97%|█████████▋| 28/29 [01:07<00:00,  1.78it/s]/usr/share/miniconda/envs/fcpy/lib/python3.10/site-packages/dask/core.py:119: RuntimeWarning: invalid value encountered in divide
  return func(*(_execute_task(a, cache) for a in args))
100%|██████████| 29/29 [01:08<00:00,  2.36s/it]

Compressor: LinQuantization

100%|██████████| 29/29 [00:13<00:00,  2.08it/s]