Limb Scatter Cloud Retrieval

Limb Scatter Cloud Retrieval#

Here we demonstrate the retrieval of cloud height, width, and optical depth. Similar to other examples we must provide information on the optical parameters.

import numpy as np
import sasktran2 as sk
import skretrieval as skr

In this example we use a Mie distribution with a fixed particle size to define our cloud optical properties.

@skr.Retrieval.register_optical_property("cloud")
def cloud_optical_prop(*args, **kwargs):
    refrac = sk.mie.refractive.H2SO4()
    dist = sk.mie.distribution.LogNormalDistribution().freeze(
        mode_width=1.6, median_radius=80.0,
    )

    db = sk.database.MieDatabase(
        dist,
        refrac,
        np.array([469.0, 750.0]),
    )

    return db

We configure the state vector for the three parameters (vertical_optical_depth, width_fwhm_m, height_m) required by the Gaussian extinction profile constituent.

# Measurement tangent altitudes
tan_alts = np.arange(10000, 45000, 1000)

# Triplets
triplets = {
    "745": {
        "wavelength": [745],
        "weights": [1],
        "altitude_range": [0, 40000],
        "normalization_range": [40000, 45000],
    },
}

# Get the wavelengths required for the triplets
wavel = np.unique(np.concatenate([triplets[t]["wavelength"] for t in triplets])).astype(float)

# Set up a simulated observation with our tangent altitudes, wavelengths, and use 1.5x the initial guess for ozone
obs = skr.observation.SimulatedLimbObservation(
    cos_sza=0.2,
    relative_azimuth=0,
    observer_altitude=200000,
    reference_latitude=20,
    reference_longitude=20,
    tangent_altitudes=tan_alts,
    sample_wavelengths=wavel,
    state_adjustment_factors={
        "cloud": {
            "vertical_optical_depth": 1.5,
            "width_fwhm_m": 1.5,
            "height_m": 1.5,
        }
    },
)

# Construct our measurement vectors
meas_vec = {}
for name, t in triplets.items():
    meas_vec[name] = skr.measvec.Triplet(
        t["wavelength"], t["weights"], t["altitude_range"], t["normalization_range"], normalize=False
    )

# Set up the retrieval object
ret = skr.Retrieval(
    obs,
    measvec=meas_vec,
    minimizer="rodgers",
    target_kwargs={"rescale_state_space": True},
    state_kwargs={
        "altitude_grid": np.arange(0, 70000, 1000),
        "aerosols": {
            "cloud": {
                "type": "gaussian_extinction_profile",
                "nominal_wavelength": 745,
                "retrieved_quantities": {
                    "vertical_optical_depth": {
                        "prior_influence": 1e0,
                        "tikh_factor": 1e-2,
                        "min_value": 1e-9,
                        "max_value": 10
                    },
                    "width_fwhm_m": {
                        "prior_influence": 1e0,
                        "tikh_factor": 1e-2,
                        "min_value": 1e-9,
                        "max_value": 20000
                    },
                    "height_m": {
                        "prior_influence": 1e0,
                        "tikh_factor": 1e-2,
                        "min_value": 0,
                        "max_value": 25000
                    },
                },
                "prior": {
                    "vertical_optical_depth": {"type": "constant", "value": 0.05},
                    "width_fwhm_m": {"type": "constant", "value": 5000},
                    "height_m": {"type": "constant", "value": 15000},
                }
            },
        },
    },
)

# Do the retrieval
results = ret.retrieve()

# Display the results
print(results["state"])

---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
Cell In[3], line 85
   },
)

# Do the retrieval
---> 85 results = ret.retrieve()

# Display the results
print(results["state"])

File ~/checkouts/readthedocs.org/user_builds/skretrieval/checkouts/stable/src/skretrieval/retrieval/processing.py:402, in Retrieval.retrieve(self, enabled_state_elements, enabled_measurement_vectors)
results["meas_l1"] = meas_l1
results["simulated_l1"] = final_l1
--> 402 results["state"] = self._state_vector.describe(min_results)
return self._construct_output(results)

File ~/checkouts/readthedocs.org/user_builds/skretrieval/checkouts/stable/src/skretrieval/retrieval/statevector/altitude.py:109, in AltitudeNativeStateVector.describe(self, rodgers_output)
def describe(self, rodgers_output: dict) -> xr.Dataset:
   """
   Returns a human readable dataset containing the state vector elements

   (...)    107     xr.Dataset
   """
--> 109     result = super().describe(rodgers_output)
   result.coords["altitude"] = self._altitude_grid
   return result

File ~/checkouts/readthedocs.org/user_builds/skretrieval/checkouts/stable/src/skretrieval/retrieval/statevector/__init__.py:129, in StateVector.describe(self, rodgers_output, **kwargs)
end = start + len(state_element.state())
s = slice(start, end)
--> 129 ds = state_element.describe(
   covariance=covar[s, s],
   averaging_kernel=averaging_kernel[s, s],
   **kwargs,
)
if ds is not None:
   all_ds.append(ds)

File ~/checkouts/readthedocs.org/user_builds/skretrieval/checkouts/stable/src/skretrieval/retrieval/statevector/constituent.py:254, in StateVectorElementConstituent.describe(self, **kwargs)
ds[self._constituent_name + "_" + property_name + "_prior"] = (
   xr.DataArray(prior_values, dims=["altitude"])
)
if end - start == 1:  # scalar property
   ds[self._constituent_name + "_" + property_name] = xr.DataArray(
--> 254         float(getattr(self._constituent, property_name))
   )
   ds[self._constituent_name + "_" + property_name + "_prior"] = float(
       self._prior[property_name].state
   )
   if "covariance" in kwargs:

TypeError: only 0-dimensional arrays can be converted to Python scalars