---
jupytext:
  text_representation:
    extension: .md
    format_name: myst
    format_version: 0.13
    jupytext_version: 1.19.5
kernelspec:
  display_name: Python 3 (ipykernel)
  language: python
  name: python3
---

# Explore HARP2 L2 V4.0 Land Products (FastMAPOL)

**Authors:** Meng Gao (NASA/SSAI), Skye Caplan (NASA/SSAI)

[edl]: https://urs.earthdata.nasa.gov/
[oci-data-access]: https://oceancolor.gsfc.nasa.gov/resources/docs/tutorials/notebooks/oci_data_access/

## Summary

This notebook explores the HARP2 Level-2 (L2) land products generated by FastMAPOL. For more detailed information about the algorithm and data products, please refer to the [FastMAPOL ATBD](https://fastmapol.github.io/mapol-doc/), including the [data product list](https://fastmapol.github.io/mapol-doc/chapters/fastmapol_product_general.html) & [land surface model](https://fastmapol.github.io/mapol-doc/chapters/fastmapol_model_land_summary.html).

You should already be familiar with exploring the FastMAPOL aerosol products from the other notebooks. Beginning with Version 4.0, FastMAPOL provides a new HARP2 aerosol and land surface product under the `MAPOL_LAND` data suite. This product is still undergoing evaluation and should be used with caution, as the retrievals may exhibit significant biases.

In this notebook, we focus on the retrieved land surface reflectance and its spectral and angular dependencies.

## How to Cite

If you use the PACE HARP2 MAPOL_LAND Version 4.0 data, please refer to the dataset information page for the latest citation and DOI:

https://www.earthdata.nasa.gov/data/catalog/ob-cloud-pace-harp2-l2-mapol-land-4.0

The dataset may be cited as:

> NASA Ocean Biology Processing Group. (2026). *PACE HARP2 Level-2 Regional Aerosol Over LAND Optical Properties, FastMAPOL Algorithm, Version 4.0*. NASA Ocean Biology Distributed Active Archive Center (OB.DAAC). DOI: 10.5067/PACE/HARP2/L2/MAPOL_LAND/4.0. Accessed on: YYYY-MM-DD.

## Learning Objectives

By the end of this notebook, you will understand:

* How to acquire HARP2 L2 data
* What land surface products are available
* How to visualize land surface reflectance
* How to examine the angular dependence of surface reflectance and BRDF correction
* How to evaluate basic retrieval quality metrics

## 1. Setup

Begin by importing all the packages used in this notebook. If your kernel uses an environment defined following the guidance in the [tutorials], the imports should be successful.

[tutorials]: https://oceancolor.gsfc.nasa.gov/resources/docs/tutorials/

```{code-cell} ipython3
import math
from pathlib import Path

import earthaccess
import numpy as np
import pandas as pd
import xarray as xr
```

```{code-cell} ipython3
auth = earthaccess.login(persist=True)
```

## 2. Get Level-2 Data

+++

HARP2 L2 data are available through both OB.DAAC and the Earthdata Cloud. Please refer to the L1C notebook for additional information on accessing data from the cloud. The following example retrieves a single HARP2 L2 `MAPOL_LAND` Version 4.0 granule.

Download the HARP2 `MAPOL_LAND` data over Rail Road Valley (38.4958, -115.6964) as one of the validaiton site part of RadCalNet. Further analysis on spectral and anglar information will be discussed in later section.

```{code-cell} ipython3
results = earthaccess.search_data(
    short_name="PACE_HARP2_L2_MAPOL_LAND",
    temporal=("2024-09-28T19:55:26", "2024-09-28T19:55:27"),
    granule_name='*V4_0*',
    count=1,
)
paths = earthaccess.open(results)
```

```{code-cell} ipython3
:tags: [remove-cell]

# this cell is tagged to be removed from HTML renders,
# but we currently want to download when we don't have direct access
if not earthaccess.__store__.in_region:
    paths = earthaccess.download(results, "./")
```

```{code-cell} ipython3
:scrolled: true

datatree = xr.open_datatree(paths[0])
datatree
```

Here, we merge all data groups for convenience in subsequent data manipulation and mark the latitude and longitude variables as coordinates.

```{code-cell} ipython3
dataset = xr.merge(datatree.to_dict().values())
dataset=dataset.set_coords(("latitude", "longitude"))
dataset
```

## 3. Understanding the HARP2 L2 Product Structure

+++

The HARP2 FastMAPOL L2 `MAPOL_LAND` product includes aerosol optical properties for both fine and coarse aerosol modes, as well as retrieved land surface properties.

Land surface product variables include:

* Land surface reflectance (`rhos*`)
* BRDF parameters based on the Ross–Li model
* Land NDVI
* Land white-sky albedo

Multi-angle land surface reflectance is derived by applying atmospheric correction independently to each HARP2 viewing angle. `rhos_angular` represents the retrieved surface reflectance before BRDF correction, while `rhos_nadir` represents the reflectance after BRDF correction and adjustment to nadir viewing geometry. Angular means and standard deviations are also provided as `rhos_nadir_mean`, `rhos_nadir_std`, `rhos_angular_mean`, and `rhos_angular_std`.

The following lists variables related to the land surface model (names beginning with `land`) and land surface reflectance (names beginning with `rhos`).

```{code-cell} ipython3
vars_keep = [
    v for v in dataset.data_vars
    if v.lower().startswith(("land", "rhos"))
]

dataset_sub = dataset[vars_keep]
dataset_sub
```

## 4. Visualize HARP2 L2 Land Surface Properties

+++

In this example, we visualize the retrieved land surface reflectance. We first read the angular means and standard deviations before and after BRDF correction.

```{code-cell} ipython3
rhos_angular_mean = dataset["rhos_angular_mean"].values
rhos_angular_std = dataset["rhos_angular_std"].values

rhos_nadir_mean = dataset["rhos_nadir_mean"].values
rhos_nadir_std = dataset["rhos_nadir_std"].values

rhos_nadir_mean.shape
```

We also read the spatial coordinates and wavelength information.

```{code-cell} ipython3
lat = dataset["latitude"].values
lon = dataset["longitude"].values
plot_range = [lon.min(), lon.max(), lat.min(), lat.max()]
wavelength = dataset["wavelength"].values
wavelength
```

### Define helper functions for visualization

+++

The following helper functions are used to visualize single-band and RGB maps, examine the angular dependence of surface reflectance, and calculate reflectance ratios between selected viewing angles.

```{code-cell} ipython3
import matplotlib.pyplot as plt
from matplotlib.colors import LogNorm
import cartopy.crs as ccrs
import cartopy.feature as cfeature


def _prepare_rgba(data, vmin=None, vmax=None):
    """Normalize RGB/RGBA data and make invalid pixels transparent."""
    data = np.asarray(data, float)

    rgb = data[..., :3].copy()
    alpha = data[..., 3].copy() if data.shape[-1] == 4 else np.ones(data.shape[:2])

    valid = np.all(np.isfinite(rgb), axis=-1) & np.isfinite(alpha)

    if vmin is not None or vmax is not None:
        vmin1 = 0 if vmin is None else vmin
        vmax1 = np.nanmax(rgb) if vmax is None else vmax
        rgb = (rgb - vmin1) / (vmax1 - vmin1)

    rgb = np.clip(np.nan_to_num(rgb), 0, 1)
    alpha = np.where(valid, alpha, 0)

    return np.dstack([rgb, np.clip(alpha, 0, 1)])


def plot_l2_product(
    lon, lat, data, plot_range, label="", title="",
    vmin=None, vmax=None, figsize=(12, 4),
    cmap="viridis", log_scale=False,
    land_color="#f2efe9", ocean_color="#dbe9f6",
):
    """Plot scalar map + histogram, or RGB/RGBA map."""

    data = np.asarray(data)

    # Remove single-band dimension
    if data.ndim == 3 and data.shape[-1] == 1:
        data = data[..., 0]

    is_rgb = data.ndim == 3 and data.shape[-1] in (3, 4)

    # ---------------- Figure ----------------
    if is_rgb:
        fig, ax = plt.subplots(
            figsize=(figsize[0] * 0.75, figsize[1]),
            subplot_kw={"projection": ccrs.PlateCarree()},
        )

    else:
        if data.ndim != 2:
            raise ValueError("data must be 2D scalar or 3/4-band RGB.")

        values = data[np.isfinite(data)]
        if log_scale:
            values = values[values > 0]

        if values.size == 0:
            raise ValueError("No valid data.")

        vmin = np.percentile(values, 2) if vmin is None else vmin
        vmax = np.percentile(values, 98) if vmax is None else vmax

        fig = plt.figure(figsize=figsize)
        gs = fig.add_gridspec(1, 2, width_ratios=[3, 1], wspace=0.3)
        ax = fig.add_subplot(gs[0], projection=ccrs.PlateCarree())

    # ---------------- Background ----------------
    ax.set_extent(plot_range, crs=ccrs.PlateCarree())
    ax.add_feature(cfeature.OCEAN, facecolor=ocean_color, zorder=0)
    ax.add_feature(cfeature.LAND, facecolor=land_color, zorder=0)

    # ---------------- RGB ----------------
    if is_rgb:
        ax.pcolormesh(
            lon, lat, data,
            shading="nearest",
            transform=ccrs.PlateCarree(),
            zorder=1,
        )

    # ---------------- Scalar ----------------
    else:
        norm = LogNorm(vmin=vmin, vmax=vmax) if log_scale else None

        pm = ax.pcolormesh(
            lon, lat, data,
            cmap=cmap,
            norm=norm,
            vmin=None if log_scale else vmin,
            vmax=None if log_scale else vmax,
            shading="auto",
            transform=ccrs.PlateCarree(),
            zorder=1,
        )

        plt.colorbar(pm, ax=ax, pad=0.1).set_label(label)

        # Histogram
        axh = fig.add_subplot(gs[1])

        bins = (
            np.logspace(np.log10(vmin), np.log10(vmax), 40)
            if log_scale else 40
        )

        axh.hist(
            values,
            bins=bins,
            range=None if log_scale else (vmin, vmax),
            color="gray",
            edgecolor="black",
        )

        if log_scale:
            axh.set_xscale("log")

        axh.set(
            xlabel=label,
            ylabel="Count",
            title=f"Histogram: N={values.size}",
        )

    # ---------------- Decoration ----------------
    ax.coastlines(resolution="110m", linewidth=0.8, zorder=3)
    ax.gridlines(draw_labels=True, zorder=4)
    ax.set_title(title)

    plt.show()


########### Plot angular view for a pixel (or within a box)

def plot_rhos_angular(
    dataset, var, iy, ix, n=2,
    title=None, ylim=None
):
    rhos = dataset[var].values[..., 0]
    angle = dataset["sensor_view_angle"].values
    wv = dataset["intensity_wavelength"].values.squeeze()

    # Spatial mean/std
    box = rhos[iy-n:iy+n+1, ix-n:ix+n+1]
    mean = np.nanmean(box, axis=(0, 1))
    std = np.nanstd(box, axis=(0, 1))
    center = rhos[iy, ix]

    # Group by nominal wavelength
    wv_group = np.round(wv)
    wv_unique = np.unique(wv_group[np.isfinite(wv_group)])

    ncols = 2
    nrows = math.ceil(len(wv_unique) / ncols)

    fig, axes = plt.subplots(
        nrows, ncols,
        figsize=(10, 3.5*nrows),
        sharex=True
    )
    axes = np.atleast_1d(axes).ravel()

    for ax, wv1 in zip(axes, wv_unique):
        ind = wv_group == wv1
        order = np.argsort(angle[ind])

        a = angle[ind][order]
        m = mean[ind][order]
        s = std[ind][order]
        c = center[ind][order]

        ax.fill_between(a, m-s, m+s, alpha=0.25)
        ax.plot(a, m, "o-", label=f"{2*n+1}×{2*n+1} mean")
        ax.plot(a, c, "o--", label="Center pixel")

        ax.set_title(f"{wv1:.0f} nm")
        ax.set_xlabel("Sensor View Angle (°)")
        ax.set_ylabel(r"$\rho_s$")
        ax.grid(alpha=0.3)

        if ylim is not None:
            ax.set_ylim(ylim)

    # Remove unused panels
    for ax in axes[len(wv_unique):]:
        ax.remove()

    axes[0].legend()

    fig.suptitle(
        title or
        f"{var}\n"
        f"({dataset['latitude'].values[iy, ix]:.4f}°, "
        f"{dataset['longitude'].values[iy, ix]:.4f}°)"
    )

    plt.tight_layout()
    plt.show()


########### Compute ratios between two angles

def rhos_angle_ratio(dataset, var, iy, ix, angle1, angle2):
    rhos = dataset[var].values[..., 0]
    angle = dataset["sensor_view_angle"].values
    wv = np.round(dataset["intensity_wavelength"].values.squeeze())

    result = []

    for wv1 in np.unique(wv[np.isfinite(wv)]):
        ind = (wv == wv1)

        # Nearest observation within this wavelength band
        j1 = np.where(ind)[0][np.argmin(np.abs(angle[ind] - angle1))]
        j2 = np.where(ind)[0][np.argmin(np.abs(angle[ind] - angle2))]

        r1 = rhos[iy, ix, j1]
        r2 = rhos[iy, ix, j2]

        result.append([
            wv1,
            angle[j1],
            angle[j2],
            r1,
            r2,
            r2 / r1
        ])

    return pd.DataFrame(
        result,
        columns=[
            "wavelength",
            "angle1",
            "angle2",
            "rhos1",
            "rhos2",
            "ratio"
        ]
    )
```

### Visualizations

+++

The following example shows land surface reflectance at a single wavelength. Note that `rhos_nadir_mean` represents the default reflectance to play with, which is the land surface reflectance after BRDF correction and nadir adjustment, followed by averaging over the available viewing angles per band.

```{code-cell} ipython3
wavelength_index = 1
title = "Land surface reflectance: " + str(wavelength[wavelength_index]) + " nm"
label = "Rhos"
data = rhos_nadir_mean[:, :, wavelength_index]

plot_l2_product(
    lon, lat, data,
    plot_range=plot_range,
    label=label,
    title=title,
    vmin=0,
    vmax=0.5,
    cmap="viridis"
)
```

The following example creates an RGB composite of the land surface reflectance. The selected bands are reordered from BGR to RGB, negative reflectances are clipped, and a gamma adjustment is applied for visualization.

```{code-cell} ipython3
# BGR -> RGB
data_rgb = rhos_nadir_mean[:, :, [2, 1, 0]].copy()

# Valid pixels
valid = np.all(np.isfinite(data_rgb), axis=-1)

# Gamma correction and scaling
data_rgb = np.clip(data_rgb, 0, None) ** 0.6
data_rgb = np.clip(data_rgb / 0.5, 0, 1)

# Add transparency
data_rgba = np.dstack([
    data_rgb,
    valid.astype(float),
])

plot_l2_product(
    lon,
    lat,
    data_rgba,
    plot_range=plot_range,
    title="Land surface reflectance (RGB)",
)
```

## 5. Spectral Surface Reflectance

+++

The HARP2 land product also provides spectral information. Here, we examine a location near Railroad Valley at approximately 38.4958°N, 115.6964°W.

+++

### Find the center pixel

+++

The nearest HARP2 pixel to the target latitude and longitude is identified from the two-dimensional geolocation arrays.

```{code-cell} ipython3
lat_target = 38.4958
lon_target = -115.6964

dist2 = (lat - lat_target)**2 + (lon - lon_target)**2

iy, ix = np.unravel_index(np.nanargmin(dist2), dist2.shape)

print("Index:", iy, ix)
print("Lat/Lon:", lat[iy, ix], lon[iy, ix])
```

### Extract the reflectance

+++

To explore the influence of pixel-to-pixel variability, we calculate the mean and standard deviation within a ±n-pixel neighborhood, corresponding to a (2n+1) × (2n+1) pixel box centered on the target location. The center-pixel reflectance is retained for comparison. The variability within the box can be used to assess the homogeneity of the surface properties.

```{code-cell} ipython3
rhos = dataset["rhos_nadir_mean"].values

# explore near by pixels
n = 1
npixel = 2*n+1

y0 = max(0, iy - n)
y1 = min(rhos.shape[0], iy + n + 1)

x0 = max(0, ix - n)
x1 = min(rhos.shape[1], ix + n + 1)

# shape: (ny, nx, wavelength)
rhos_box = rhos[y0:y1, x0:x1, :]

# Center pixel
rhos_target = rhos[iy, ix, :]

print("center pixel rhos:", rhos_target)

# Spatial statistics over the 5x5 neighborhood
rhos_mean = np.nanmean(rhos_box, axis=(0, 1))
rhos_std  = np.nanstd(rhos_box, axis=(0, 1))

print("average pixel rhos:", rhos_mean)
```

### Plot the reflectance spectrum

+++

The figure compares the center-pixel spectrum with the 5 × 5 spatial mean. The shaded region represents ±1 standard deviation within the neighborhood.

```{code-cell} ipython3
plt.figure(figsize=(7, 4))

# Variation among neighboring pixels
plt.fill_between(
    wavelength,
    rhos_mean - rhos_std,
    rhos_mean + rhos_std,
    alpha=0.25,
    label=r"$\pm$1 STD"
)

# Neighborhood mean
plt.plot(
    wavelength, rhos_mean,
    "o-", lw=2,
    label=f"{npixel}×{npixel} mean"
)

# Center pixel
plt.plot(
    wavelength, rhos_target,
    "o--", lw=1.5,
    label="Center pixel"
)

plt.xlabel("Wavelength (nm)")
plt.ylabel(r"Surface Reflectance, $\rho_s$")
plt.title(
    f"Land Surface Reflectance\n"
    f"({lat[iy, ix]:.4f}°, {lon[iy, ix]:.4f}°)"
)

plt.ylim(bottom=0)
plt.grid(alpha=0.3)
plt.legend()
plt.tight_layout()
plt.show()
```

## 6. Analyze the Angular Dependence

+++

The HARP2 land product preserves multi-angle information, allowing the angular dependence of the retrieved surface reflectance to be examined directly. Here, we again use the location near Railroad Valley.

`rhos_angular` represents surface reflectance derived independently at each viewing angle before BRDF correction. The figure below shows the angular variation for each spectral band, together with the spatial variability within the 5 × 5 neighborhood.

```{code-cell} ipython3
plot_rhos_angular(
    dataset,
    "rhos_angular",
    iy, ix,
    title="Angular Surface Reflectance (Before BRDF Correction)",
    ylim=(0, 0.4)
)
```

We can also compare the surface reflectance at two selected viewing angles. The ratio between the two angles provides a simple measure of angular variation. Comparing this ratio among spectral bands can help assess whether the angular dependence is spectrally consistent.

The ratio reported below is defined as the reflectance at `angle2` divided by the reflectance at `angle1`.

```{code-cell} ipython3
angle1, angle2 = 0, 30

ratio = rhos_angle_ratio(
    dataset,
    "rhos_angular",
    iy, ix,
    angle1=angle1,
    angle2=angle2
)

ratio
```

The reverse ratio can be calculated by exchanging the two viewing angles.

```{code-cell} ipython3
angle1, angle2 = 30, 0

ratio = rhos_angle_ratio(
    dataset,
    "rhos_angular",
    iy, ix,
    angle1=angle1,
    angle2=angle2
)

ratio
```

## 7. Angular Correction

+++

Angular correction is applied to `rhos_angular` to produce the nadir-adjusted surface reflectance, `rhos_nadir`. This remains an experimental product under evaluation.

Comparing `rhos_nadir` with `rhos_angular` provides a direct way to examine how the angular correction reduces the angular dependence of the retrieved surface reflectance.

```{code-cell} ipython3
plot_rhos_angular(
    dataset,
    "rhos_nadir",
    iy, ix,
    title="Angular Surface Reflectance (After BRDF Correction)",
    ylim=(0, 0.4)
)
```

The reflectance ratios between selected viewing angles can also be recalculated using `rhos_nadir`. Ratios closer to unity after correction indicate reduced angular dependence.

```{code-cell} ipython3
angle1, angle2 = 0, 30

ratio = rhos_angle_ratio(
    dataset,
    "rhos_nadir",
    iy, ix,
    angle1=angle1,
    angle2=angle2
)

ratio
```

```{code-cell} ipython3
angle1, angle2 = 30, 0

ratio = rhos_angle_ratio(
    dataset,
    "rhos_nadir",
    iy, ix,
    angle1=angle1,
    angle2=angle2
)

ratio
```

For this example, the BRDF correction reduces a substantial portion of the angular variation, although residual angular dependence remains. The performance of this correction is still under evaluation.

+++

## 8. Advanced: Quality Assessment

+++

As with the aerosol products, retrieval quality metrics are important for evaluating the land products, particularly when analyzing multi-angle information.

The primary metrics examined here are the retrieval cost function (`chi2`), the number of retained reflectance measurements (`nv_ref`), the number of retained DoLP (degree of linear polarization) measurements (`nv_dolp`), and the overall `quality_flag`. Please see the [aerosol product tutorial](https://nasa.github.io/oceandata-notebooks/sections/cloud-atmosphere.html) or [ATBD](https://fastmapol.github.io/mapol-doc/chapters/fastmapol_product_quality.html) for more details.

```{code-cell} ipython3
chi2 = dataset["chi2"].values
nv_ref = dataset["nv_ref"].values
nv_dolp = dataset["nv_dolp"].values
quality_flag = dataset["quality_flag"].values
```

The retrieval cost function, $\chi^2$, provides an overall measure of the agreement between the observations and the fitted forward model.

```{code-cell} ipython3
title = r"Retrieval cost function: $\chi^2$"
label = r"$\chi^2$"
data = chi2

plot_l2_product(
    lon, lat, data,
    plot_range=plot_range,
    label=label,
    title=title,
    vmin=0,
    vmax=3,
    cmap="viridis"
)
```

`nv_ref` gives the number of reflectance measurements retained in the retrieval.

```{code-cell} ipython3
title = r"Total number of reflectance measurements"
label = r"$N_{ref}$"
data = nv_ref

plot_l2_product(
    lon, lat, data,
    plot_range=plot_range,
    label=label,
    title=title,
    vmin=0,
    vmax=90,
    cmap="viridis"
)
```

`nv_dolp` gives the number of DoLP measurements retained in the retrieval.

```{code-cell} ipython3
title = r"Total number of DoLP measurements"
label = r"$N_{dolp}$"
data = nv_dolp

plot_l2_product(
    lon, lat, data,
    plot_range=plot_range,
    label=label,
    title=title,
    vmin=0,
    vmax=90,
    cmap="viridis"
)
```

The following provides simple granule-level averages of these retrieval metrics.

```{code-cell} ipython3
np.nanmean(chi2), np.nanmean(nv_ref), np.nanmean(nv_dolp)
```

The overall retrieval quality is summarized by `quality_flag`.

+++

The min and max value of the quality flags:

```{code-cell} ipython3
np.nanmin(data), np.nanmax(data)
```

```{code-cell} ipython3
title = "Retrieval quality flag"
label = "quality_flag"
data = quality_flag

plot_l2_product(
    lon, lat, data,
    plot_range=plot_range,
    label=label,
    title=title,
    vmin=0,
    vmax=6,
    cmap="viridis"
)
```

The `quality_flag` is determined using retrieval metrics including $\chi^2$ and the number of retained measurements. In this example, only a relatively small fraction of the retrievals, primarily near the center of the swath, meet the criteria for the highest-quality category (`quality_flag = 0`). Future improvements in instrument calibration and retrieval performance may increase the fraction of retrievals meeting the highest-quality criteria.

+++

## 9. Advanced: Multi-Angle Data Mask

+++

To better understand the fitting and screening of individual angular measurements, we can examine the adaptive data masks.

`mask_ref` indicates which reflectance measurements are retained or excluded during the retrieval, while `mask_dolp` provides the corresponding information for DoLP measurements.

```{code-cell} ipython3
mask_ref = dataset["mask_ref"].values
mask_dolp = dataset["mask_dolp"].values

mask_ref.shape, mask_dolp.shape
```

The following example shows the reflectance data mask for one viewing-angle index.

```{code-cell} ipython3
angle_index = 5
title = "Adaptive data mask on reflectance: angle index " + str(angle_index)
label = "mask_ref"
data = mask_ref[:, :, angle_index, 0]

plot_l2_product(
    lon, lat, data,
    plot_range=plot_range,
    label=label,
    title=title,
    vmin=0,
    vmax=1,
    cmap="viridis"
)
```

The corresponding DoLP data mask can be examined in the same way.

```{code-cell} ipython3
angle_index = 5
title = "Adaptive data mask on DoLP: angle index " + str(angle_index)
label = "mask_DOLP"
data = mask_dolp[:, :, angle_index, 0]

plot_l2_product(
    lon, lat, data,
    plot_range=plot_range,
    label=label,
    title=title,
    vmin=0,
    vmax=1,
    cmap="viridis"
)
```

Examining these masks for individual viewing angles can help identify where angular observations have been screened from the inversion.

+++

## 10. References

* [FastMAPOL ATBD](https://fastmapol.github.io/mapol-doc/).

* [**Algorithm**] Gao, M., Franz, B. A., Zhai, P.-W., Knobelspiesse, K., Sayer, A. M., Xu, X., Martins, J. V., Cairns, B., Castellanos, P., Fu, G., Hannadige, N., Hasekamp, O., Hu, Y., Ibrahim, A., Patt, F., Puthukkudy, A., and Werdell, P. J.: Simultaneous retrieval of aerosol and ocean properties from PACE HARP2 with uncertainty assessment using cascading neural network radiative transfer models, *Atmos. Meas. Tech.*, **16**, 5863–5881, https://doi.org/10.5194/amt-16-5863-2023, 2023.