Skip to content

19 - Multispectral and hyperspectral imaging

In this chapter one maize field is imaged in RGB, as Sentinel-2 and hyperspectrally, and a vegetation index and the red-edge position are computed.

The scene

import numpy as np
import less

scene = less.Scene()
scene.size = 10.0
scene.terrain = less.Terrain(property=less.Lambertian(reflectance=0.15))
scene.illumination = less.Illumination(
    source=less.Sun(zenith=30, azimuth=150),
    atmosphere=less.SimpleSpectralAtmosphere(turbidity=2.0))

maize = less.Object("maize", mesh=less.asset_path("maize.obj"))
maize.set_property(less.Prospect(cab=45, car=10, cw=0.012, cm=0.006, N=1.55))
positions = less.place.grid(scene, spacing=(0.75, 0.40), jitter=0.04, margin=0.3, seed=42)
rng = np.random.default_rng(42)
scene.add(maize, positions=positions,
          scales=rng.uniform(0.9, 1.1, len(positions)),
          rotations=rng.uniform(0, 360, len(positions)))

RGB

rgb = scene.simulate(less.OpticalImager(
    less.Orthographic(image_size=128), bands=[650, 550, 450], quality=16))
rgb.save("maize_rgb.png")

Satellite multispectral

LESS includes the spectral response functions of Sentinel-2A MSI and Landsat 8 OLI:

bands = less.Sentinel2A()
print(bands.band_names)      # B1 … B12
print(bands.band_centers)    # band centres (nm)

s2 = scene.simulate(less.OpticalImager(
    less.Orthographic(image_size=128), bands=bands,
    spectral_resolution=20, quality=16))
s2.save("maize_sentinel2.tif")

Each band is integrated across its spectral response, which matters especially for broad bands. spectral_resolution is the spectral sampling interval of the integration (nm): 20 nm keeps this example fast; for quantitative work reduce it to 10 nm and 5 nm until the results stop changing. Smaller intervals mean more wavelengths to compute.

Find the red and near-infrared bands by name and compute NDVI:

data = s2.to_brf().data
red = data[:, :, bands.band_names.index("B4")]
nir = data[:, :, bands.band_names.index("B8")]
ndvi = (nir - red) / (nir + red)
print(f"mean NDVI: {ndvi.mean():.3f}")

Typical NDVI ranges: bare soil 0.05–0.15, sparse vegetation 0.15–0.4, dense vegetation 0.4–0.9.

less.Landsat8_OLI() is used the same way.

Hyperspectral

Hyperspectral creates evenly spaced contiguous bands, including both ends:

hyper_bands = less.Hyperspectral(start=400, stop=1000, step=10)
hyper = scene.simulate(less.OpticalImager(
    less.Orthographic(image_size=64), bands=hyper_bands, quality=16))
hyper.save("maize_hyperspectral.tif")

wavelengths = hyper_bands.get_band_centers()
mean_spectrum = hyper.to_brf().data.mean(axis=(0, 1))
print(len(wavelengths), "bands")

More bands take more time and memory.

Red-edge position

The red edge (about 700–750 nm) is the steep rise of vegetation spectra from strong red absorption to high near-infrared reflectance. Its position is often taken as the wavelength of the largest first derivative:

edge_bands = less.Hyperspectral(start=650, stop=800, step=2)
edge = scene.simulate(less.OpticalImager(
    less.Orthographic(image_size=64), bands=edge_bands, quality=32))

edge_wavelengths = edge_bands.get_band_centers()
spectrum = edge.to_brf().data.mean(axis=(0, 1))
derivative = np.gradient(spectrum, edge_wavelengths)
print(f"red-edge position: {edge_wavelengths[np.argmax(derivative)]:.0f} nm")

Illumination and atmosphere

This chapter uses the simple spectral atmosphere, which covers 300–2500 nm. For radiance at airborne or satellite altitude use less.Atmosphere with observation_altitude_km or observation_level="toa" (see chapter 16).

Companion script: scripts/19_multispectral_imaging.py

Next