Skip to content

rslearn.train.transforms.terrain

terrain

Transforms deriving terrain products from an elevation band.

ElevationToSlopeAspect

Bases: Transform

Derive slope and aspect from a single-band elevation image.

The input must be a single-channel elevation raster in meters. Because the transform runs after materialization, the image is already in the window's projection, so the pixel spacing is uniform and must be provided in meters via pixel_size_m (a RasterImage does not carry projection metadata). Pixels are assumed to be square, which holds for rslearn windows.

Nodata is assumed to be represented as NaN in the input elevation. NaN pixels produce NaN slope and aspect, and since the gradients read neighboring pixels, each NaN also invalidates the slope and aspect of its four orthogonally adjacent pixels.

Gradients use central differences in the interior and one-sided differences on the image border, so slope and aspect on the outermost pixel ring are less accurate than the interior.

Source code in rslearn/train/transforms/terrain.py 
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
class ElevationToSlopeAspect(Transform):
    """Derive slope and aspect from a single-band elevation image.

    The input must be a single-channel elevation raster in meters. Because the
    transform runs after materialization, the image is already in the window's
    projection, so the pixel spacing is uniform and must be provided in meters via
    ``pixel_size_m`` (a ``RasterImage`` does not carry projection metadata). Pixels
    are assumed to be square, which holds for rslearn windows.

    Nodata is assumed to be represented as NaN in the input elevation. NaN pixels
    produce NaN slope and aspect, and since the gradients read neighboring pixels,
    each NaN also invalidates the slope and aspect of its four orthogonally adjacent
    pixels.

    Gradients use central differences in the interior and one-sided differences on
    the image border, so slope and aspect on the outermost pixel ring are less
    accurate than the interior.
    """

    def __init__(
        self,
        pixel_size_m: float,
        input_selector: str = "elevation",
        output_selector: str | None = None,
        bands: list[str] = list(SUPPORTED_BANDS),
        flat_aspect_value: float = -1.0,
        skip_missing: bool = False,
    ):
        """Initialize a new ElevationToSlopeAspect.

        Args:
            pixel_size_m: the size of a pixel in meters in the window's projection.
            input_selector: the selector containing the elevation image.
            output_selector: the selector to write the result to. Defaults to
                ``input_selector``.
            bands: the bands to emit, in order. Choose from "elevation", "slope",
                and "aspect".
            flat_aspect_value: the aspect value to use for flat pixels, where the
                compass direction of steepest descent is undefined.
            skip_missing: if True, skip the transform when the input selector does
                not exist in the input/target dicts.
        """
        super().__init__(skip_missing=skip_missing)
        if pixel_size_m <= 0:
            raise ValueError("pixel_size_m must be positive")
        if not bands:
            raise ValueError("expected at least one band")
        for band in bands:
            if band not in SUPPORTED_BANDS:
                raise ValueError(
                    f"unsupported band '{band}'; choose from {SUPPORTED_BANDS}"
                )

        self.pixel_size_m = pixel_size_m
        self.input_selector = input_selector
        self.output_selector = (
            output_selector if output_selector is not None else input_selector
        )
        self.bands = list(bands)
        self.flat_aspect_value = flat_aspect_value

    def compute(self, elevation: torch.Tensor) -> torch.Tensor:
        """Compute the configured terrain bands from an elevation tensor.

        Args:
            elevation: THW tensor of elevations in meters. Rows run north to south,
                following the standard raster convention.

        Returns:
            a CTHW tensor with one channel per configured band.
        """
        dem = elevation.to(torch.float32)
        band_arrays: dict[str, torch.Tensor] = {ELEVATION_BAND: dem}

        if SLOPE_BAND in self.bands or ASPECT_BAND in self.bands:
            # torch.gradient returns derivatives per unit spacing, so passing the
            # metric pixel size directly yields elevation change per meter.
            dz_drow = torch.gradient(dem, spacing=self.pixel_size_m, dim=-2)[0]
            dz_dcol = torch.gradient(dem, spacing=self.pixel_size_m, dim=-1)[0]

            # Rows increase southward, so negate to get the northward derivative.
            dz_north = -dz_drow
            dz_east = dz_dcol

            slope = torch.rad2deg(torch.atan(torch.hypot(dz_east, dz_north)))
            band_arrays[SLOPE_BAND] = slope

            if ASPECT_BAND in self.bands:
                # atan2(east, north) is the clockwise bearing of steepest ascent;
                # add 180 degrees for the direction of steepest descent.
                aspect = (torch.rad2deg(torch.atan2(dz_east, dz_north)) + 180.0) % 360.0
                aspect = torch.where(
                    slope < FLAT_SLOPE_EPSILON,
                    torch.full_like(aspect, self.flat_aspect_value),
                    aspect,
                )
                band_arrays[ASPECT_BAND] = aspect

        # Central differences do not read the center pixel, so a pixel with unknown
        # elevation would otherwise get a slope derived purely from its neighbors.
        nan_mask = torch.isnan(dem)
        for band in (SLOPE_BAND, ASPECT_BAND):
            if band in band_arrays:
                band_arrays[band] = torch.where(
                    nan_mask,
                    torch.full_like(band_arrays[band], float("nan")),
                    band_arrays[band],
                )

        return torch.stack([band_arrays[band] for band in self.bands], dim=0)

    def apply_image(self, image: RasterImage) -> RasterImage:
        """Derive terrain bands from the given elevation image.

        Args:
            image: a single-channel CTHW elevation image.

        Returns:
            a CTHW image with one channel per configured band.
        """
        if image.image.shape[0] != 1:
            raise ValueError(
                "expected a single-channel elevation image, got "
                f"{image.image.shape[0]} channels"
            )
        return RasterImage(
            self.compute(image.image[0]),
            timestamps=image.timestamps,
        )

    def forward(
        self, input_dict: dict[str, Any], target_dict: dict[str, Any]
    ) -> tuple[dict[str, Any], dict[str, Any]]:
        """Apply the transform over the inputs and targets.

        Args:
            input_dict: the input
            target_dict: the target

        Returns:
            (input_dict, target_dict) where output_selector contains the terrain
            bands.
        """
        if self.skip_missing and not selector_exists(
            input_dict, target_dict, self.input_selector
        ):
            return input_dict, target_dict

        image = read_selector(input_dict, target_dict, self.input_selector)
        write_selector(
            input_dict, target_dict, self.output_selector, self.apply_image(image)
        )
        return input_dict, target_dict

compute 

compute(elevation: Tensor) -> Tensor

Compute the configured terrain bands from an elevation tensor.

Parameters:

Name Type Description Default
elevation Tensor

THW tensor of elevations in meters. Rows run north to south, following the standard raster convention.

 required

Returns:

Type Description
Tensor

a CTHW tensor with one channel per configured band.
Source code in rslearn/train/transforms/terrain.py 
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
def compute(self, elevation: torch.Tensor) -> torch.Tensor:
    """Compute the configured terrain bands from an elevation tensor.

    Args:
        elevation: THW tensor of elevations in meters. Rows run north to south,
            following the standard raster convention.

    Returns:
        a CTHW tensor with one channel per configured band.
    """
    dem = elevation.to(torch.float32)
    band_arrays: dict[str, torch.Tensor] = {ELEVATION_BAND: dem}

    if SLOPE_BAND in self.bands or ASPECT_BAND in self.bands:
        # torch.gradient returns derivatives per unit spacing, so passing the
        # metric pixel size directly yields elevation change per meter.
        dz_drow = torch.gradient(dem, spacing=self.pixel_size_m, dim=-2)[0]
        dz_dcol = torch.gradient(dem, spacing=self.pixel_size_m, dim=-1)[0]

        # Rows increase southward, so negate to get the northward derivative.
        dz_north = -dz_drow
        dz_east = dz_dcol

        slope = torch.rad2deg(torch.atan(torch.hypot(dz_east, dz_north)))
        band_arrays[SLOPE_BAND] = slope

        if ASPECT_BAND in self.bands:
            # atan2(east, north) is the clockwise bearing of steepest ascent;
            # add 180 degrees for the direction of steepest descent.
            aspect = (torch.rad2deg(torch.atan2(dz_east, dz_north)) + 180.0) % 360.0
            aspect = torch.where(
                slope < FLAT_SLOPE_EPSILON,
                torch.full_like(aspect, self.flat_aspect_value),
                aspect,
            )
            band_arrays[ASPECT_BAND] = aspect

    # Central differences do not read the center pixel, so a pixel with unknown
    # elevation would otherwise get a slope derived purely from its neighbors.
    nan_mask = torch.isnan(dem)
    for band in (SLOPE_BAND, ASPECT_BAND):
        if band in band_arrays:
            band_arrays[band] = torch.where(
                nan_mask,
                torch.full_like(band_arrays[band], float("nan")),
                band_arrays[band],
            )

    return torch.stack([band_arrays[band] for band in self.bands], dim=0)

apply_image 

apply_image(image: RasterImage) -> RasterImage

Derive terrain bands from the given elevation image.

Parameters:

Name Type Description Default
image RasterImage

a single-channel CTHW elevation image.

 required

Returns:

Type Description
RasterImage

a CTHW image with one channel per configured band.
Source code in rslearn/train/transforms/terrain.py 
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
def apply_image(self, image: RasterImage) -> RasterImage:
    """Derive terrain bands from the given elevation image.

    Args:
        image: a single-channel CTHW elevation image.

    Returns:
        a CTHW image with one channel per configured band.
    """
    if image.image.shape[0] != 1:
        raise ValueError(
            "expected a single-channel elevation image, got "
            f"{image.image.shape[0]} channels"
        )
    return RasterImage(
        self.compute(image.image[0]),
        timestamps=image.timestamps,
    )

forward 

forward(input_dict: dict[str, Any], target_dict: dict[str, Any]) -> tuple[dict[str, Any], dict[str, Any]]

Apply the transform over the inputs and targets.

Parameters:

Name Type Description Default
input_dict dict[str, Any]

the input

 required
target_dict dict[str, Any]

the target

 required

Returns:

Type Description
dict[str, Any]

(input_dict, target_dict) where output_selector contains the terrain
dict[str, Any]

bands.
Source code in rslearn/train/transforms/terrain.py 
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
def forward(
    self, input_dict: dict[str, Any], target_dict: dict[str, Any]
) -> tuple[dict[str, Any], dict[str, Any]]:
    """Apply the transform over the inputs and targets.

    Args:
        input_dict: the input
        target_dict: the target

    Returns:
        (input_dict, target_dict) where output_selector contains the terrain
        bands.
    """
    if self.skip_missing and not selector_exists(
        input_dict, target_dict, self.input_selector
    ):
        return input_dict, target_dict

    image = read_selector(input_dict, target_dict, self.input_selector)
    write_selector(
        input_dict, target_dict, self.output_selector, self.apply_image(image)
    )
    return input_dict, target_dict