# Import the following libraries
# For fetching from the Raster API
import requests
# For making maps
import folium
import folium.plugins
from folium import Map, TileLayer
# For talking to the STAC API
from pystac_client import Client
# For working with data
import pandas as pd
# For making time series
import matplotlib.pyplot as plt
# For formatting date/time data
import datetime
# Custom functions for working with GHGC data via the API
import ghgc_utils
CASA-GFED3 Land Carbon Flux
Run this notebook
You can launch this notebook in the US GHG Center JupyterHub by clicking the link below.
Approach
- Identify available dates and temporal frequency of observations for a given collection using the GHGC API
/stac
endpoint. The collection processed in this notebook is the Land-Atmosphere Carbon Flux data product. - Pass the STAC item into the raster API
/collections/{collection_id}/items/{item_id}/tilejson.json
endpoint. - Using
folium.plugins.DualMap
, visualize two tiles (side-by-side), allowing time point comparison. - After the visualization, perform zonal statistics for a given polygon.
About the Data
CASA-GFED3
This dataset presents a variety of carbon flux parameters derived from the Carnegie-Ames-Stanford-Approach – Global Fire Emissions Database version 3 (CASA-GFED3) model. The model’s input data includes air temperature, precipitation, incident solar radiation, a soil classification map, and a number of satellite derived products. All model calculations are driven by analyzed meteorological data from NASA’s Modern-Era Retrospective analysis for Research and Application, Version 2 (MERRA-2). The resulting product provides monthly, global data at 0.5 degree resolution from January 2003 through December 2017. It includes the following carbon flux variables expressed in units of kilograms of carbon per square meter per month (kg Carbon m²/mon) from the following sources: net primary production (NPP), net ecosystem exchange (NEE), heterotrophic respiration (Rh), wildfire emissions (FIRE), and fuel wood burning emissions (FUEL). This product and earlier versions of MERRA-driven CASA-GFED carbon fluxes have been used in a number of atmospheric CO₂ transport studies, and through the support of NASA’s Carbon Monitoring System (CMS), it helps characterize, quantify, understand and predict the evolution of global carbon sources and sinks.
Terminology
Navigating data via the GHGC API, you will encounter terminology that is different from browsing in a typical filesystem. We’ll define some terms here which are used throughout this notebook. - catalog
: All datasets available at the /stac
endpoint - collection
: A specific dataset, e.g. CarbonTracker-CH₄ Isotopic Methane Inverse Fluxes - item
: One granule in the dataset, e.g. one monthly file of methane inverse fluxes - asset
: A variable available within the granule, e.g. microbial, fossil, or pyrogenic methane fluxes - STAC API
: SpatioTemporal Asset Catalogs - Endpoint for fetching metadata about available datasets - Raster API
: Endpoint for fetching data itself, for imagery and statistics
Installing the Required Libraries
Required libraries are pre-installed on the GHG Center Hub. If you need to run this notebook elsewhere, please install them with this line in a code cell:
%pip install requests, folium, rasterstats, pystac_client, pandas, matplotlib
Query the STAC API
Now, you must fetch the dataset from the STAC API by defining its associated STAC API collection ID as a variable. The collection ID, also known as the collection name, for the CASA-GFED Land-Atmosphere Carbon Flux dataset is casagfed-carbonflux-monthgrid-v3.*
**You can find the collection name of any dataset on the GHGC data portal by navigating to the dataset landing page within the data catalog. The collection name is the last portion of the dataset landing page’s URL, and is also listed in the pop-up box after clicking “ACCESS DATA.”*
# Provide STAC and RASTER API endpoints
= "https://earth.gov/ghgcenter/api/stac"
STAC_API_URL = "https://earth.gov/ghgcenter/api/raster"
RASTER_API_URL
# Name of the collection for CASA GFED Land-Atmosphere Carbon Flux monthly emissions.
= "casagfed-carbonflux-monthgrid-v3" collection_name
= Client.open(STAC_API_URL)
catalog = catalog.get_collection(collection_name)
collection collection
- type "Collection"
- id "casagfed-carbonflux-monthgrid-v3"
- stac_version "1.0.0"
- description "This dataset presents a variety of carbon flux parameters derived from the Carnegie-Ames-Stanford-Approach – Global Fire Emissions Database version 3 (CASA-GFED3) model. All model calculations are driven by analyzed meteorological data from NASA's Modern-Era Retrospective analysis for Research and Application, Version 2 (MERRA-2). The resulting model output provides monthly, global data at 0.5 degree resolution from January 2003 through December 2017. It includes the following carbon flux variables expressed in units of kilograms of carbon per square meter per month (kg Carbon/m2/mon): net primary production (NPP), net ecosystem exchange (NEE), heterotrophic respiration (Rh), wildfire emissions (FIRE), and fuel wood burning emissions (FUEL). This product and earlier versions of MERRA-driven CASA-GFED carbon fluxes have been used in a number of atmospheric carbon dioxide (CO₂) transport studies, and through the support of NASA's Carbon Monitoring System (CMS), it helps characterize, quantify, understand and predict the evolution of global carbon sources and sinks. The source dataset can be found at https://doi.org/10.5067/03147VMJE8J9. As of April 2024, this dataset has been replaced by an updated version in the US GHG Center titled MiCASA Land Carbon Flux v1 (STAC ids: micasa-carbonflux-daygrid-v1 and micasa-carbonflux-monthgrid-v1)."
links[] 4 items
0
- rel "items"
- href "https://earth.gov/ghgcenter/api/stac/collections/casagfed-carbonflux-monthgrid-v3/items"
- type "application/geo+json"
1
- rel "parent"
- href "https://earth.gov/ghgcenter/api/stac/"
- type "application/json"
2
- rel "root"
- href "https://earth.gov/ghgcenter/api/stac"
- type "application/json"
- title "US GHG Center STAC API"
3
- rel "self"
- href "https://earth.gov/ghgcenter/api/stac/collections/casagfed-carbonflux-monthgrid-v3"
- type "application/json"
stac_extensions[] 2 items
- 0 "https://stac-extensions.github.io/render/v1.0.0/schema.json"
- 1 "https://stac-extensions.github.io/item-assets/v1.0.0/schema.json"
renders
rh
assets[] 1 items
- 0 "rh"
rescale[] 1 items
0[] 2 items
- 0 0
- 1 0.3
- colormap_name "purd"
nee
assets[] 1 items
- 0 "nee"
rescale[] 1 items
0[] 2 items
- 0 -0.1
- 1 0.1
- colormap_name "coolwarm"
npp
assets[] 1 items
- 0 "npp"
rescale[] 1 items
0[] 2 items
- 0 0
- 1 0.3
- colormap_name "purd"
fire
assets[] 1 items
- 0 "fire"
rescale[] 1 items
0[] 2 items
- 0 0
- 1 0.3
- colormap_name "purd"
fuel
assets[] 1 items
- 0 "fuel"
rescale[] 1 items
0[] 2 items
- 0 0
- 1 0.03
- colormap_name "bupu"
item_assets
rh
- type "image/tiff; application=geotiff; profile=cloud-optimized"
roles[] 2 items
- 0 "data"
- 1 "layer"
- title "Heterotrophic Respiration (Rh)"
- description "Model-estimated heterotrophic respiration (Rh), which is the flux of carbon from the soil to the atmosphere."
nee
- type "image/tiff; application=geotiff; profile=cloud-optimized"
roles[] 2 items
- 0 "data"
- 1 "layer"
- title "Net Ecosystem Exchange (NEE)"
- description "Model-estimated net ecosystem exchange (NEE), which is the net carbon flux to the atmosphere."
npp
- type "image/tiff; application=geotiff; profile=cloud-optimized"
roles[] 2 items
- 0 "data"
- 1 "layer"
- title "Net Primary Production (NPP)"
- description "Model-estimated net primary production (NPP), which is the amount of carbon available from plants."
fire
- type "image/tiff; application=geotiff; profile=cloud-optimized"
roles[] 2 items
- 0 "data"
- 1 "layer"
- title "Fire Emissions (FIRE)"
- description "Model-estimated flux of carbon to the atmosphere from wildfires."
fuel
- type "image/tiff; application=geotiff; profile=cloud-optimized"
roles[] 2 items
- 0 "data"
- 1 "layer"
- title "Wood Fuel Emissions (FUEL)"
- description "Model-estimated flux of carbon to the atmosphere from wood burned for fuel."
- dashboard:is_periodic True
- dashboard:time_density "month"
- title "CASA-GFED3 Land Carbon Flux v3"
extent
spatial
bbox[] 1 items
0[] 4 items
- 0 -180.0
- 1 -90.0
- 2 180.0
- 3 90.0
temporal
interval[] 1 items
0[] 2 items
- 0 "2003-01-01T00:00:00Z"
- 1 "2017-12-31T00:00:00Z"
- license "CC0-1.0"
summaries
datetime[] 2 items
- 0 "2003-01-01T00:00:00Z"
- 1 "2017-12-31T00:00:00Z"
= list(collection.get_items()) # Convert the iterator to a list
items print(f"Found {len(items)} items")
Found 180 items
= catalog.search(
search =collection_name,
collections=['2010-01-01T00:00:00Z','2010-12-31T00:00:00Z']
datetime
)# Take a look at the items we found
for item in search.item_collection():
print(item)
<Item id=casagfed-carbonflux-monthgrid-v3-201012>
<Item id=casagfed-carbonflux-monthgrid-v3-201011>
<Item id=casagfed-carbonflux-monthgrid-v3-201010>
<Item id=casagfed-carbonflux-monthgrid-v3-201009>
<Item id=casagfed-carbonflux-monthgrid-v3-201008>
<Item id=casagfed-carbonflux-monthgrid-v3-201007>
<Item id=casagfed-carbonflux-monthgrid-v3-201006>
<Item id=casagfed-carbonflux-monthgrid-v3-201005>
<Item id=casagfed-carbonflux-monthgrid-v3-201004>
<Item id=casagfed-carbonflux-monthgrid-v3-201003>
<Item id=casagfed-carbonflux-monthgrid-v3-201002>
<Item id=casagfed-carbonflux-monthgrid-v3-201001>
Examining the contents of our collection
under the temporal
variable, we see that the data is available from January 2003 to December 2017. By looking at the dashboard:time density
, we observe that the periodic frequency of these observations is monthly.
# Examine the first item in the collection
0] items[
- type "Feature"
- stac_version "1.0.0"
stac_extensions[] 2 items
- 0 "https://stac-extensions.github.io/raster/v1.1.0/schema.json"
- 1 "https://stac-extensions.github.io/projection/v1.1.0/schema.json"
- id "casagfed-carbonflux-monthgrid-v3-201712"
geometry
- type "Polygon"
coordinates[] 1 items
0[] 5 items
0[] 2 items
- 0 -180
- 1 -90
1[] 2 items
- 0 180
- 1 -90
2[] 2 items
- 0 180
- 1 90
3[] 2 items
- 0 -180
- 1 90
4[] 2 items
- 0 -180
- 1 -90
bbox[] 4 items
- 0 -180.0
- 1 -90.0
- 2 180.0
- 3 90.0
properties
- end_datetime "2017-12-31T00:00:00+00:00"
- start_datetime "2017-12-01T00:00:00+00:00"
- datetime None
links[] 4 items
0
- rel "collection"
- href "https://earth.gov/ghgcenter/api/stac/collections/casagfed-carbonflux-monthgrid-v3"
- type "application/json"
1
- rel "parent"
- href "https://earth.gov/ghgcenter/api/stac/collections/casagfed-carbonflux-monthgrid-v3"
- type "application/json"
2
- rel "root"
- href "https://earth.gov/ghgcenter/api/stac/"
- type "application/json"
3
- rel "self"
- href "https://earth.gov/ghgcenter/api/stac/collections/casagfed-carbonflux-monthgrid-v3/items/casagfed-carbonflux-monthgrid-v3-201712"
- type "application/geo+json"
assets
rh
- href "s3://ghgc-data-store/casagfed-carbonflux-monthgrid-v3/GEOSCarb_CASAGFED3v3_Rh_Flux_Monthly_x720_y360_201712.tif"
- type "image/tiff; application=geotiff"
- title "Heterotrophic Respiration (Rh)"
- description "Model-estimated heterotrophic respiration (Rh), which is the flux of carbon from the soil to the atmosphere."
proj:bbox[] 4 items
- 0 -180.0
- 1 -90.0
- 2 180.0
- 3 90.0
- proj:epsg 4326.0
proj:shape[] 2 items
- 0 360.0
- 1 720.0
raster:bands[] 1 items
0
- scale 1.0
- offset 0.0
- sampling "area"
- data_type "float32"
histogram
- max 0.6039900183677673
- min 0.0
- count 11.0
buckets[] 10 items
- 0 249101.0
- 1 7375.0
- 2 2429.0
- 3 252.0
- 4 32.0
- 5 5.0
- 6 2.0
- 7 2.0
- 8 0.0
- 9 2.0
statistics
- mean 0.006758838426321745
- stddev 0.022668374702334404
- maximum 0.6039900183677673
- minimum 0.0
- valid_percent 0.0003858024691358025
proj:geometry
- type "Polygon"
coordinates[] 1 items
0[] 5 items
0[] 2 items
- 0 -180.0
- 1 -90.0
1[] 2 items
- 0 180.0
- 1 -90.0
2[] 2 items
- 0 180.0
- 1 90.0
3[] 2 items
- 0 -180.0
- 1 90.0
4[] 2 items
- 0 -180.0
- 1 -90.0
proj:projjson
id
- code 4326.0
- authority "EPSG"
- name "WGS 84"
- type "GeographicCRS"
datum
- name "World Geodetic System 1984"
- type "GeodeticReferenceFrame"
ellipsoid
- name "WGS 84"
- semi_major_axis 6378137.0
- inverse_flattening 298.257223563
- $schema "https://proj.org/schemas/v0.4/projjson.schema.json"
coordinate_system
axis[] 2 items
0
- name "Geodetic latitude"
- unit "degree"
- direction "north"
- abbreviation "Lat"
1
- name "Geodetic longitude"
- unit "degree"
- direction "east"
- abbreviation "Lon"
- subtype "ellipsoidal"
proj:transform[] 9 items
- 0 0.5
- 1 0.0
- 2 -180.0
- 3 0.0
- 4 -0.5
- 5 90.0
- 6 0.0
- 7 0.0
- 8 1.0
roles[] 2 items
- 0 "data"
- 1 "layer"
nee
- href "s3://ghgc-data-store/casagfed-carbonflux-monthgrid-v3/GEOSCarb_CASAGFED3v3_NEE_Flux_Monthly_x720_y360_201712.tif"
- type "image/tiff; application=geotiff"
- title "Net Ecosystem Exchange (NEE)"
- description "Model-estimated net ecosystem exchange (NEE), which is the net carbon flux to the atmosphere."
proj:bbox[] 4 items
- 0 -180.0
- 1 -90.0
- 2 180.0
- 3 90.0
- proj:epsg 4326.0
proj:shape[] 2 items
- 0 360.0
- 1 720.0
raster:bands[] 1 items
0
- scale 1.0
- offset 0.0
- sampling "area"
- data_type "float32"
histogram
- max 0.48997998237609863
- min -0.11027999967336655
- count 11.0
buckets[] 10 items
- 0 663.0
- 1 234393.0
- 2 23809.0
- 3 282.0
- 4 37.0
- 5 10.0
- 6 4.0
- 7 0.0
- 8 0.0
- 9 2.0
statistics
- mean 0.0015448036137968302
- stddev 0.00977976992726326
- maximum 0.48997998237609863
- minimum -0.11027999967336655
- valid_percent 0.0003858024691358025
proj:geometry
- type "Polygon"
coordinates[] 1 items
0[] 5 items
0[] 2 items
- 0 -180.0
- 1 -90.0
1[] 2 items
- 0 180.0
- 1 -90.0
2[] 2 items
- 0 180.0
- 1 90.0
3[] 2 items
- 0 -180.0
- 1 90.0
4[] 2 items
- 0 -180.0
- 1 -90.0
proj:projjson
id
- code 4326.0
- authority "EPSG"
- name "WGS 84"
- type "GeographicCRS"
datum
- name "World Geodetic System 1984"
- type "GeodeticReferenceFrame"
ellipsoid
- name "WGS 84"
- semi_major_axis 6378137.0
- inverse_flattening 298.257223563
- $schema "https://proj.org/schemas/v0.4/projjson.schema.json"
coordinate_system
axis[] 2 items
0
- name "Geodetic latitude"
- unit "degree"
- direction "north"
- abbreviation "Lat"
1
- name "Geodetic longitude"
- unit "degree"
- direction "east"
- abbreviation "Lon"
- subtype "ellipsoidal"
proj:transform[] 9 items
- 0 0.5
- 1 0.0
- 2 -180.0
- 3 0.0
- 4 -0.5
- 5 90.0
- 6 0.0
- 7 0.0
- 8 1.0
roles[] 2 items
- 0 "data"
- 1 "layer"
npp
- href "s3://ghgc-data-store/casagfed-carbonflux-monthgrid-v3/GEOSCarb_CASAGFED3v3_NPP_Flux_Monthly_x720_y360_201712.tif"
- type "image/tiff; application=geotiff"
- title "npp"
- description "Model-estimated net primary production (NPP), which is the amount of carbon available from plants."
proj:bbox[] 4 items
- 0 -180.0
- 1 -90.0
- 2 180.0
- 3 90.0
- proj:epsg 4326.0
proj:shape[] 2 items
- 0 360.0
- 1 720.0
raster:bands[] 1 items
0
- scale 1.0
- offset 0.0
- sampling "area"
- data_type "float32"
histogram
- max 0.23635999858379364
- min 0.0
- count 11.0
buckets[] 10 items
- 0 244636.0
- 1 3051.0
- 2 1928.0
- 3 2634.0
- 4 4088.0
- 5 2211.0
- 6 428.0
- 7 156.0
- 8 59.0
- 9 9.0
statistics
- mean 0.005214035045355558
- stddev 0.021809572353959084
- maximum 0.23635999858379364
- minimum 0.0
- valid_percent 0.0003858024691358025
proj:geometry
- type "Polygon"
coordinates[] 1 items
0[] 5 items
0[] 2 items
- 0 -180.0
- 1 -90.0
1[] 2 items
- 0 180.0
- 1 -90.0
2[] 2 items
- 0 180.0
- 1 90.0
3[] 2 items
- 0 -180.0
- 1 90.0
4[] 2 items
- 0 -180.0
- 1 -90.0
proj:projjson
id
- code 4326.0
- authority "EPSG"
- name "WGS 84"
- type "GeographicCRS"
datum
- name "World Geodetic System 1984"
- type "GeodeticReferenceFrame"
ellipsoid
- name "WGS 84"
- semi_major_axis 6378137.0
- inverse_flattening 298.257223563
- $schema "https://proj.org/schemas/v0.4/projjson.schema.json"
coordinate_system
axis[] 2 items
0
- name "Geodetic latitude"
- unit "degree"
- direction "north"
- abbreviation "Lat"
1
- name "Geodetic longitude"
- unit "degree"
- direction "east"
- abbreviation "Lon"
- subtype "ellipsoidal"
proj:transform[] 9 items
- 0 0.5
- 1 0.0
- 2 -180.0
- 3 0.0
- 4 -0.5
- 5 90.0
- 6 0.0
- 7 0.0
- 8 1.0
roles[] 2 items
- 0 "data"
- 1 "layer"
fire
- href "s3://ghgc-data-store/casagfed-carbonflux-monthgrid-v3/GEOSCarb_CASAGFED3v3_FIRE_Flux_Monthly_x720_y360_201712.tif"
- type "image/tiff; application=geotiff"
- title "Fire Emissions (FIRE)"
- description "Fire emissions"
proj:bbox[] 4 items
- 0 -180.0
- 1 -90.0
- 2 180.0
- 3 90.0
- proj:epsg 4326.0
proj:shape[] 2 items
- 0 360.0
- 1 720.0
raster:bands[] 1 items
0
- scale 1.0
- offset 0.0
- sampling "area"
- data_type "float32"
histogram
- max 0.7556899785995483
- min 0.0
- count 11.0
buckets[] 10 items
- 0 258952.0
- 1 161.0
- 2 53.0
- 3 22.0
- 4 11.0
- 5 0.0
- 6 0.0
- 7 0.0
- 8 0.0
- 9 1.0
statistics
- mean 0.00025634843041189015
- stddev 0.005492232274264097
- maximum 0.7556899785995483
- minimum 0.0
- valid_percent 0.0003858024691358025
proj:geometry
- type "Polygon"
coordinates[] 1 items
0[] 5 items
0[] 2 items
- 0 -180.0
- 1 -90.0
1[] 2 items
- 0 180.0
- 1 -90.0
2[] 2 items
- 0 180.0
- 1 90.0
3[] 2 items
- 0 -180.0
- 1 90.0
4[] 2 items
- 0 -180.0
- 1 -90.0
proj:projjson
id
- code 4326.0
- authority "EPSG"
- name "WGS 84"
- type "GeographicCRS"
datum
- name "World Geodetic System 1984"
- type "GeodeticReferenceFrame"
ellipsoid
- name "WGS 84"
- semi_major_axis 6378137.0
- inverse_flattening 298.257223563
- $schema "https://proj.org/schemas/v0.4/projjson.schema.json"
coordinate_system
axis[] 2 items
0
- name "Geodetic latitude"
- unit "degree"
- direction "north"
- abbreviation "Lat"
1
- name "Geodetic longitude"
- unit "degree"
- direction "east"
- abbreviation "Lon"
- subtype "ellipsoidal"
proj:transform[] 9 items
- 0 0.5
- 1 0.0
- 2 -180.0
- 3 0.0
- 4 -0.5
- 5 90.0
- 6 0.0
- 7 0.0
- 8 1.0
roles[] 2 items
- 0 "data"
- 1 "layer"
fuel
- href "s3://ghgc-data-store/casagfed-carbonflux-monthgrid-v3/GEOSCarb_CASAGFED3v3_FUEL_Flux_Monthly_x720_y360_201712.tif"
- type "image/tiff; application=geotiff"
- title "Wood Fuel Emissions (FUEL)"
- description "Fuel emissions"
proj:bbox[] 4 items
- 0 -180.0
- 1 -90.0
- 2 180.0
- 3 90.0
- proj:epsg 4326.0
proj:shape[] 2 items
- 0 360.0
- 1 720.0
raster:bands[] 1 items
0
- scale 1.0
- offset 0.0
- sampling "area"
- data_type "float32"
histogram
- max 0.020759999752044678
- min 0.0
- count 11.0
buckets[] 10 items
- 0 257568.0
- 1 1150.0
- 2 284.0
- 3 115.0
- 4 47.0
- 5 21.0
- 6 5.0
- 7 6.0
- 8 3.0
- 9 1.0
statistics
- mean 5.057307134848088e-05
- stddev 0.0003876804548781365
- maximum 0.020759999752044678
- minimum 0.0
- valid_percent 0.0003858024691358025
proj:geometry
- type "Polygon"
coordinates[] 1 items
0[] 5 items
0[] 2 items
- 0 -180.0
- 1 -90.0
1[] 2 items
- 0 180.0
- 1 -90.0
2[] 2 items
- 0 180.0
- 1 90.0
3[] 2 items
- 0 -180.0
- 1 90.0
4[] 2 items
- 0 -180.0
- 1 -90.0
proj:projjson
id
- code 4326.0
- authority "EPSG"
- name "WGS 84"
- type "GeographicCRS"
datum
- name "World Geodetic System 1984"
- type "GeodeticReferenceFrame"
ellipsoid
- name "WGS 84"
- semi_major_axis 6378137.0
- inverse_flattening 298.257223563
- $schema "https://proj.org/schemas/v0.4/projjson.schema.json"
coordinate_system
axis[] 2 items
0
- name "Geodetic latitude"
- unit "degree"
- direction "north"
- abbreviation "Lat"
1
- name "Geodetic longitude"
- unit "degree"
- direction "east"
- abbreviation "Lon"
- subtype "ellipsoidal"
proj:transform[] 9 items
- 0 0.5
- 1 0.0
- 2 -180.0
- 3 0.0
- 4 -0.5
- 5 90.0
- 6 0.0
- 7 0.0
- 8 1.0
roles[] 2 items
- 0 "data"
- 1 "layer"
- collection "casagfed-carbonflux-monthgrid-v3"
# Restructure the items into a dictionary where keys are derived from the datetime items; we can then query more easily by date/time, e.g. "2003-07"
= {item.properties["start_datetime"][:7]: item for item in collection.get_items()} items_dict
# Before we go further, let's pick which asset to focus on for the remainder of the notebook.
# 'npp' = net primary production
= "npp" asset_name
Creating Maps using Folium
We will explore differences in the land atmosphere Carbon flux Net Primary Productivity between two different dates/times. We’ll then visualize the outputs on a map using folium
.
Fetch Imagery Using the Raster API
Here we get information from the Raster API which we will add to our map in the next section.
# Specify which two date/times you would like to visualize, using the format of items_dict.keys()
= ["2010-01","2003-01"] dates
Below, we use some statistics of the raster data to set upper and lower limits for our color bar. These are saved as the rescale_values
, and will be passed to the Raster API in the following step(s).
# Extract collection name and item ID for the first date
= items_dict[dates[0]]
first_date = first_date.collection_id
collection_id = first_date.id
item_id # Select relevant asset (NPP)
object = first_date.assets[asset_name]
= object.extra_fields.get("raster:bands", [{}])
raster_bands # Print the raster bands' information
raster_bands
[{'scale': 1.0,
'offset': 0.0,
'sampling': 'area',
'data_type': 'float32',
'histogram': {'max': 0.23026999831199646,
'min': 0.0,
'count': 11.0,
'buckets': [244259.0,
3221.0,
2076.0,
2056.0,
3230.0,
3103.0,
996.0,
195.0,
53.0,
11.0]},
'statistics': {'mean': 0.00534185953438282,
'stddev': 0.022541318088769913,
'maximum': 0.23026999831199646,
'minimum': 0.0,
'valid_percent': 0.0003858024691358025}}]
# Use mean, scaled stddev, and minimum to generate an appropriate color bar range.
= {
rescale_values #"max": raster_bands[0]['statistics']['mean'] + 4*raster_bands[0]['statistics']['stddev'],
"max": raster_bands[0]['statistics']['maximum'],
"min": raster_bands[0]['statistics']['minimum'],
}
print(rescale_values)
{'max': 0.23026999831199646, 'min': 0.0}
Now, you will pass the item id
, collection name
, asset name
, and the rescale values
to the Raster API endpoint, along with a colormap. This step is done twice, one for each date/time you will visualize, and tells the Raster API which collection, item, and asset you want to view, specifying the colormap and colorbar ranges to use for visualization. The API returns a JSON with information about the requested image. Each image will be referred to as a tile.
# Choose a color for displaying the data
# For more information on Colormaps in Matplotlib, please visit https://matplotlib.org/stable/users/explain/colors/colormaps.html
= "viridis" color_map
# Make a GET request to retrieve information for your first date/time
= requests.get(
date_1_tile f"{RASTER_API_URL}/collections/{collection_id}/items/{item_id}/tilejson.json?"
f"&assets={asset_name}"
f"&color_formula=gamma+r+1.05&colormap_name={color_map}"
f"&rescale={rescale_values['min']},{rescale_values['max']}"
).json()
# Print the properties of the retrieved granule to the console
date_1_tile
{'tilejson': '2.2.0',
'version': '1.0.0',
'scheme': 'xyz',
'tiles': ['https://earth.gov/ghgcenter/api/raster/collections/casagfed-carbonflux-monthgrid-v3/items/casagfed-carbonflux-monthgrid-v3-201001/tiles/WebMercatorQuad/{z}/{x}/{y}@1x?assets=npp&color_formula=gamma+r+1.05&colormap_name=viridis&rescale=0.0%2C0.23026999831199646'],
'minzoom': 0,
'maxzoom': 24,
'bounds': [-180.0, -90.0, 180.0, 90.0],
'center': [0.0, 0.0, 0]}
# Repeat the above for your second date/time
# Note that we do not calculate new rescale_values for this tile, because we dates tiles 1 and 2 to have the same colorbar range for best visual comparison.
= items_dict[dates[1]]
second_date # Extract collection name and item ID
= second_date.collection_id
collection_id = second_date.id
item_id
= requests.get(
date_2_tile f"{RASTER_API_URL}/collections/{collection_id}/items/{item_id}/tilejson.json?"
f"&assets={asset_name}"
f"&color_formula=gamma+r+1.05&colormap_name={color_map}"
f"&rescale={rescale_values['min']},{rescale_values['max']}"
).json()
date_2_tile
{'tilejson': '2.2.0',
'version': '1.0.0',
'scheme': 'xyz',
'tiles': ['https://earth.gov/ghgcenter/api/raster/collections/casagfed-carbonflux-monthgrid-v3/items/casagfed-carbonflux-monthgrid-v3-200301/tiles/WebMercatorQuad/{z}/{x}/{y}@1x?assets=npp&color_formula=gamma+r+1.05&colormap_name=viridis&rescale=0.0%2C0.23026999831199646'],
'minzoom': 0,
'maxzoom': 24,
'bounds': [-180.0, -90.0, 180.0, 90.0],
'center': [0.0, 0.0, 0]}
Generate Map
First, we’ll define the Area of Interest (AOI) as a GEOJSON. This will be visualized as a filled polygon on the map.
# The AOI is currently set to the Amazon rainforest in South America
= {
aoi "type": "Feature",
"properties": {},
"geometry": {
"coordinates": [
[# [longitude, latitude]
-74.0, -3.0], # Southwest Bounding Coordinate
[-74.0, 5.0], # Southeast Bounding Coordinate
[-60.0, 5.0], # Northeast Bounding Coordinate
[-60.0, -3.0], # Northwest Bounding Coordinate
[-74.0, -3.0] # Closing the polygon at the Southwest Bounding Coordinate
[
]
],"type": "Polygon",
}, }
# Initialize the map, specifying the center of the map and the starting zoom level.
# 'folium.plugins' allows mapping side-by-side via 'DualMap'
# Map is centered on the position specified by "location=(lat,lon)"
= folium.plugins.DualMap(location=(0, -66), zoom_start=5)
map_
# Define the first map layer using the tile fetched for the first date
# The TileLayer library helps in manipulating and displaying raster layers on a map
= TileLayer(
map_layer_1 =date_1_tile["tiles"][0],
tiles="GHG",
attr=0.8,
opacity=f"{dates[0]} NPP",
name= True,
overlay= True
legendEnabled
)# Add the first layer to the Dual Map
# This will appear on the left side, specified by 'm1'
map_layer_1.add_to(map_.m1)
# Define the second map layer using the tile fetched for the second date
= TileLayer(
map_layer_2 =date_2_tile["tiles"][0],
tiles="GHG",
attr=0.8,
opacity=f"{dates[1]} NPP",
name= True,
overlay= True
legendEnabled
)# Add the second layer to the Dual Map
# This will appear on the right side, specified by 'm2'
map_layer_2.add_to(map_.m2)
# Add data markers to both sides of map
0, -66), tooltip="Amazon Rainforest").add_to(map_)
folium.Marker((# Add AOI to both sides of map
="Amazon Rainforest, South America",
folium.GeoJson(aoi, name=lambda feature: {
style_function"fillColor": "none",
}).add_to(map_)# Add controls to turn different elements on/off, for both sides of map
=False).add_to(map_)
folium.LayerControl(collapsed
# Add a colorbar
# For this, use one of our custom 'ghgc_utils' functions to create an HTML colorbar representation.
= ghgc_utils.generate_html_colorbar(color_map, rescale_values, label='NPP (kg Carbon/m2/month)')
legend_html # Now add colorbar to the map
map_.get_root().html.add_child(folium.Element(legend_html))
# Visualizing the map
map_
Calculating Zonal Statistics
To perform zonal statistics, first we need to create a polygon. In this case we are focusing on an Area of Interest (AOI) in the Amazon Rainforest, South America.
# Give the AOI a name to be used in your time series plot later on.
= 'Amazon Rainforest'
aoi_name # The AOI is defined as a GEOJSON
= {
aoi "type": "Feature",
"properties": {},
"geometry": {
"coordinates": [
[# [longitude, latitude]
-74.0, -3.0], # Southwest Bounding Coordinate
[-74.0, 5.0], # Southeast Bounding Coordinate
[-60.0, 5.0], # Northeast Bounding Coordinate
[-60.0, -3.0], # Northwest Bounding Coordinate
[-74.0, -3.0] # Closing the polygon at the Southwest Bounding Coordinate
[
]
],"type": "Polygon",
}, }
# Use one of the custom 'ghgc_utils' functions to generate statistics over your AOI using the Raster API
# This step may take a minute.
= ghgc_utils.generate_stats(items,aoi,url=RASTER_API_URL,asset=asset_name)
df # Print the first 5 lines of our statistics
5) df.head(
Generating stats...
Visualizing the Data as a Time Series
We can now explore the NPP time series for the Amazon Rainforest, South America area. We can plot the data set using the code below:
# Plot data
= plt.figure(figsize=(10,5)) # Set the size of the figure
fig
# Change "which_stat" if you would rather look at something like mean, median, or standard deviation.
= 'max'
which_stat
plt.plot("date"][:],
df[
df[which_stat],="purple",
color="-",
linestyle=1.5,
linewidth
)
# Add x labels at desired positions (for example, every 6 months)
plt.xticks("date"][::6],
df[=45, # Rotate labels to avoid overlap
rotation="right"
ha
)
# Labels and title
"Month")
plt.xlabel("kg Carbon/m2/month")
plt.ylabel(f"{which_stat.capitalize()} Monthly {asset_name.upper()} for {aoi_name} (2003-2017)")
plt.title(
# Add data citation
plt.text("date"][:].min(), # X-coordinate of the text (first datetime value)
df[min(), # Y-coordinate of the text (minimum CO2 value)
df[which_stat].
# Text to be displayed
f"Source: {collection.title}",
=9, # Font size
fontsize="left", # Horizontal alignment
horizontalalignment="bottom", # Vertical alignment
verticalalignment="blue", # Text color
color
)
plt.show()
# Fetch the third granule in the collection and set the color scheme and rescale values.
= 2
n = requests.get(
october_tile f"{RASTER_API_URL}/collections/{collection_name}/items/{items[n].id}/tilejson.json?"
f"&assets={asset_name}"
f"&color_formula=gamma+r+1.05&colormap_name={color_map}"
f"&rescale={rescale_values['min']},{rescale_values['max']}",
).json() october_tile
# Map the NPP level for the Congo area for the chosen tile
= Map(
aoi_map_bbox ="OpenStreetMap",
tiles=[
location1, # latitude
17, # longitude
],=5,
zoom_start
)
= TileLayer(
map_layer =october_tile["tiles"][0],
tiles="GHG", opacity = 0.7, name=f"{items[n].properties['start_datetime'][0:7]} {asset_name}", overlay= True, legendEnabled = True
attr
)
map_layer.add_to(aoi_map_bbox)
# Display data marker (title) on the map
1,17), tooltip="Dallas").add_to(aoi_map_bbox)
folium.Marker((=False).add_to(aoi_map_bbox)
folium.LayerControl(collapsed
# Add a colorbar
# Add a colorbar
# For this, use one of the custom 'ghgc_utils' functions to create an HTML colorbar representation.
= ghgc_utils.generate_html_colorbar(color_map, rescale_values, label='NPP (kg Carbon/m2/month)')
legend_html
# Add colorbar to the map
aoi_map_bbox.get_root().html.add_child(folium.Element(legend_html))
aoi_map_bbox
Summary
In this notebook we have successfully completed the following steps for the STAC collection for CASA GFED Land-Atmosphere Carbon Flux data: 1. Install and import the necessary libraries 2. Fetch the collection from STAC collections using the appropriate endpoints 3. Count the number of existing granules within the collection 4. Map and compare the Net Primary Production (NPP) levels over the Amazon Rainforest, South America area for two distinctive years 5. Create a table that displays the minimum, maximum, and sum of the NPP values for a specified region 6. Generate a time-series graph of the NPP values for a specified region
If you have any questions regarding this user notebook, please contact us using the feedback form.