# -*- coding: utf-8 -*-
# I2 Bab 2: turunan DEM (kemiringan, arah lereng, bayangan bukit). Penulis: Badar Mubarok Yogaswara
import os
import processing
from qgis.core import QgsRasterLayer

PAKET_I2 = os.environ.get("I2_PAKET", r"C:/KPH_Contoh/paket-i2")     # folder paket-i2
HASIL = r"C:/temp/hasil_i2"
os.makedirs(HASIL, exist_ok=True)

dem = QgsRasterLayer(os.path.join(PAKET_I2, "DEM_Lereng.tif"), "DEM_Lereng")
s = dem.dataProvider().bandStatistics(1)
print("DEM min, maks, rata-rata:", round(s.minimumValue, 1), round(s.maximumValue, 1), round(s.mean, 1))

kemiringan = processing.run("native:slope", {"INPUT": dem, "Z_FACTOR": 1.0, "OUTPUT": os.path.join(HASIL, "slope_derajat.tif")})["OUTPUT"]
arah = processing.run("native:aspect", {"INPUT": dem, "Z_FACTOR": 1.0, "OUTPUT": os.path.join(HASIL, "aspect.tif")})["OUTPUT"]
bayangan = processing.run("native:hillshade", {"INPUT": dem, "Z_FACTOR": 1.0, "AZIMUTH": 315.0, "V_ANGLE": 45.0,
                                               "OUTPUT": os.path.join(HASIL, "hillshade.tif")})["OUTPUT"]

# kemiringan dalam persen = tan(derajat) x 100
lyr = QgsRasterLayer(kemiringan, "slope_derajat")
persen = processing.run("native:rastercalc", {
    "EXPRESSION": 'tan("slope_derajat@1" * 3.14159265358979 / 180) * 100', "LAYERS": [lyr],
    "CELL_SIZE": dem.rasterUnitsPerPixelX(), "EXTENT": dem.extent(), "CRS": dem.crs(),
    "OUTPUT": os.path.join(HASIL, "slope_persen.tif")})["OUTPUT"]

for nama, jalur in (("Slope (derajat)", kemiringan), ("Aspect (derajat dari utara)", arah), ("Hillshade (0-255)", bayangan), ("Slope (persen)", persen)):
    lapisan = QgsRasterLayer(jalur, "x")      # simpan di variabel agar tidak dibuang
    st = lapisan.dataProvider().bandStatistics(1)
    print("%-28s min %.2f  maks %.2f  rata-rata %.2f" % (nama, st.minimumValue, st.maximumValue, st.mean))