Skip to content

Raster

Raster analysis and parsing module.

Info

Grabs raster information.

Source code in teledetection/upload/raster.py 
 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
class Info:
    """Grabs raster information."""

    def __init__(self, raster_file):
        """Init Info class.

        Args:
            raster_file: str

        Returns:
            Info class

        """
        self.raster_file = raster_file
        with ropen(self.raster_file) as src:
            bbox = src.bounds
            geom = bbox_to_geom(bbox)
            # Reproject the geometry to "epsg:4326"
            geom = warp.transform_geom(src.crs, EPSG_4326, geom)
            self.bbox = features.bounds(geom)
            self.geom = bbox_to_geom(self.bbox)
            self.meta = src.meta
            self.gsd = src.res[0]
            self.proj_ext_info = get_projection_info(src)
            self.nodata = src.nodata
            self.area_or_point = src.tags().get("AREA_OR_POINT", "").lower()
            self.bands = src.indexes
        self.stats = None

    def band_info(self, band: int):
        """Get band info.

        Args:
            band: band index

        Returns:
            band metadata and band statistics

        """
        if band <= 0:
            raise ValueError('The "band" parameter value starts at 1')

        with ropen(self.raster_file) as src_dst:
            md = {
                "data_type": src_dst.dtypes[band - 1],
                "scale": src_dst.scales[band - 1],
                "offset": src_dst.offsets[band - 1],
            }
            if self.area_or_point:
                md["sampling"] = self.area_or_point

            # If the Nodata is not set we don't forward it.
            if src_dst.nodata is not None:
                if numpy.isnan(src_dst.nodata):
                    md["nodata"] = "nan"
                elif numpy.isposinf(src_dst.nodata):
                    md["nodata"] = "inf"
                elif numpy.isneginf(src_dst.nodata):
                    md["nodata"] = "-inf"
                else:
                    md["nodata"] = src_dst.nodata

            if src_dst.units[band - 1] is not None:
                md["unit"] = src_dst.units[band - 1]

            stats = {}
            try:
                if not self.stats:
                    self.stats = src_dst.stats(approx=True)
                statistics = self.stats[band - 1]
                stats.update(
                    {
                        "mean": statistics.mean,
                        "minimum": statistics.min,
                        "maximum": statistics.max,
                        "stddev": statistics.std,
                    }
                )
            except rasterio.errors.StatisticsError as e:
                logger.warning("Unable to compute relevant statistics: %s", e)

            return md, stats

__init__(raster_file)

Init Info class.

Parameters:

Name Type Description Default
raster_file

str

 required

Returns:

Type Description

Info class
Source code in teledetection/upload/raster.py 
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
def __init__(self, raster_file):
    """Init Info class.

    Args:
        raster_file: str

    Returns:
        Info class

    """
    self.raster_file = raster_file
    with ropen(self.raster_file) as src:
        bbox = src.bounds
        geom = bbox_to_geom(bbox)
        # Reproject the geometry to "epsg:4326"
        geom = warp.transform_geom(src.crs, EPSG_4326, geom)
        self.bbox = features.bounds(geom)
        self.geom = bbox_to_geom(self.bbox)
        self.meta = src.meta
        self.gsd = src.res[0]
        self.proj_ext_info = get_projection_info(src)
        self.nodata = src.nodata
        self.area_or_point = src.tags().get("AREA_OR_POINT", "").lower()
        self.bands = src.indexes
    self.stats = None

band_info(band)

Get band info.

Parameters:

Name Type Description Default
band int

band index

 required

Returns:

Type Description

band metadata and band statistics
Source code in teledetection/upload/raster.py 
 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
def band_info(self, band: int):
    """Get band info.

    Args:
        band: band index

    Returns:
        band metadata and band statistics

    """
    if band <= 0:
        raise ValueError('The "band" parameter value starts at 1')

    with ropen(self.raster_file) as src_dst:
        md = {
            "data_type": src_dst.dtypes[band - 1],
            "scale": src_dst.scales[band - 1],
            "offset": src_dst.offsets[band - 1],
        }
        if self.area_or_point:
            md["sampling"] = self.area_or_point

        # If the Nodata is not set we don't forward it.
        if src_dst.nodata is not None:
            if numpy.isnan(src_dst.nodata):
                md["nodata"] = "nan"
            elif numpy.isposinf(src_dst.nodata):
                md["nodata"] = "inf"
            elif numpy.isneginf(src_dst.nodata):
                md["nodata"] = "-inf"
            else:
                md["nodata"] = src_dst.nodata

        if src_dst.units[band - 1] is not None:
            md["unit"] = src_dst.units[band - 1]

        stats = {}
        try:
            if not self.stats:
                self.stats = src_dst.stats(approx=True)
            statistics = self.stats[band - 1]
            stats.update(
                {
                    "mean": statistics.mean,
                    "minimum": statistics.min,
                    "maximum": statistics.max,
                    "stddev": statistics.std,
                }
            )
        except rasterio.errors.StatisticsError as e:
            logger.warning("Unable to compute relevant statistics: %s", e)

        return md, stats

NoSpatialLayerException

Bases: Exception

No spatial layer exception.

Source code in teledetection/upload/raster.py 
24
25
class NoSpatialLayerException(Exception):
    """No spatial layer exception."""

OffsetsInputFormatException

Bases: Exception

No spatial layer exception.

Source code in teledetection/upload/raster.py 
32
33
class OffsetsInputFormatException(Exception):
    """No spatial layer exception."""

ScalesInputFormatException

Bases: Exception

No spatial layer exception.

Source code in teledetection/upload/raster.py 
28
29
class ScalesInputFormatException(Exception):
    """No spatial layer exception."""

apply_created_metadata(item)

Apply time metadata to item.

Source code in teledetection/upload/raster.py 
332
333
334
335
336
def apply_created_metadata(item: Item):
    """Apply time metadata to item."""
    # Timestamps
    logger.debug("Apply created metadata to item")
    item.common_metadata.created = datetime.now()

apply_proj_extension(asset)

Apply projection extension.

Parameters:

Name Type Description Default
asset Asset

cog raster asset

 required

Returns:

Type Description

None
Source code in teledetection/upload/raster.py 
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
def apply_proj_extension(asset: Asset):
    """Apply projection extension.

    Args:
        asset: cog raster asset

    Returns:
        None
    """
    # Projection
    proj_ext = AssetProjectionExtension.ext(asset, add_if_missing=True)
    try:
        info = Info(raster_file=asset.href)
    except NoSpatialLayerException:
        logger.warning("Failed to retrieve spatial info for %s", asset.href)
    proj_ext_args = info.proj_ext_info
    logger.debug("Projection extension args: %s", proj_ext_args)
    if proj_ext_args is not None:
        logger.debug("Applying projection extension")
        proj_ext.apply(**proj_ext_args)

apply_raster_extension(asset)

Apply raster extension.

Parameters:

Name Type Description Default
asset Asset

cog raster asset

 required

Returns:

Type Description

None
Source code in teledetection/upload/raster.py 
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
def apply_raster_extension(asset: Asset):
    """Apply raster extension.

    Args:
        asset: cog raster asset

    Returns:
        None
    """
    # Raster
    # For the raster extension, we merge RasterBand arguments from
    # (1) what is implemented in the `task`, and (2) stacflow built-in
    # properties (stats, nodata, etc)
    raster_ext_args = get_args_for_raster_ext(asset.href)
    logger.debug("Built-in raster extension args: %s", raster_ext_args)
    if raster_ext_args is not None:
        logger.debug("Applying raster extension")
        try:
            raster_ext = RasterExtension.ext(asset)
            if raster_ext.bands:
                builtin_bands = raster_ext_args["bands"]

                # When a raster extension has already been applied to
                # the asset, we want to keep the information that is
                # already there.
                logger.debug(
                    "Found an existing raster extension for this asset. Merging %s and %s",
                    [band.to_dict() for band in builtin_bands],
                    [band.to_dict() for band in raster_ext.bands],
                )
                raster_ext_args = {
                    "bands": [
                        _merge_raster_bands(b1, b2)
                        for b1, b2 in zip(
                            builtin_bands,  # built-in
                            raster_ext.bands,  # in `task`
                        )
                    ]
                }
                logger.debug("Result: %s", raster_ext_args)
            else:
                logger.debug("No raster extension properties for asset")
        except ExtensionNotImplemented:
            logger.debug("No existing raster extension found. Creating.")
            raster_ext = RasterExtension.ext(asset, add_if_missing=True)

        raster_ext.apply(**raster_ext_args)

convert_to_cog(local_filename, keep_cog_dir='')

Convert raster to COG in a tmp directory.

Parameters:

Name Type Description Default
local_filename str

input file path

 required
keep_cog_dir str

path to the directory to keep COG files (not used when "")

 ''

Returns:

Type Description
str

file path of the converted raster
Source code in teledetection/upload/raster.py 
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
def convert_to_cog(local_filename: str, keep_cog_dir: str = "") -> str:
    """Convert raster to COG in a tmp directory.

    Args:
        local_filename: input file path
        keep_cog_dir: path to the directory to keep COG files (not used when "")

    Returns:
        file path of the converted raster
    """
    logger.debug(
        "Input raster file: %s (keep COG dir = %s)", local_filename, keep_cog_dir
    )
    cog_dir = keep_cog_dir or os.path.join(os.path.dirname(local_filename), "TMPCOG")
    logger.debug("Output directory for COG file: %s", cog_dir)
    pathlib.Path(cog_dir).mkdir(parents=True, exist_ok=True)
    cog_file = os.path.join(cog_dir, os.path.basename(local_filename))
    cog_file, _ = os.path.splitext(cog_file)
    cog_file += ".tif"
    if is_raster(cog_file) and is_cog(cog_file) and keep_cog_dir:
        return cog_file
    raster2cog(local_filename, cog_file)
    return cog_file

get_args_for_raster_ext(asset_path)

This function returns the arguments, as a dict, for the raster extension.

Parameters:

Name Type Description Default
asset_path str

asset path

 required

Returns:

Type Description

Dict of arguments for the RasterExtension.ext(...).apply()

function.
Source code in teledetection/upload/raster.py 
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
def get_args_for_raster_ext(asset_path: str):
    """This function returns the arguments, as a dict, for the raster extension.

    Args:
        asset_path: asset path

    Returns:
        Dict of arguments for the `RasterExtension.ext(...).apply()`
        function.
    """
    try:
        raster_info = Info(asset_path)
        bands = []
        for band in raster_info.bands:
            md, stats = raster_info.band_info(band=band)

            raster_stats = Statistics.create(**stats)
            bands.append(RasterBand.create(statistics=raster_stats, **md))
        return {"bands": bands}
    except (rasterio.errors.RasterioIOError, rasterio.errors.CRSError) as err:
        logger.warning(
            "Failed to retrieve raster information for %s (%s). "
            "Maybe it is expected because the file is not a raster.",
            err,
            asset_path,
        )
    return None

is_cog(src_raster)

Check if the raster is in COG format.

Parameters:

Name Type Description Default
src_raster str

source raster file path

 required

Returns:

Type Description
bool

True if the raster is a COG, else False
Source code in teledetection/upload/raster.py 
160
161
162
163
164
165
166
167
168
169
170
def is_cog(src_raster: str) -> bool:
    """Check if the raster is in COG format.

    Args:
        src_raster: source raster file path

    Returns:
        True if the raster is a COG, else False
    """
    is_valid, _, _ = cog_validate(src_raster, quiet=True, strict=True)
    return is_valid

is_raster(src_filepath)

Check if the file is a raster.

Parameters:

Name Type Description Default
src_filepath str

input file path

 required

Returns:

Type Description
bool

True if the provided file is a raster
Source code in teledetection/upload/raster.py 
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
def is_raster(src_filepath: str) -> bool:
    """Check if the file is a raster.

    Args:
        src_filepath: input file path

    Returns:
        True if the provided file is a raster
    """
    try:
        Info(src_filepath)
        return True
    except (errors.RasterioIOError, errors.CRSError):
        pass
    return False

raster2cog(src_raster, dst_raster)

Convert a raster to Cloud Optimized Geotiff.

Parameters:

Name Type Description Default
src_raster str

source raster file path

 required
dst_raster str

destination raster file path

 required
Source code in teledetection/upload/raster.py 
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
def raster2cog(src_raster: str, dst_raster: str):
    """Convert a raster to Cloud Optimized Geotiff.

    Args:
        src_raster: source raster file path
        dst_raster: destination raster file path
    """
    profile = {
        "driver": "GTiff",
        "interleave": "pixel",
        "tiled": True,
        "blockxsize": 512,
        "blockysize": 512,
        "compress": "DEFLATE",
        "BIGTIFF": "IF_SAFER",
    }

    config = {
        "GDAL_NUM_THREADS": os.environ.get("GDAL_NUM_THREADS") or "ALL_CPUS",
        "GDAL_TIFF_INTERNAL_MASK": True,
        "GDAL_TIFF_OVR_BLOCKSIZE": "512",
    }

    assert dst_raster.lower().endswith(".tif"), (
        "Destination COG file extension is not .tif"
    )

    logger.debug("Converting %s to COG file %s", src_raster, dst_raster)
    cog_translate(
        source=src_raster,
        dst_path=dst_raster,
        dst_kwargs=profile,
        config=config,
        web_optimized=False,
        progress_out=None,
        in_memory=False,
        allow_intermediate_compression=True,
    )