Snow Depth and Snow Cover Data Exploration¶

This tutorial demonstrates how to access and compare coincident snow data across in-situ, airborne, and satellite platforms from NASA's SnowEx, ASO, and MODIS data sets, respectively. All data are available from the NASA National Snow and Ice Data Center Distributed Active Archive Center, or NSIDC DAAC.

Here are the steps you will learn in this snow data notebook:¶

1. Explore the coverage and structure of select NSIDC DAAC snow data products, as well as available resources to search and access data.
2. Search and download spatiotemporally coincident data across in-situ, airborne, and satellite observations.
3. Subset and reformat MODIS data using the NSIDC DAAC API.
4. Read CSV and GeoTIFF formatted data using geopandas and rasterio libraries.
5. Subset data based on buffered area.
6. Extract and visualize raster values at point locations.
7. Save output as shapefile for further GIS analysis.

Explore snow products and resources¶

NSIDC introduction¶

The National Snow and Ice Data Center provides over 1100 data sets covering the Earth's cryosphere and more, all of which are available to the public free of charge. Beyond providing these data, NSIDC creates tools for data access, supports data users, performs scientific research, and educates the public about the cryosphere.

Select Data Resources¶

• NSIDC Data Search
  • Search NSIDC snow data
• NSIDC Data Update Announcements
  • News and tips for data users
• NASA Earthdata Search
  • Search and access data across the NASA Earthdata
• NASA Worldview
  • Interactive interface for browsing full-resolution, global, daily satellite images

Snow Today¶

Snow Today, a collaboration with the University of Colorado's Institute of Alpine and Arctic Research (INSTAAR), provides near-real-time snow analysis for the western United States and regular reports on conditions during the winter season. Snow Today is funded by NASA Hydrological Sciences Program and utilizes data from the Moderate Resolution Imaging Spectroradiometer (MODIS) instrument and snow station data from the Snow Telemetry (SNOTEL) network by the Natural Resources Conservation Service (NRCS), United States Department of Agriculture (USDA) and the California Department of Water Resources: www.wcc.nrcs.usda.gov/snow.

Snow-related missions and data sets used in the following steps:¶

• SnowEx
  • SnowEx17 Ground Penetrating Radar, Version 2: https://doi.org/10.5067/G21LGCNLFSC5
• ASO
  • ASO L4 Lidar Snow Depth 3m UTM Grid, Version 1: https://doi.org/10.5067/KIE9QNVG7HP0
• MODIS
  • MODIS/Terra Snow Cover Daily L3 Global 500m SIN Grid, Version 6: https://doi.org/10.5067/MODIS/MOD10A1.006

Other relevant snow products:¶

• VIIRS Snow Data

• AMSR-E and AMSR-E/AMSR2 Unified Snow Data

• MEaSUREs Snow Data

• Near-Real-Time SSM/I-SSMIS EASE-Grid Daily Global Ice Concentration and Snow Extent (NISE), Version 5: https://doi.org/10.5067/3KB2JPLFPK3R

Import Packages¶

Get started by importing packages needed to run the following code blocks, including the tutorial_helper_functions module provided within this repository.

Data Discovery¶

Start by identifying your study area and exploring coincident data over the same time and area.

NASA Earthdata Search can be used to visualize file coverage over mulitple data sets and to access the same data you will be working with below: https://search.earthdata.nasa.gov/projects?projectId=5366449248

Identify area and time of interest¶

Since our focus is on the Grand Mesa study site of the NASA SnowEx campaign, we'll use that area to search for coincident data across other data products. From the SnowEx17 Ground Penetrating Radar Version 2 landing page, you can find the rectangular spatial coverage under the Overview tab, or you can draw a polygon over your area of interest in the map under the Download Data tab and export the shape as a geojson file using the Export Polygon icon shown below. An example polygon geojson file is provided in the /Data folder of this repository.

Create polygon coordinate string¶

Read in the geojson file as a GeoDataFrame object and simplify and reorder using the shapely package. This will be converted back to a dictionary to be applied as our polygon search parameter.

Set time range¶

We are interested in accessing files within each data set over the same time range, so we'll start by searching all of 2017.

Create data dictionary¶

Create a nested dictionary with each data set shortname and version, as well as shared temporal range and polygonal area of interest. Data set shortnames, or IDs, as well as version numbers, are located at the top of every NSIDC landing page.

Determine how many files exist over this time and area of interest, as well as the average size and total volume of those files¶

We will use the granule_info function to query metadata about each data set and associated files using the Common Metadata Repository (CMR), which is a high-performance, high-quality, continuously evolving metadata system that catalogs Earth Science data and associated service metadata records. Note that not all NSIDC data can be searched at the file level using CMR, particularly those outside of the NASA DAAC program.

Find coincident data¶

The function above tells us the size of data available for each data set over our time and area of interest, but we want to go a step further and determine what time ranges are coincident based on our bounding box. This time_overlap helper function returns a dataframe with file names, dataset_id, start date, and end date for all files that overlap in temporal range across all data sets of interest.

Data Access¶

The number of files has been greatly reduced to only those needed to compare data across these data sets. This CMR query will collect the data file URLs, including the associated quality and metadata files if available.

Additional data access and subsetting services¶

API Access¶

Data can be accessed directly from our HTTPS file system through the URLs collected above, or through our Application Programming Interface (API). Our API offers you the ability to order data using specific temporal and spatial filters, as well as subset, reformat, and reproject select data sets. The same subsetting, reformatting, and reprojection services available on select data sets through NASA Earthdata Search can also be applied using this API. These options can be requested in a single access command without the need to script against our data directory structure. See our programmatic access guide for more information on API options.

Add service request options for MODIS data¶

According to https://nsidc.org/support/faq/what-data-subsetting-reformatting-and-reprojection-services-are-available-for-MODIS-data, we can see that spatial subsetting and GeoTIFF reformatting are available for MOD10A1 so those options are requested below. The area subset must be described as a bounding box, which can be created based on the polygon bounds above. We will also add GeoTIFF reformatting to the MOD10A1 data dictionary and the temporal range will be set based on the range of MOD10A1 files in the dataframe above. These new parameters will be added to the API request below.

Create the data request API endpoint¶

Programmatic API requests are formatted as HTTPS URLs that contain key-value-pairs specifying the service operations that we specified above. We will first create a string of key-value-pairs from our data dictionary and we'll feed those into our API endpoint. This API endpoint can be executed via command line, a web browser, or in Python below.

Download options¶

The following functions will return the file URLs and the MOD10A1 API request. For demonstration purposes, these functions have been commented out, and instead the data utilized in the following steps will be accessed from a staged directory. *Note that if you are running this notebook in Binder, the memory may not be sufficient to download these files. Please use the Docker or local Conda options provided in the README if you are interested in downloading all files.*

Read in SnowEx data and buffer points around Snotel location¶

This SnowEx data set is provided in CSV. A Pandas DataFrame is used to easily read in data. For these next steps, just one day's worth of data will be selected from this file and the coincident ASO and MODIS data will be selected.

Convert to time values and extract a single day¶

The collection date needs to be extracted from the collection value and a new dataframe will be generated as a subset of the original based on a single day:

Convert to Geopandas dataframe to provide point geometry¶

According to the SnowEx documentation, the data are available in UTM Zone 12N so we'll set to this projection so that we can buffer in meters in the next step:

Buffer data around SNOTEL site¶

We can further subset the SnowEx snow depth data to get within a 500 m radius of the SNOTEL Mesa Lakes site.

First we'll create a new geodataframe with the SNOTEL site location, set to our SnowEx UTM coordinate reference system, and create a 500 meter buffer around this point. Then we'll subset the SnowEx points to the buffer and convert back to the WGS84 CRS:

Read in Airborne Snow Observatory data and clip to SNOTEL buffer¶

Snow depth data from the ASO L4 Lidar Snow Depth 3m UTM Grid data set were calculated from surface elevation measured by the Riegl LMS-Q1560 airborne laser scanner (ALS). The data are provided in GeoTIFF format, so we'll use the Rasterio library to read in the data.

Clip data to SNOTEL buffer¶

In order to reduce the data volume to the buffered region of interest, we can subset this GeoTIFF to the same SNOTEL buffer:

Read in MODIS Snow Cover data¶

We are interested in the Normalized Difference Snow Index (NDSI) snow cover value from the MOD10A1 data set, which is an index that is related to the presence of snow in a pixel. According to the MOD10A1 FAQ, snow cover is detected using the NDSI ratio of the difference in visible reflectance (VIS) and shortwave infrared reflectance (SWIR), where NDSI = ((band 4-band 6) / (band 4 + band 6)).

Note that you may need to change this filename output below if you download the data outside of the staged bucket, as the output names may vary per request.

Add ASO and MODIS data to GeoPandas dataframe¶

In order to add data from these ASO and MODIS gridded data sets, we need to define the geometry parameters for theses, as well as the SnowEx data. The SnowEx geometry is set using the longitude and latitude values of the geodataframe:

With ASO and MODIS data on regular grids, we can create area definitions for these using projection and extent metadata:

Interpolate ASO and MODIS values onto SnowEx points¶

To interpolate ASO snow depth and MODIS snow cover data to SnowEx snow depth points, we can use the pyresample library. The radius_of_influence parameter determines maximum radius to look for nearest neighbor interpolation.

Visualize data and export for further GIS analysis¶

The rasterio plot module allows you to directly plot GeoTIFFs objects. The SnowEx Thickness values are plotted against the clipped ASO snow depth raster.

We can do the same for MOD10A1: This was subsetted to the entire Grand Mesa region defined by the SnowEx data set coverage.

Additional data imagery services¶

NASA Worldview and the Global Browse Imagery Service¶

NASA’s EOSDIS Worldview mapping application provides the capability to interactively browse over 900 global, full-resolution satellite imagery layers and then download the underlying data. Many of the available imagery layers are updated within three hours of observation, essentially showing the entire Earth as it looks “right now."

According to the MOD10A1 landing page, snow cover imagery layers from this data set are available through NASA Worldview. This layer can be downloaded as various image files including GeoTIFF using the snapshot feature at the top right of the page. This link presents the MOD10A1 NDSI layer over our time and area of interest: https://go.nasa.gov/35CgYMd.

Additionally, the NASA Global Browse Imagery Service provides up to date, full resolution imagery for select NSIDC DAAC data sets as web services including WMTS, WMS, KML, and more. These layers can be accessed in GIS applications following guidance on the GIBS documentation pages.

Export dataframe to Shapefile¶

Finally, the dataframe can be exported as an Esri shapefile for further analysis in GIS: