Click here to download this notebook.

Editings

from casys.report import ReportAltiData

from casys.editing import Clip, ClipCondition, EditingComponent, StatisticsByPass
from casys.readers import CLSTableReader

from casys import CasysPlot, DataParams, DateHandler, Field, NadirData, PlotParams

NadirData.enable_loginfo()

# Reader definition
table_name = "TABLE_C_J3_B_GDRD"
orf_name = "C_J3_GDRD"
cycle_number = 122
start = DateHandler.from_orf(orf_name, cycle_number, 1, pos="first")
end = DateHandler.from_orf(orf_name, cycle_number, 154, pos="last")

reader = CLSTableReader(
    name=table_name,
    date_start=start,
    date_end=end,
    orf=orf_name,
    time="time",
    longitude="LONGITUDE",
    latitude="LATITUDE",
)

# Data container definition
ad = NadirData(source=reader)

var_wind = Field(
    name="wind", source="IIF(FLAG_VAL.ALTI==0, WIND_SPEED.ALTI, DV)", unit="m/s"
)

Define the invalidity conditions

e0 = EditingComponent(
    name="Invalid values",
    invalidity_conditions=[
        ClipCondition(name="Wind speed is NaN", clip=Clip("EQ_DV(WIND_SPEED.ALTI)"))
    ],
    value=1,
)

e1 = EditingComponent(
    name="Quality_indicator",
    invalidity_conditions=[
        ClipCondition(
            name="Sea ice",
            clip=Clip("IIF(FLAG_ICE.ALTI, IIF(FLAG_SURFACE_TYPE :< 2, 1, 0), 0) :!= 0"),
        ),
        ClipCondition(
            name="Land and enclosed seas except caspian",
            clip=Clip(
                "IIF(FLAG_SURFACE_TYPE :== 0 || (FLAG_SURFACE_TYPE :== 1 && OCEANIC_BASIN_NUMBER :== 22), 0, 1) :!= 0"
            ),
        ),
    ],
    value=2,
)

e2 = EditingComponent(
    name="Wind > 5",
    invalidity_conditions=[
        ClipCondition(name="Wind speed out of range", clip=Clip("WIND_SPEED.ALTI :> 5"))
    ],
    value=3,
)

e3_clip = (
    "IIF(BATHYMETRY.MODEL :< -1000 && DISTANCE_SHORELINE.MODEL.ABSOLUTE :> 100000 && OCEANIC_VARIABILITY :< 0.3,"
    "SEA_SURFACE_HEIGHT.ALTI.INTERP - SEA_STATE_BIAS.ALTI.NON_PARAMETRIC - SOLID_EARTH_TIDE_HEIGHT.MODEL.CARTWRIGHT_TAYLER_71"
    "- POLE_TIDE_HEIGHT.MODEL.WAHR_85 - DYNAMICAL_ATMOSPHERIC_CORRECTION.MODEL.MOG2D_HR - DRY_TROPOSPHERIC_CORRECTION.MODEL.ECMWF_GAUSS"
    "- OCEAN_TIDE_HEIGHT.MODEL.GOT4V8 - WET_TROPOSPHERIC_CORRECTION.RAD - IONOSPHERIC_CORRECTION.ALTI - MEAN_SEA_SURFACE.MODEL.CNESCLS11"
    "- SEA_LEVEL_ANOMALY_BIAS.ALTI, DV)"
)

# Invalidate passes where mean or std statistics exceed the thresholds.
e3 = EditingComponent(
    name="Per pass statistics",
    invalidity_conditions=[
        StatisticsByPass(
            name="e3_1",
            clip=Clip(e3_clip),
            nbr_min_pts=3,
            threshold={"mean": 0.3, "std": 0.4},
            orf=ad.orf,
        ),
        StatisticsByPass(
            name="e3_2",
            clip=Clip(e3_clip),
            nbr_min_pts=200,
            threshold={"mean": 0.15, "std": 0.2},
            orf=ad.orf,
        ),
    ],
    value=4,
)

Using the add_editing method:

var_edt = Field(name="EDITING.CUSTOM", description="My personal editing")

ad.add_editing(name="My editing", field=var_edt, components=[e0, e1, e2, e3])

ad.add_raw_data(name="Wind speed", field=var_wind)

ad.compute()

2025-05-14 10:56:23 INFO    Reading ['LONGITUDE', 'LATITUDE', 'wind', 'time', 'DYNAMICAL_ATMOSPHERIC_CORRECTION.MODEL.MOG2D_HR', 'FLAG_ICE.ALTI', 'SOLID_EARTH_TIDE_HEIGHT.MODEL.CARTWRIGHT_TAYLER_71', 'BATHYMETRY.MODEL', 'POLE_TIDE_HEIGHT.MODEL.WAHR_85', 'OCEANIC_BASIN_NUMBER', 'SEA_SURFACE_HEIGHT.ALTI.INTERP', 'OCEANIC_VARIABILITY', 'OCEAN_TIDE_HEIGHT.MODEL.GOT4V8', 'SEA_STATE_BIAS.ALTI.NON_PARAMETRIC', 'DRY_TROPOSPHERIC_CORRECTION.MODEL.ECMWF_GAUSS', 'WET_TROPOSPHERIC_CORRECTION.RAD', 'SEA_LEVEL_ANOMALY_BIAS.ALTI', 'FLAG_SURFACE_TYPE', 'IONOSPHERIC_CORRECTION.ALTI', 'WIND_SPEED.ALTI', 'DISTANCE_SHORELINE.MODEL.ABSOLUTE', 'MEAN_SEA_SURFACE.MODEL.CNESCLS11']
2025-05-14 10:56:31 INFO    Computing My editing
2025-05-14 10:56:32 INFO    Computing done.

Plot

Editings can be visualized in 3 different ways:

• A graphical plot along time using no additional parameter (or plot="time")

• A map plot using the plot="map" parameter

• A graphical view showing which points were invalidated by which condition using the plot="editing" parameter

Each component can be visualized separately using the group parameter.

• group="all": shows all group

• group="valid": shows valid, unedited values

• group="component name": shows only the data flagged by the provided component name

• group=["name1", "name2"]: shows data flagged by the provided components name.

Editings map

cp = CasysPlot(
    data=ad,
    data_name="My editing",
    plot="map",
    group=["Per pass statistics", "Wind > 5"],
)

cp.show()

cp = CasysPlot(data=ad, data_name="My editing", plot="map", group="valid")

cp.show()

Editings graphic

cp = CasysPlot(data=ad, data_name="My editing", plot="editing", group="all")

cp.show()

Merging editing with raw data

cp = CasysPlot(
    data=ad,
    data_name="My editing",
    plot="editing",
    plot_params=PlotParams(
        x_limits=(
            DateHandler("2019/06/02T07:00:00"),
            DateHandler("2019/06/02T08:30:00"),
        ),
        y_limits=(-0.5, 8),
    ),
)

cp2 = CasysPlot(
    data=ad,
    data_name="Wind speed",
    plot_params=PlotParams(y_limits=(-6.4, 21.18)),
    data_params=DataParams(remove_nan=True),
)

cp.add_plot(cp2, shared_ax="x")
cp.show()