This notebook is from EMIT-Data-Resources

Source: Generating Methane Spectral Fingerprint

Imported on: 2025-01-22

Generating Methane Spectral Fingerprint

Summary

In this notebook, we’ll examine a EMIT L1B At-Sensor Calibrated Radiance (EMITL1BRAD) scene to visualize the spectra for a sample point inside a methane plume, and two sample points outside of methane plumes. To highlight the spectral differences we will look at an in-plume divided by out-of-plume ratio to get an example of the spectral signature or fingerprint of methane and visually compare this to the target methane signature.

Background

Methane is a major contributor to atmospheric radiative forcing because its more efficient at trapping radiation than other greenhouse gases, like carbon dioxide. Because the lifetime of methane in the atmosphere is only about 10 years, reducing methane emissions offers an effective way to curb anthropogenic contributions to atmospheric radiative forcing.

The EMIT instrument is an imaging spectrometer that measures light in visible and infrared wavelengths. These measurements display unique spectral signatures that correspond to chemical composition. Although the primary mission focuses on mapping the mineral composition of Earth’s surface, the EMIT instrument has also been used to successfully map methane point source emissions within its target mask using these unique spectral signatures. More details about EMIT and its associated products can be found in the README.md and on the EMIT website.

The EMITL1BRAD product provides at-sensor calibrated radiance values along with observation data in a spatially raw, non-orthocorrected format. This product is used to estimate the EMIT L2B Methane Enhancement Data (EMITL2BCH4ENH) product. The EMITL2BCH4ENH product is created using an adaptive matched filter to search each pixel’s radiance spectrum for deviations that are characteristic of methane’s absorption spectrum, and is only delivered to the LP DAAC for scenes where plume complexes have been identified. To reduce the risk of false positives, all EMITL2BCH4ENH data undergo a manual review process before being designated as a plume complex. For more information on the identification and manual review process, see Section 4.2.2 of the EMIT GHG Algorithm Theoretical Basis Document (ATBD).

References

Andrew K. Thorpe et al., Attribution of individual methane and carbon dioxide emission sources using EMIT observations from space. Sci. Adv.9, eadh2391 (2023). DOI:10.1126/sciadv.adh2391

Requirements - Set up Python Environment - See setup_instructions.md in the /setup/ folder - NASA Earthdata Login Account. Sign Up

Data Used - EMIT L1B At-Sensor Calibrated Radiance and Geolocation Data (EMITL1BRAD) - EMIT L2B Methane Enhancement Data (EMITL2BCH4ENH)

Learning Objectives - Open and orthorectify an EMITL1BRAD scene - Open and visualize the scene and point data - Extract point data using coordinates in a .csv file - Visualize spectral fingerprint of methane and compare it to the modeled methane signature

Tutorial Outline

1. Setup

Import the necessary Python libraries.

import os
import sys
import numpy as np
import pandas as pd
from osgeo import gdal
import earthaccess
import rasterio as rio
import rioxarray as rxr
import holoviews as hv
import hvplot
import hvplot.xarray
import hvplot.pandas
from skimage import exposure

sys.path.append('../modules/')
from emit_tools import emit_xarray

Download the EMIT L1B Radiance and L2B Methane Enhancement products for the scene we’re going to look at.

urls = ['https://data.lpdaac.earthdatacloud.nasa.gov/lp-prod-protected/EMITL1BRAD.001/EMIT_L1B_RAD_001_20220815T042838_2222703_003/EMIT_L1B_RAD_001_20220815T042838_2222703_003.nc',
        'https://data.lpdaac.earthdatacloud.nasa.gov/lp-prod-protected/EMITL2BCH4ENH.001/EMIT_L2B_CH4ENH_001_20220815T042838_2222703_003/EMIT_L2B_CH4ENH_001_20220815T042838_2222703_003.tif']

# Authenticate and create an https session
earthaccess.login(persist=True)
fs = earthaccess.get_requests_https_session()

for url in urls:
# Retrieve granule asset ID from URL (to maintain existing naming convention)
    granule_asset_id = url.split('/')[-1]
    # Define Local Filepath
    fp = f'../../data/{granule_asset_id}'
    # Download the Granule Asset if it doesn't exist
    if not os.path.isfile(fp):
        with fs.get(url,stream=True) as src:
            with open(fp,'wb') as dst:
                for chunk in src.iter_content(chunk_size=64*1024*1024):
                    dst.write(chunk)

Set the filepath for the radiance and methane enhancement files.

rad_fp = '../../data/EMIT_L1B_RAD_001_20220815T042838_2222703_003.nc'
enh_fp = '../../data/EMIT_L2B_CH4ENH_001_20220815T042838_2222703_003.tif'

2. Opening EMIT Data

The EMIT L1B At-Sensor Radiance data is distributed in a non-orthorectified spatially raw netCDF4 (.nc) format consisting of the data and its associated metadata. Inside the L1B file, there are 3 groups.

The root group that can be considered the main dataset contains the radiance data described by the downtrack, crosstrack, and bands dimensions.
The sensor_band_parameters group containing the wavelength center and the full-width half maximum (FWHM) of each band.
The location group contains latitude and longitude values at the center of each pixel described by the crosstrack and downtrack dimensions, as well as a geometry lookup table (GLT) described by the ortho_x and ortho_y dimensions. The GLT is an orthorectified image (EPSG:4326) consisting of 2 layers containing downtrack and crosstrack indices. These index positions allow us to quickly project the raw data onto this geographic grid.

This data can be opened using the netCDF4 and xarray libraries, or utilizing functions within the emit_tools.py module to organize them into a flattened (no groups) xarray.Dataset object. For this notebook, we will use functions from emit_tools.py to simplify working with the data. For more about the structure and use of netCDF4 and xarray please see the Exploring EMIT L2A Reflectance Jupyter Notebook.

Open the radiance file using the emit_xarray function from the emit_tools.py module and orthorecitify it.

rad = emit_xarray(rad_fp, ortho=True)
rad

The EMIT L2B Methane Enhancement Data represent an enhancement above background methane concentration for a 1 meter layer in parts-per-million (ppm) meter (m). These units are used rather than ppm because we are unable to measure the vertical extent of plumes. This data is distributed as a single band cloud-optimized geotiff (COG) and it has been orthocorrected. We can open this using the rioxarray library to place the data in an xarray.DataArray object.

Open the methane enhancement geotiff file using rioxarray and squeeze the band dimension to remove the extra dimension so our array will only have 2 dimensions.

enh = rxr.open_rasterio(enh_fp).squeeze('band',drop=True) 
enh

3. Extracting Point Data

Open the .csv file included in the /data/ directory as a pandas.DataFrame. This file contains the latitude and longitude coordinates for three points of interest. Two that are outside of a methane plume, and one that is inside.

# Define our Points for In-plume and out-of-plume
points = pd.read_csv('../../data/methane_tutorial/methane_inout_points.csv')
points

Set the index as the ID column.

points = points.set_index(['ID'])
points

Now we can use the sel function from xarray to extract the radiance data at each point in our dataframe.

# Extract target spectra from dataset
point_ds = rad.sel(latitude=points.to_xarray().latitude, longitude=points.to_xarray().longitude, method='nearest')

After extracting, we can convert the data to a pandas.DataFrame and join it with our ‘in-plume’ column from the original dataframe.

point_df = point_ds.to_dataframe().join(points['in-plume'],on=['ID'])
point_df

At this point, we can save the data to a .csv file for future use.

# point_df.to_csv('../../data/methane_tutorial/point_df.csv')

Next, visualize these points on an rgb image generated from the radiance file, and the methane enhancement data, just to get a better idea of where these points are located.

To do this, first create an RGB data array from the radiance file using the sel function to select the bands nearest to the desired wavelengths.

# Create an RGB from Radiance
rgb = rad.sel(wavelengths=[650,560,470], method='nearest')

Next, use a function to rescale the brightness, so this image is easier to see.

rgb.radiance.data[rgb.radiance.data == -9999] = 0
rgb.radiance.data = exposure.rescale_intensity(rgb.radiance.data, in_range='image', out_range=(0,1))

Now that we have the necessary pieces to visualize we can make some spatial visualizations using hvplot. For the enhancement data, set the fill value of -9999 to np.nan to make it transparent.

# Create RGB Plot
rgb_map = rgb.hvplot.rgb(x='longitude',y='latitude',bands='wavelengths',title='RGB Radiance', geo=True, crs='EPSG:4326')

# Create Methane Enhancement Plot
enh.data[enh.data == -9999] = np.nan
methane_map = enh.hvplot.image(x='x',y='y',cmap='viridis', geo=True, crs='EPSG:4326', clim=(0,enh.data.max()), title='Methane Enhancement', clabel='ppm m', xlabel='Longitude', ylabel='Latitude')

point_map = point_df.hvplot.points(x='longitude',y='latitude', color='in-plume', cmap='HighContrast', geo=True, crs='EPSG:4326', hover=False, colorbar=False)

We can combine these plots with an * operator to overlay them on the same plot. This does require that all are in the same CRS to display properly.

(rgb_map*point_map) + (methane_map*point_map)

Now lets look at the spectra.

point_df.hvplot.line(x='wavelengths',y='radiance', by=['ID'], color=hv.Cycle('Dark2'), frame_height=400, frame_width=600, title = 'Radiance Spectra, ID 0 is in-plume' , xlabel='Wavelength (nm)', ylabel='Radiance (W/m^2/sr/nm)')

4. Methane Spectral Signature

Let’s open a file containing the modeled methane signature and visualize it. This is the spectral fingerprint that we are looking for in EMIT data to identify methane plumes. Open our absorption coefficient file using pandas, add some column names and set an index.

# Open file 
ch4_ac = pd.read_csv('../../data/methane_tutorial/emit20220815t042838_ch4_target', sep='\s+', header=None)
# Add Column Names
ch4_ac.columns = ['index','wavelength','value']
# Set Index
ch4_ac.set_index('index', inplace=True)
ch4_ac

Create a figure for the absorption coefficient.

ac_plot = ch4_ac.hvplot(x='wavelength',y='value', frame_height=400, frame_width=400, line_color='black', line_width=2, xlim=(2150,2450), ylim=(-1.5,0), xlabel='Wavelength (nm)', title='Methane Absorption Coefficient', ylabel='')

ac_plot

We can visualize a similar curve by creating a band ratio of in-plume to out-of-plume spectra. dividing the radiance for a pixel with methane by radiance for a pixel outside the plume to generate a diagnostic spectral fingerprint. This spectral fingerprint is confirmation that EMIT is observing methane enhancements with characteristic methane absorption features, which agree well with the modeled methane signature.

This agreement is strong in the example we selected due to the similarity of the spectra without the contribution of methane.

To create the band ratio, first separate our data into in-plume and out-of-plume dataframes.

in_plume = point_df.loc[point_df['in-plume'] == 1].copy()
out_plume = point_df.loc[point_df['in-plume'] == 0].copy()

Next, add a column for the in/out band ratio using the in-plume divided by out-of-plume radiance.

out_plume['band_ratio'] = (in_plume.loc[0,'radiance']/out_plume['radiance'])

Create an hvplot object for our band ratio.

in_out_plot = out_plume.hvplot(x='wavelengths',y='band_ratio', by=['ID'], color=hv.Cycle('Dark2'), frame_height=400, frame_width=400, xlim=(2150,2450), ylim=(0.85,1.05), ylabel='In Plume/Out of Plume Ratio', xlabel='Wavelength (nm)', title='In Plume/Out of Plume Ratio')

Overlay our absorption coefficient to show similarity between the two within the 2150 and 2450 nanometers (nm) range where methane spectral features are present. We can do this by setting our figure xlim.

from bokeh.models import GlyphRenderer, LinearAxis, LinearScale, Range1d

def overlay_hook(plot, element):
    # Adds right y-axis
    p = plot.handles["plot"]
    p.extra_y_scales = {"right": LinearScale()}
    p.extra_y_ranges = {"right": Range1d(-1.5,0)}
    p.add_layout(LinearAxis(y_range_name="right"), "right")

   # find the last line and set it to right
    lines = [p for p in p.renderers if isinstance(p, GlyphRenderer)]
    lines[-1].y_range_name = "right"

# Create Figure
(in_out_plot.opts(ylim=(0.85,0.95)) * ac_plot.opts(color="k")).opts(hooks=[overlay_hook]).opts(title='In Plume/Out of Plume and Absorption Coefficient')

We can contrast this with our out-of-plume/out-of-plume band ratios to show how similar the in-plume and out-of-plume ratio is to this spectral signature. Create a column in our dataframe for the out-of-plume/out-of-plume band ratio and visualize it overlayed with our methane absorption feature using hvplot.

out_plume['out-out'] = (out_plume['radiance'] / out_plume.loc[2,'radiance'])

out_out_plot = out_plume.hvplot(x='wavelengths',y='out-out', by=['ID'], color=hv.Cycle('Dark2'), frame_height=400, frame_width=400, xlim=(2150,2450), ylim=(0.85,1.05), ylabel='Out/Out 2 Ratio', xlabel='Wavelength (nm)', title='Out of Plume/Out of Plume 2 Ratio')

from bokeh.models import GlyphRenderer, LinearAxis, LinearScale, Range1d

def overlay_hook(plot, element):
    # Adds right y-axis
    p = plot.handles["plot"]
    p.extra_y_scales = {"right": LinearScale()}
    p.extra_y_ranges = {"right": Range1d(-1.5,0)}
    p.add_layout(LinearAxis(y_range_name="right"), "right")

   # find the last line and set it to right
    lines = [p for p in p.renderers if isinstance(p, GlyphRenderer)]
    lines[-1].y_range_name = "right"

# Create Figure
(out_out_plot * ac_plot.opts(color="k")).opts(hooks=[overlay_hook]).opts(title='Out of Plume/Out of Plume and Absorption Coefficient')

This highlights the similarity of the in-plume/out-of-plume band ratio to the modeled methane signature, as opposed to a band ratio where methane is not present.

Contact Info:

Email: LPDAAC@usgs.gov
Voice: +1-866-573-3222
Organization: Land Processes Distributed Active Archive Center (LP DAAC)¹
Website: https://lpdaac.usgs.gov/
Date last modified: 03-13-2024

¹Work performed under USGS contract G15PD00467 for NASA contract NNG14HH33I.