ArcPy zu Rasterio: Vollständiger Migrationsleitfaden für die Raster-Verarbeitung

import arcpy
from arcpy.sa import *

# Environment einrichten
arcpy.env.workspace = "C:/data"
arcpy.CheckOutExtension("Spatial")

# Raster lesen - gibt Raster-Objekt zurück
dem = Raster("elevation.tif")

# Eigenschaften abrufen
print(f"Cell size: {dem.meanCellWidth}")
print(f"Extent: {dem.extent}")
print(f"Spatial ref: {dem.spatialReference.name}")

import rasterio
import numpy as np

# Raster öffnen - Context Manager behandelt Cleanup
with rasterio.open("elevation.tif") as src:
    # Als NumPy Array lesen (Band 1)
    dem = src.read(1)

    # Eigenschaften über Profil abrufen
    print(f"Cell size: {src.res}")
    print(f"Bounds: {src.bounds}")
    print(f"CRS: {src.crs}")

    # dem ist jetzt ein Standard-NumPy Array
    print(f"Shape: {dem.shape}")
    print(f"Dtype: {dem.dtype}")

from arcpy.sa import Slope

# Neigung in Grad berechnen
dem = Raster("elevation.tif")
slope = Slope(dem, "DEGREE")
slope.save("slope.tif")

import rasterio
import numpy as np

with rasterio.open("elevation.tif") as src:
    dem = src.read(1).astype(float)
    cellsize = src.res[0]  # Annahme: quadratische Zellen

    # Gradient mit NumPy berechnen
    dy, dx = np.gradient(dem, cellsize)
    slope_rad = np.arctan(np.sqrt(dx**2 + dy**2))
    slope_deg = np.degrees(slope_rad)

    # Ausgabe schreiben
    profile = src.profile.copy()
    profile.update(dtype=rasterio.float32)
    with rasterio.open("slope.tif", "w", **profile) as dst:
        dst.write(slope_deg.astype(np.float32), 1)

from arcpy.sa import Con, Reclassify, RemapRange

# Bedingt: Werte < 0 auf NoData setzen
dem = Raster("elevation.tif")
dem_positive = Con(dem >= 0, dem)

# Höhe in Kategorien neuklassifizieren
reclass = Reclassify(dem, "VALUE",
    RemapRange([[0, 100, 1],
                [100, 500, 2],
                [500, 1000, 3],
                [1000, 9999, 4]])
)
reclass.save("elevation_classes.tif")

import rasterio
import numpy as np

with rasterio.open("elevation.tif") as src:
    dem = src.read(1).astype(float)
    nodata = src.nodata

    # Bedingt: Negative Werte auf NoData setzen
    dem[dem < 0] = np.nan

    # Mit np.digitize oder np.select neuklassifizieren
    conditions = [
        (dem >= 0) & (dem < 100),
        (dem >= 100) & (dem < 500),
        (dem >= 500) & (dem < 1000),
        dem >= 1000
    ]
    classes = np.select(conditions, [1, 2, 3, 4], default=0)

    # Ausgabe schreiben
    profile = src.profile.copy()
    profile.update(dtype=rasterio.uint8, nodata=0)
    with rasterio.open("elevation_classes.tif", "w", **profile) as dst:
        dst.write(classes.astype(np.uint8), 1)

import rasterio

# Standard GeoTIFF schreiben
with rasterio.open(
    "output.tif", "w",
    driver="GTiff",
    height=data.shape[0],
    width=data.shape[1],
    count=1,
    dtype=data.dtype,
    crs=src.crs,
    transform=src.transform
) as dst:
    dst.write(data, 1)

# Problem: Kein internes Tiling
# Problem: Keine Überblicke
# Problem: Ganze Datei zum Lesen herunterladen

import rasterio
from rasterio.enums import Resampling

# Cloud-Optimized GeoTIFF mit COG-Treiber
with rasterio.open(
    "output_cog.tif", "w",
    driver="COG",
    height=data.shape[0],
    width=data.shape[1],
    count=1,
    dtype=data.dtype,
    crs=src.crs,
    transform=src.transform,
    compress="deflate",        # Kompression
    blocksize=512,             # Größe der internen Fliesen
    overview_resampling="average"
) as dst:
    dst.write(data, 1)

# Oder rio-cogeo CLI verwenden:
# rio cogeo create input.tif output_cog.tif

# NICHT: Metadaten hartcodieren
with rasterio.open(
    "output.tif", "w",
    driver="GTiff",
    height=1000,
    width=1000,
    count=1,
    dtype="float32",
    crs="EPSG:4326",  # Könnte nicht mit Input übereinstimmen!
    transform=Affine(0.1, 0, 0, 0, -0.1, 0)  # Falsch!
) as dst:
    dst.write(data, 1)

# Probleme:
# - CRS könnte nicht mit Input übereinstimmen
# - Transform ist völlig falsch
# - Keine Kompression (riesige Dateien)
# - nodata nicht gesetzt

# TUN: Profil kopieren und aktualisieren
with rasterio.open("input.tif") as src:
    data = src.read(1)
    profile = src.profile.copy()

    # Daten verarbeiten
    result = some_processing(data)

    # Nur das Geänderte aktualisieren
    profile.update(
        driver="COG",           # Cloud-optimiert
        dtype=result.dtype,     # Match output type
        compress="deflate",     # Kompression
        nodata=-9999           # Explizites nodata
    )

    with rasterio.open("output.tif", "w", **profile) as dst:
        dst.write(result, 1)