EMIT¶
Install package with EMIT dependecies¶
EMIT requires the netcdf4 package to read the products. We will also use the pysolar package to convert the EMIT L1B radiances to reflectances.
pip install georeader-spaceml netcdf4 pysolar
Download an EMIT image¶
from georeader.readers import emit
from georeader.readers import download_utils
import os
dir_emit_files = "emit_database/raw"
os.makedirs(dir_emit_files, exist_ok=True)
# link = 'https://data.lpdaac.earthdatacloud.nasa.gov/lp-prod-protected/EMITL1BRAD.001/EMIT_L1B_RAD_001_20220828T051941_2224004_006/EMIT_L1B_RAD_001_20220828T051941_2224004_006.nc'
link = 'https://data.lpdaac.earthdatacloud.nasa.gov/lp-prod-protected/EMITL1BRAD.001/EMIT_L1B_RAD_001_20220827T060753_2223904_013/EMIT_L1B_RAD_001_20220827T060753_2223904_013.nc'
file_save = os.path.join(dir_emit_files, os.path.basename(link))
emit.download_product(link, file_save)
Create and inspec emit object¶
EMIT objects let you open an EMIT nc file without loading the content of the file in memory. The object in the cell below has 285 spectral bands. For the API description of the EMIT class see the emit module. Since the object follows the API of georeader, you can read from the rasters using the functions from the read module.
# file_save = "emit_database/EMIT_L1B_RAD_001_20220827T060753_2223904_013.nc"
ei = emit.EMITImage(file_save)
ei
ei.nc_ds
ei.wavelengths
ei.time_coverage_start, ei.time_coverage_end
import geopandas as gpd
gpd.GeoDataFrame(geometry=[ei.footprint()],
crs=ei.crs).explore()
Load RGB¶
Select the RGB bands, we see that the raster has only 3 channels now (see Shape in the output of the cell below). The ei object has an attribute called wavelengths with the central wavelength of the hyperspectral band.
import numpy as np
wavelengths_read = np.array([640, 550, 460])
bands_read = np.argmin(np.abs(wavelengths_read[:, np.newaxis] - ei.wavelengths), axis=1).tolist()
ei_rgb = ei.read_from_bands(bands_read)
ei_rgb
Normalize radiance to reflectance¶
import pandas as pd
import matplotlib.pyplot as plt
from georeader import reflectance
thuiller = reflectance.load_thuillier_irradiance() # (dataframe with 8191 rows, 2 colums -> Nanometer, Radiance(mW/m2/nm)
ei_rgb.wavelengths, ei_rgb.fwhm # (K,) vectors with the center wavelength and the FWHM
response = reflectance.srf(ei_rgb.wavelengths, ei_rgb.fwhm, thuiller["Nanometer"].values) # (8191, K)
colors = ["red","green", "blue"]
fig, ax = plt.subplots(1,1,figsize=(12,4))
ax.plot(thuiller["Nanometer"].values,
thuiller["Radiance(mW/m2/nm)"].values,
label="Thuiller irradiance")
ax.set_ylabel("Solar irradiance (mW/m$^2$/nm)")
ax = ax.twinx()
for k in range(3):
ax.plot(thuiller["Nanometer"].values, response[:, k],
label=f"{wavelengths_read[k]}nm",c=colors[k])
ax.set_xlabel("Wavelength nm")
ax.set_ylabel("SRF")
ax.set_xlim(410,700)
ax.legend(loc="upper left")
# solar_irradiance_norm = thuiller["Radiance(mW/m2/nm)"].values.dot(response) # mW/m$^2$/nm
# solar_irradiance_norm/=1_000 # W/m$^2$/nm
# solar_irradiance_norm
# ei_rgb_local = ei_rgb.load(as_reflectance=False)
# ei_rgb_local = reflectance.radiance_to_reflectance(ei_rgb_local, solar_irradiance_norm,
# ei.time_coverage_start,units=units)
ei_rgb_local = ei_rgb.load(as_reflectance=True)
from georeader.plot import show
show((ei_rgb_local).clip(0,1),
mask=ei_rgb_local.values == ei_rgb_local.fill_value_default)
Reproject to UTM¶
import georeader
crs_utm = georeader.get_utm_epsg(ei.footprint("EPSG:4326"))
emit_image_utm = ei.to_crs(crs_utm)
emit_image_utm_rgb = emit_image_utm.read_from_bands(bands_read)
emit_image_utm_rgb_local = emit_image_utm_rgb.load(as_reflectance=True)
emit_image_utm_rgb_local
emit_image_utm_rgb.observation_date_correction_factor
show((emit_image_utm_rgb_local).clip(0,1),
mask=emit_image_utm_rgb_local.values == emit_image_utm_rgb_local.fill_value_default,
add_scalebar=True)
show(emit_image_utm.elevation(), add_colorbar_next_to=True, add_scalebar=True,
mask=True,title="Elevation")
Subset EMIT image¶
from georeader import read
point_tup = (61.28, 36.21)
ei_subset = read.read_from_center_coords(emit_image_utm_rgb, point_tup,
shape=(200,200),crs_center_coords="EPSG:4326")
ei_subset_local = ei_subset.load(as_reflectance=True)
ei_subset_local
from georeader.plot import add_shape_to_plot
from shapely.geometry import Point
ax = show((ei_subset_local).clip(0,1),
mask=ei_subset_local.values == ei_rgb_local.fill_value_default,
add_scalebar=True)
add_shape_to_plot(Point(*point_tup),crs_shape="EPSG:4326",
crs_plot=ei_subset_local.crs, ax=ax)
import matplotlib.pyplot as plt
ei_rgb_raw = ei_rgb.load_raw(transpose=False)
plt.imshow((ei_rgb_raw/12).clip(0,1))
Subset¶
ei_rgb_subset = ei_subset.load_raw(transpose=False)
plt.imshow((ei_rgb_subset/12).clip(0,1))
Load L2AMask¶
The L2AMask is stored in a separated file. The mask array has the following variables:
emit_image_utm.mask_bands
# mask filtering cloudy pixels
show(emit_image_utm.validmask(), add_scalebar=True,vmin=0, vmax=1,
mask=True,title="Valid Mask")
Load metadata¶
Metadata is stored in a separated file (with suffix _OBS_ instead of _RAD_). In metadata we have the following variables:
emit_image_utm.observation_bands
show(emit_image_utm.sza(), add_colorbar_next_to=True, add_scalebar=True,
mask=True,title="Solar Zenith Angle")
show(emit_image_utm.vza(), add_colorbar_next_to=True, add_scalebar=True,
mask=True,title="Viewing Zenith Angle")
show(emit_image_utm.observation('Path length (sensor-to-ground in meters)'),
add_colorbar_next_to=True, add_scalebar=True,
mask=True,title='Path length (sensor-to-ground in meters)')
show(emit_image_utm.observation('Slope (local surface slope as derived from DEM in degrees)'),
add_colorbar_next_to=True, add_scalebar=True,
mask=True,title='Slope (local surface slope as derived from DEM in degrees)')
Licence¶
The georeader package is published under a GNU Lesser GPL v3 licence
georeader tutorials and notebooks are released under a Creative Commons non-commercial licence.
If you find this work useful please cite:
@article{ruzicka_starcop_2023,
title = {Semantic segmentation of methane plumes with hyperspectral machine learning models},
volume = {13},
issn = {2045-2322},
url = {https://www.nature.com/articles/s41598-023-44918-6},
doi = {10.1038/s41598-023-44918-6},
number = {1},
journal = {Scientific Reports},
author = {Růžička, Vít and Mateo-Garcia, Gonzalo and Gómez-Chova, Luis and Vaughan, Anna, and Guanter, Luis and Markham, Andrew},
month = nov,
year = {2023},
pages = {19999},
}