# Skrip 4.1: memotong raster dengan poligon (clip by mask). Penulis: Badar Mubarok Yogaswara
import os
import processing
from qgis.core import QgsVectorLayer

DATA = os.environ.get("B3_DATA", r"C:/KPH_Contoh/paket-b3")
OUT = os.environ.get("B3_OUT", r"C:/kerja_b3")
os.makedirs(OUT, exist_ok=True)

# 1) Potong DEM dengan Area_Kajian (poligon tidak beraturan)
hasil = processing.run("gdal:cliprasterbymasklayer", {
    "INPUT": DATA + "/DEM_KPH.tif",
    "MASK": DATA + "/Area_Kajian.gpkg",
    "NODATA": -9999,
    "CROP_TO_CUTLINE": True,
    "KEEP_RESOLUTION": True,
    "OUTPUT": OUT + "/DEM_kajian.tif"})
u = processing.run("native:rasterlayeruniquevaluesreport", {"INPUT": hasil["OUTPUT"], "BAND": 1})
isi = u["TOTAL_PIXEL_COUNT"] - u["NODATA_PIXEL_COUNT"]
print("Piksel total:", u["TOTAL_PIXEL_COUNT"], "| NoData:", u["NODATA_PIXEL_COUNT"])
print("Piksel berisi nilai:", isi, "= luas (ha):", round(isi * 4 / 10000, 2))

# 2) Potong dengan satu KPH saja (filter lapisan batas lebih dulu)
batas = QgsVectorLayer(DATA + "/Batas_KPH.gpkg|layername=Batas_KPH", "Batas_KPH", "ogr")
batas.setSubsetString("\"NAMA_KPH\" = 'KPH Alpha'")
h2 = processing.run("gdal:cliprasterbymasklayer", {
    "INPUT": DATA + "/DEM_KPH.tif", "MASK": batas, "NODATA": -9999,
    "CROP_TO_CUTLINE": True, "KEEP_RESOLUTION": True, "OUTPUT": OUT + "/DEM_Alpha.tif"})
u2 = processing.run("native:rasterlayeruniquevaluesreport", {"INPUT": h2["OUTPUT"], "BAND": 1})
print("KPH Alpha: piksel total", u2["TOTAL_PIXEL_COUNT"], "| NoData", u2["NODATA_PIXEL_COUNT"])