Summary

This example is based on the OceanWatch tutorial meterial edited with Great Lakes satellite data.In this example you will see how to extract Great Lakes ice concentration data from the ERDDAP server and make a ice concentration map, and caculate the monthly

The example demonstrates the following techniques:

• Loading Great Lakes ice concentration data from Great Lakes ERDDAP data server.
• Create a map of ice concentration.
• Compute the daily mean over the selected region.

Datesets used:

• Great Lakes ice concentration: Great Lakes ice concentration product.

Import the required Python modules

import xarray as xr
import pandas as pd
import numpy as np
import urllib.request

from matplotlib import pyplot as plt
from matplotlib.colors import LinearSegmentedColormap
import warnings
#warnings.filterwarnings('ignore')

Downlading data from ERDDAP server

Because ERDDAP includes RESTful services, you can download data listed on any ERDDAP platform from Python using the URL structure. For example, the following link allows you to subset daily ice concentration data from the dataset GL_Ice_Concentration_GCS

In this specific example, we will get the SST data from 2023-06-01 to 2023-06-30. the URL we generated is :

https://apps.glerl.noaa.gov/erddap/griddap/GL_Ice_Concentration_GCS.nc?ice_concentration%5B(1995-01-01T12:00:00Z):1:(2024-05-01T12:00:00Z)%5D%5B(41.38):1:(42.10)%5D%5B(-83.59):1:(-82.5)%5D

we extract the data in csv format due to the nc library not available.

url="https://apps.glerl.noaa.gov/erddap/griddap/GL_Ice_Concentration_GCS.nc?ice_concentration%5B(1995-01-01T12:00:00Z):1:(2024-05-01T12:00:00Z)%5D%5B(41.38):1:(42.10)%5D%5B(-83.59):1:(-82.5)%5D"
urllib.request.urlretrieve(url, "w_e_ice_concentration.nc")

('w_e_ice_concentration.nc', <http.client.HTTPMessage at 0x1b687bf22d0>)

Importing NetCDF4 data in Python

Now that we’ve downloaded the data locally, we can import it and extract our variables of interest.

The xarray package makes it very convenient to work with NetCDF files. Documentation is available here: http://xarray.pydata.org/en/stable/why-xarray.html

import xarray as xr
import netCDF4 as nc

Open the file and load it as an xarray dataset:

ds = xr.open_dataset('w_e_ice_concentration.nc',decode_cf=False)
#ds = xr.open_dataset('e_sst.nc')

Examine the data structure:

ds

<xarray.Dataset>
Dimensions:            (time: 7146, latitude: 52, longitude: 79)
Coordinates:
  * time               (time) float64 7.89e+08 7.89e+08 ... 1.714e+09 1.715e+09
  * latitude           (latitude) float64 41.38 41.4 41.41 ... 42.07 42.08 42.1
  * longitude          (longitude) float64 -83.59 -83.58 -83.56 ... -82.51 -82.5
Data variables:
    ice_concentration  (time, latitude, longitude) float32 ...
Attributes: (12/28)
    cdm_data_type:              Grid
    Conventions:                CF-1.6, COARDS, ACDD-1.3
    Easternmost_Easting:        -82.496964466524
    GDAL:                       GDAL 3.4.3, released 2022/04/22
    geospatial_lat_max:         42.0985561800016
    geospatial_lat_min:         41.3837647963109
    ...                         ...
    standard_name_vocabulary:   CF Standard Name Table v29
    summary:                    Ice Concentration from Great Lakes Surface En...
    time_coverage_end:          2024-05-01T12:00:00Z
    time_coverage_start:        1995-01-01T12:00:00Z
    title:                      Ice Concentration from Great Lakes Surface En...
    Westernmost_Easting:        -83.590174818051

xarray.Dataset

• Dimensions:
  • time: 7146
  • latitude: 52
  • longitude: 79
• Coordinates: (3)
  • time
    (time)
    float64
    7.89e+08 7.89e+08 ... 1.715e+09

    _CoordinateAxisType :

    Time

    actual_range :

    [7.8896160e+08 1.7145648e+09]

    axis :

    T

    calendar :

    Gregorian

    ioos_category :

    Time

    long_name :

    Time

    standard_name :

    time

    time_origin :

    01-JAN-1970 00:00:00

    units :

    seconds since 1970-01-01T00:00:00Z

    array([7.889616e+08, 7.890480e+08, 7.891344e+08, ..., 1.714392e+09,
           1.714478e+09, 1.714565e+09])

  • latitude
    (latitude)
    float64
    41.38 41.4 41.41 ... 42.08 42.1

    _CoordinateAxisType :

    Lat

    actual_range :

    [41.3837648 42.09855618]

    axis :

    Y

    ioos_category :

    Location

    long_name :

    Latitude

    standard_name :

    latitude

    units :

    degrees_north

    array([41.383765, 41.39778 , 41.411796, 41.425811, 41.439827, 41.453842,
           41.467858, 41.481873, 41.495889, 41.509904, 41.52392 , 41.537935,
           41.551951, 41.565967, 41.579982, 41.593998, 41.608013, 41.622029,
           41.636044, 41.65006 , 41.664075, 41.678091, 41.692106, 41.706122,
           41.720137, 41.734153, 41.748168, 41.762184, 41.776199, 41.790215,
           41.80423 , 41.818246, 41.832261, 41.846277, 41.860292, 41.874308,
           41.888323, 41.902339, 41.916354, 41.93037 , 41.944385, 41.958401,
           41.972417, 41.986432, 42.000448, 42.014463, 42.028479, 42.042494,
           42.05651 , 42.070525, 42.084541, 42.098556])

  • longitude
    (longitude)
    float64
    -83.59 -83.58 ... -82.51 -82.5

    _CoordinateAxisType :

    Lon

    actual_range :

    [-83.59017482 -82.49696447]

    axis :

    X

    ioos_category :

    Location

    long_name :

    Longitude

    standard_name :

    longitude

    units :

    degrees_east

    array([-83.590175, -83.576159, -83.562144, -83.548128, -83.534113, -83.520097,
           -83.506082, -83.492066, -83.478051, -83.464035, -83.45002 , -83.436004,
           -83.421989, -83.407973, -83.393958, -83.379942, -83.365927, -83.351911,
           -83.337896, -83.32388 , -83.309864, -83.295849, -83.281833, -83.267818,
           -83.253802, -83.239787, -83.225771, -83.211756, -83.19774 , -83.183725,
           -83.169709, -83.155694, -83.141678, -83.127663, -83.113647, -83.099632,
           -83.085616, -83.071601, -83.057585, -83.04357 , -83.029554, -83.015539,
           -83.001523, -82.987508, -82.973492, -82.959477, -82.945461, -82.931446,
           -82.91743 , -82.903414, -82.889399, -82.875383, -82.861368, -82.847352,
           -82.833337, -82.819321, -82.805306, -82.79129 , -82.777275, -82.763259,
           -82.749244, -82.735228, -82.721213, -82.707197, -82.693182, -82.679166,
           -82.665151, -82.651135, -82.63712 , -82.623104, -82.609089, -82.595073,
           -82.581058, -82.567042, -82.553027, -82.539011, -82.524996, -82.51098 ,
           -82.496964])

• Data variables: (1)
  • ice_concentration
    (time, latitude, longitude)
    float32
    ...

    _FillValue :

    -99999.0

    colorBarMaximum :

    100.0

    colorBarMinimum :

    0.0

    colorBarPalette :

    WhiteBlack

    grid_mapping :

    crs

    ioos_category :

    Ocean Color

    long_name :

    Ice Concentration

    standard_name :

    ice_concentration

    units :

    percent

    [29355768 values with dtype=float32]

• Indexes: (3)
  • time
    PandasIndex

    PandasIndex(Index([ 788961600.0,  789048000.0,  789134400.0,  789220800.0,  789307200.0,
            789393600.0,  789480000.0,  789566400.0,  789652800.0,  789739200.0,
           ...
           1713787200.0, 1713873600.0, 1713960000.0, 1714046400.0, 1714132800.0,
           1714219200.0, 1714305600.0, 1714392000.0, 1714478400.0, 1714564800.0],
          dtype='float64', name='time', length=7146))

  • latitude
    PandasIndex

    PandasIndex(Index([41.3837647963109, 41.3977803136382, 41.4117958309654, 41.4258113482927,
             41.43982686562, 41.4538423829472, 41.4678579002745, 41.4818734176018,
            41.495888934929, 41.5099044522563, 41.5239199695836, 41.5379354869108,
           41.5519510042381, 41.5659665215654, 41.5799820388927, 41.5939975562199,
           41.6080130735472, 41.6220285908745, 41.6360441082017,  41.650059625529,
           41.6640751428563, 41.6780906601835, 41.6921061775108, 41.7061216948381,
           41.7201372121653, 41.7341527294926, 41.7481682468199, 41.7621837641471,
           41.7761992814744, 41.7902147988017,  41.804230316129, 41.8182458334562,
           41.8322613507835, 41.8462768681108,  41.860292385438, 41.8743079027653,
           41.8883234200926, 41.9023389374198, 41.9163544547471, 41.9303699720744,
           41.9443854894016, 41.9584010067289, 41.9724165240562, 41.9864320413835,
           42.0004475587107,  42.014463076038, 42.0284785933653, 42.0424941106925,
           42.0565096280198, 42.0705251453471, 42.0845406626743, 42.0985561800016],
          dtype='float64', name='latitude'))

  • longitude
    PandasIndex

    PandasIndex(Index([ -83.590174818051, -83.5761593007237, -83.5621437833965,
           -83.5481282660692, -83.5341127487419, -83.5200972314146,
           -83.5060817140874, -83.4920661967601, -83.4780506794328,
           -83.4640351621056, -83.4500196447783,  -83.436004127451,
           -83.4219886101238, -83.4079730927965, -83.3939575754692,
           -83.3799420581419, -83.3659265408147, -83.3519110234874,
           -83.3378955061601, -83.3238799888329, -83.3098644715056,
           -83.2958489541783, -83.2818334368511, -83.2678179195238,
           -83.2538024021965, -83.2397868848693,  -83.225771367542,
           -83.2117558502147, -83.1977403328874, -83.1837248155602,
           -83.1697092982329, -83.1556937809057, -83.1416782635784,
           -83.1276627462511, -83.1136472289238, -83.0996317115966,
           -83.0856161942693,  -83.071600676942, -83.0575851596148,
           -83.0435696422875, -83.0295541249602,  -83.015538607633,
           -83.0015230903057, -82.9875075729784, -82.9734920556511,
           -82.9594765383239, -82.9454610209966, -82.9314455036693,
           -82.9174299863421, -82.9034144690148, -82.8893989516875,
           -82.8753834343603,  -82.861367917033, -82.8473523997057,
           -82.8333368823785, -82.8193213650512, -82.8053058477239,
           -82.7912903303966, -82.7772748130694, -82.7632592957421,
           -82.7492437784149, -82.7352282610876, -82.7212127437603,
            -82.707197226433, -82.6931817091058, -82.6791661917785,
           -82.6651506744512,  -82.651135157124, -82.6371196397967,
           -82.6231041224694, -82.6090886051422, -82.5950730878149,
           -82.5810575704876, -82.5670420531603, -82.5530265358331,
           -82.5390110185058, -82.5249955011785, -82.5109799838513,
            -82.496964466524],
          dtype='float64', name='longitude'))

• Attributes: (28)

  cdm_data_type :

  Grid

  Conventions :

  CF-1.6, COARDS, ACDD-1.3

  Easternmost_Easting :

  -82.496964466524

  GDAL :

  GDAL 3.4.3, released 2022/04/22

  geospatial_lat_max :

  42.0985561800016

  geospatial_lat_min :

  41.3837647963109

  geospatial_lat_resolution :

  0.014015517327269056

  geospatial_lat_units :

  degrees_north

  geospatial_lon_max :

  -82.496964466524

  geospatial_lon_min :

  -83.590174818051

  geospatial_lon_resolution :

  0.01401551732726889

  geospatial_lon_units :

  degrees_east

  history :

  Ice concentration from Great Lakes Surface Environmental Analysis (GLSEA) asc format to nc fromat 2024-09-18T16:53:43Z (local files) 2024-09-18T16:53:43Z https://apps.glerl.noaa.gov/erddap/griddap/GL_Ice_Concentration_GCS.nc?ice_concentration%5B(1995-01-01T12:00:00Z):1:(2024-05-01T12:00:00Z)%5D%5B(41.38):1:(42.10)%5D%5B(-83.59):1:(-82.5)%5D

  infoUrl :

  https://coastwatch.glerl.noaa.gov/glsea/glsea.html

  institution :

  CoastWatch Great Lakes Node

  keywords :

  analysis, concentration, data, distribution, environmental, glsea, great, great lakes, ice, ice distribution, ice_concentration, lakes, nic, surface, time

  license :

  The data may be used and redistributed for free but is not intended for legal use, since it may contain inaccuracies. Neither the data Contributor, ERD, NOAA, nor the United States Government, nor any of their employees or contractors, makes any warranty, express or implied, including warranties of merchantability and fitness for a particular purpose, or assumes any legal liability for the accuracy, completeness, or usefulness, of this information.

  NCO :

  netCDF Operators version 5.0.7 (Homepage = http://nco.sf.net, Code = https://github.com/nco/nco)

  Northernmost_Northing :

  42.0985561800016

  source :

  Ice Concentration from Great Lakes Surface Environmental Analysis (GLSEA) and NIC, Geodetic coordinate system (LAT, LON), 1995-present

  sourceUrl :

  (local files)

  Southernmost_Northing :

  41.3837647963109

  standard_name_vocabulary :

  CF Standard Name Table v29

  summary :

  Ice Concentration from Great Lakes Surface Environmental Analysis (GLSEA) and NIC

  time_coverage_end :

  2024-05-01T12:00:00Z

  time_coverage_start :

  1995-01-01T12:00:00Z

  title :

  Ice Concentration from Great Lakes Surface Environmental Analysis (GLSEA) and NIC, Geodetic coordinate system (LAT, LON), 1995-present

  Westernmost_Easting :

  -83.590174818051

Examine which coordinates and variables are included in the dataset:

 ds.dims

Frozen({'time': 7146, 'latitude': 52, 'longitude': 79})

ds.coords

Coordinates:
  * time       (time) float64 7.89e+08 7.89e+08 ... 1.714e+09 1.715e+09
  * latitude   (latitude) float64 41.38 41.4 41.41 41.43 ... 42.07 42.08 42.1
  * longitude  (longitude) float64 -83.59 -83.58 -83.56 ... -82.52 -82.51 -82.5

ds.data_vars

Data variables:
    ice_concentration  (time, latitude, longitude) float32 ...

ds.attrs

{'cdm_data_type': 'Grid',
 'Conventions': 'CF-1.6, COARDS, ACDD-1.3',
 'Easternmost_Easting': -82.496964466524,
 'GDAL': 'GDAL 3.4.3, released 2022/04/22',
 'geospatial_lat_max': 42.0985561800016,
 'geospatial_lat_min': 41.3837647963109,
 'geospatial_lat_resolution': 0.014015517327269056,
 'geospatial_lat_units': 'degrees_north',
 'geospatial_lon_max': -82.496964466524,
 'geospatial_lon_min': -83.590174818051,
 'geospatial_lon_resolution': 0.01401551732726889,
 'geospatial_lon_units': 'degrees_east',
 'history': 'Ice concentration from Great Lakes Surface Environmental Analysis (GLSEA) asc format to nc fromat\n2024-09-18T16:53:43Z (local files)\n2024-09-18T16:53:43Z https://apps.glerl.noaa.gov/erddap/griddap/GL_Ice_Concentration_GCS.nc?ice_concentration%5B(1995-01-01T12:00:00Z):1:(2024-05-01T12:00:00Z)%5D%5B(41.38):1:(42.10)%5D%5B(-83.59):1:(-82.5)%5D',
 'infoUrl': 'https://coastwatch.glerl.noaa.gov/glsea/glsea.html',
 'institution': 'CoastWatch Great Lakes Node',
 'keywords': 'analysis, concentration, data, distribution, environmental, glsea, great, great lakes, ice, ice distribution, ice_concentration, lakes, nic, surface, time',
 'license': 'The data may be used and redistributed for free but is not intended\nfor legal use, since it may contain inaccuracies. Neither the data\nContributor, ERD, NOAA, nor the United States Government, nor any\nof their employees or contractors, makes any warranty, express or\nimplied, including warranties of merchantability and fitness for a\nparticular purpose, or assumes any legal liability for the accuracy,\ncompleteness, or usefulness, of this information.',
 'NCO': 'netCDF Operators version 5.0.7 (Homepage = http://nco.sf.net, Code = https://github.com/nco/nco)',
 'Northernmost_Northing': 42.0985561800016,
 'source': 'Ice Concentration from Great Lakes Surface Environmental Analysis (GLSEA) and NIC, Geodetic coordinate system (LAT, LON), 1995-present',
 'sourceUrl': '(local files)',
 'Southernmost_Northing': 41.3837647963109,
 'standard_name_vocabulary': 'CF Standard Name Table v29',
 'summary': 'Ice Concentration from Great Lakes Surface Environmental Analysis (GLSEA) and NIC',
 'time_coverage_end': '2024-05-01T12:00:00Z',
 'time_coverage_start': '1995-01-01T12:00:00Z',
 'title': 'Ice Concentration from Great Lakes Surface Environmental Analysis (GLSEA) and NIC, Geodetic coordinate system (LAT, LON), 1995-present',
 'Westernmost_Easting': -83.590174818051}

Examine the structure of ice concentration:

ds.ice_concentration.shape

(7146, 52, 79)

Our dataset is a 3-D array with 52 rows corresponding to latitudes and 79 columns corresponding to longitudes, for each of the 7146 time steps. #### Get the dates for each time step:

ds.time

<xarray.DataArray 'time' (time: 7146)>
array([7.889616e+08, 7.890480e+08, 7.891344e+08, ..., 1.714392e+09,
       1.714478e+09, 1.714565e+09])
Coordinates:
  * time     (time) float64 7.89e+08 7.89e+08 7.891e+08 ... 1.714e+09 1.715e+09
Attributes:
    _CoordinateAxisType:  Time
    actual_range:         [7.8896160e+08 1.7145648e+09]
    axis:                 T
    calendar:             Gregorian
    ioos_category:        Time
    long_name:            Time
    standard_name:        time
    time_origin:          01-JAN-1970 00:00:00
    units:                seconds since 1970-01-01T00:00:00Z

xarray.DataArray

'time'

• time: 7146

• 7.89e+08 7.89e+08 7.891e+08 ... 1.714e+09 1.714e+09 1.715e+09

  array([7.889616e+08, 7.890480e+08, 7.891344e+08, ..., 1.714392e+09,
         1.714478e+09, 1.714565e+09])

• Coordinates: (1)
  • time
    (time)
    float64
    7.89e+08 7.89e+08 ... 1.715e+09

    _CoordinateAxisType :

    Time

    actual_range :

    [7.8896160e+08 1.7145648e+09]

    axis :

    T

    calendar :

    Gregorian

    ioos_category :

    Time

    long_name :

    Time

    standard_name :

    time

    time_origin :

    01-JAN-1970 00:00:00

    units :

    seconds since 1970-01-01T00:00:00Z

    array([7.889616e+08, 7.890480e+08, 7.891344e+08, ..., 1.714392e+09,
           1.714478e+09, 1.714565e+09])

• Indexes: (1)
  • time
    PandasIndex

    PandasIndex(Index([ 788961600.0,  789048000.0,  789134400.0,  789220800.0,  789307200.0,
            789393600.0,  789480000.0,  789566400.0,  789652800.0,  789739200.0,
           ...
           1713787200.0, 1713873600.0, 1713960000.0, 1714046400.0, 1714132800.0,
           1714219200.0, 1714305600.0, 1714392000.0, 1714478400.0, 1714564800.0],
          dtype='float64', name='time', length=7146))

• Attributes: (9)

  _CoordinateAxisType :

  Time

  actual_range :

  [7.8896160e+08 1.7145648e+09]

  axis :

  T

  calendar :

  Gregorian

  ioos_category :

  Time

  long_name :

  Time

  standard_name :

  time

  time_origin :

  01-JAN-1970 00:00:00

  units :

  seconds since 1970-01-01T00:00:00Z

ds.time.attrs

{'_CoordinateAxisType': 'Time',
 'actual_range': array([7.8896160e+08, 1.7145648e+09]),
 'axis': 'T',
 'calendar': 'Gregorian',
 'ioos_category': 'Time',
 'long_name': 'Time',
 'standard_name': 'time',
 'time_origin': '01-JAN-1970 00:00:00',
 'units': 'seconds since 1970-01-01T00:00:00Z'}

print(ds.time)

<xarray.DataArray 'time' (time: 7146)>
array([7.889616e+08, 7.890480e+08, 7.891344e+08, ..., 1.714392e+09,
       1.714478e+09, 1.714565e+09])
Coordinates:
  * time     (time) float64 7.89e+08 7.89e+08 7.891e+08 ... 1.714e+09 1.715e+09
Attributes:
    _CoordinateAxisType:  Time
    actual_range:         [7.8896160e+08 1.7145648e+09]
    axis:                 T
    calendar:             Gregorian
    ioos_category:        Time
    long_name:            Time
    standard_name:        time
    time_origin:          01-JAN-1970 00:00:00
    units:                seconds since 1970-01-01T00:00:00Z

The time units is seconds, we need to convert the seconds to dates.

dates=nc.num2date(ds.time,ds.time.units,only_use_cftime_datetimes=False, 
                        only_use_python_datetimes=True )
dates

array([real_datetime(1995, 1, 1, 12, 0), real_datetime(1995, 1, 2, 12, 0),
       real_datetime(1995, 1, 3, 12, 0), ...,
       real_datetime(2024, 4, 29, 12, 0),
       real_datetime(2024, 4, 30, 12, 0),
       real_datetime(2024, 5, 1, 12, 0)], dtype=object)

Find the index of dates for 2019-03-01

for i, date in enumerate(dates):
    if date.strftime("%Y-%m-%d") == "2019-03-01":
        print(i, date)

5872 2019-03-01 12:00:00

5872 is the index on the array dates for 2019-03-01.

Create a map of ice concentration for March 1, 2019 (our 5872th time step).

import pandas as pd
import numpy as np
from matplotlib import pyplot as plt
from matplotlib.colors import LinearSegmentedColormap
#np.warnings.filterwarnings('ignore')

#### Examine the values of ice concentration:

print(ds.ice_concentration.values)
print(ds.ice_concentration.shape)

[[[-99999. -99999. -99999. ... -99999. -99999. -99999.]
  [-99999. -99999. -99999. ... -99999. -99999.      0.]
  [-99999. -99999. -99999. ...      0.      0.      0.]
  ...
  [-99999. -99999. -99999. ... -99999. -99999. -99999.]
  [-99999. -99999. -99999. ... -99999. -99999. -99999.]
  [-99999. -99999. -99999. ... -99999. -99999. -99999.]]

 [[-99999. -99999. -99999. ... -99999. -99999. -99999.]
  [-99999. -99999. -99999. ... -99999. -99999.      0.]
  [-99999. -99999. -99999. ...      0.      0.      0.]
  ...
  [-99999. -99999. -99999. ... -99999. -99999. -99999.]
  [-99999. -99999. -99999. ... -99999. -99999. -99999.]
  [-99999. -99999. -99999. ... -99999. -99999. -99999.]]

 [[-99999. -99999. -99999. ... -99999. -99999. -99999.]
  [-99999. -99999. -99999. ... -99999. -99999.      0.]
  [-99999. -99999. -99999. ...      0.      0.      0.]
  ...
  [-99999. -99999. -99999. ... -99999. -99999. -99999.]
  [-99999. -99999. -99999. ... -99999. -99999. -99999.]
  [-99999. -99999. -99999. ... -99999. -99999. -99999.]]

 ...

 [[-99999. -99999. -99999. ... -99999. -99999. -99999.]
  [-99999. -99999. -99999. ... -99999. -99999.      0.]
  [-99999. -99999. -99999. ...      0.      0.      0.]
  ...
  [-99999. -99999. -99999. ... -99999. -99999. -99999.]
  [-99999. -99999. -99999. ... -99999. -99999. -99999.]
  [-99999. -99999. -99999. ... -99999. -99999. -99999.]]

 [[-99999. -99999. -99999. ... -99999. -99999. -99999.]
  [-99999. -99999. -99999. ... -99999. -99999.      0.]
  [-99999. -99999. -99999. ...      0.      0.      0.]
  ...
  [-99999. -99999. -99999. ... -99999. -99999. -99999.]
  [-99999. -99999. -99999. ... -99999. -99999. -99999.]
  [-99999. -99999. -99999. ... -99999. -99999. -99999.]]

 [[-99999. -99999. -99999. ... -99999. -99999. -99999.]
  [-99999. -99999. -99999. ... -99999. -99999.      0.]
  [-99999. -99999. -99999. ...      0.      0.      0.]
  ...
  [-99999. -99999. -99999. ... -99999. -99999. -99999.]
  [-99999. -99999. -99999. ... -99999. -99999. -99999.]
  [-99999. -99999. -99999. ... -99999. -99999. -99999.]]]
(7146, 52, 79)

ds.ice_concentration.attrs

{'_FillValue': -99999.0,
 'colorBarMaximum': 100.0,
 'colorBarMinimum': 0.0,
 'colorBarPalette': 'WhiteBlack',
 'grid_mapping': 'crs',
 'ioos_category': 'Ocean Color',
 'long_name': 'Ice Concentration',
 'standard_name': 'ice_concentration',
 'units': 'percent'}

ds.ice_concentration.attrs['_FillValue']

-99999.0

Make a new ice concentration DataArray and replace _fillValue with NaN

nan_ice_concentration = ds.ice_concentration.where(ds.ice_concentration.values != ds.ice_concentration.attrs['_FillValue'])

print(nan_ice_concentration)

<xarray.DataArray 'ice_concentration' (time: 7146, latitude: 52, longitude: 79)>
array([[[nan, nan, nan, ..., nan, nan, nan],
        [nan, nan, nan, ..., nan, nan,  0.],
        [nan, nan, nan, ...,  0.,  0.,  0.],
        ...,
        [nan, nan, nan, ..., nan, nan, nan],
        [nan, nan, nan, ..., nan, nan, nan],
        [nan, nan, nan, ..., nan, nan, nan]],

       [[nan, nan, nan, ..., nan, nan, nan],
        [nan, nan, nan, ..., nan, nan,  0.],
        [nan, nan, nan, ...,  0.,  0.,  0.],
        ...,
        [nan, nan, nan, ..., nan, nan, nan],
        [nan, nan, nan, ..., nan, nan, nan],
        [nan, nan, nan, ..., nan, nan, nan]],

       [[nan, nan, nan, ..., nan, nan, nan],
        [nan, nan, nan, ..., nan, nan,  0.],
        [nan, nan, nan, ...,  0.,  0.,  0.],
        ...,
...
        ...,
        [nan, nan, nan, ..., nan, nan, nan],
        [nan, nan, nan, ..., nan, nan, nan],
        [nan, nan, nan, ..., nan, nan, nan]],

       [[nan, nan, nan, ..., nan, nan, nan],
        [nan, nan, nan, ..., nan, nan,  0.],
        [nan, nan, nan, ...,  0.,  0.,  0.],
        ...,
        [nan, nan, nan, ..., nan, nan, nan],
        [nan, nan, nan, ..., nan, nan, nan],
        [nan, nan, nan, ..., nan, nan, nan]],

       [[nan, nan, nan, ..., nan, nan, nan],
        [nan, nan, nan, ..., nan, nan,  0.],
        [nan, nan, nan, ...,  0.,  0.,  0.],
        ...,
        [nan, nan, nan, ..., nan, nan, nan],
        [nan, nan, nan, ..., nan, nan, nan],
        [nan, nan, nan, ..., nan, nan, nan]]], dtype=float32)
Coordinates:
  * time       (time) float64 7.89e+08 7.89e+08 ... 1.714e+09 1.715e+09
  * latitude   (latitude) float64 41.38 41.4 41.41 41.43 ... 42.07 42.08 42.1
  * longitude  (longitude) float64 -83.59 -83.58 -83.56 ... -82.52 -82.51 -82.5
Attributes:
    _FillValue:       -99999.0
    colorBarMaximum:  100.0
    colorBarMinimum:  0.0
    colorBarPalette:  WhiteBlack
    grid_mapping:     crs
    ioos_category:    Ocean Color
    long_name:        Ice Concentration
    standard_name:    ice_concentration
    units:            percent

Set some color breaks

# find min value in man_sst
np.nanmin(nan_ice_concentration)

0.0

np.nanmax(nan_ice_concentration)

99.99847

levs = np.arange(0, 101, 10)
len(levs)

11

Define a color palette

# init a color list
jet=["blue", "#007FFF", "cyan","#7FFF7F", "yellow", "#FF7F00", "red", "#7F0000"]

Set color scale using the jet palette

cm = LinearSegmentedColormap.from_list('my_jet', jet, N=len(levs))

plot the ice_concentration map

plt.subplots(figsize=(10, 5))

#plot 5872th ice concentration image: nan_ice_concentration[5872 ,:,:]
plt.contourf(nan_ice_concentration.longitude, nan_ice_concentration.latitude, nan_ice_concentration[5872,:,:], levs,cmap=cm)

#plot the color scale
plt.colorbar()

#example of how to add points to the map
#plt.scatter(np.linspace(-82,-80.5,num=4),np.repeat(42,4),c='black')

#example of how to add a contour line
#step = np.arange(9,26, 1)

#plt.contour(ds.longitude, ds.latitude, ds.sst[0,:,:],levels=step,linewidths=1)

#plot title
plt.title("West Lake Erie Ice Concentration - " + dates[5872].strftime('%b %d, %Y'))
plt.show()

Let’s compute the daily mean over the west Lake Erie region:

res=np.nanmean(nan_ice_concentration,axis=(1,2))
res

array([0., 0., 0., ..., 0., 0., 0.], dtype=float32)

Let’s plot the time-series (from 2019-03-01 to 2019-03-31):

for i, date in enumerate(dates):
    if date.strftime("%Y-%m-%d") == "2019-03-01":
        print(i, date)
    
    if date.strftime("%Y-%m-%d") == "2019-03-31":
        print(i, date)

5872 2019-03-01 12:00:00
5902 2019-03-31 12:00:00

print(res.shape)

(7146,)

plt.figure(figsize=(8,4))

plt.scatter(dates[5872:5902+1],res[5872:5902+1])

#degree_sign = u"\N{DEGREE SIGN}"
plt.ylabel('Ice Concentration (%)')

plt.xlim(dates[5872], dates[5902+1])

plt.xticks(dates[5872:5902+1],rotation=70, fontsize=10 )
plt.show()

!jupyter nbconvert --to html GL_python_tutorial4.ipynb

[NbConvertApp] Converting notebook GL_python_tutorial4.ipynb to html
[NbConvertApp] Writing 313369 bytes to GL_python_tutorial4.html