BAB 16: Kisi Piksel yang Bisa Dihitung: Raster dan Citra dengan Skrip

Studi kasus: "Seberapa hijau KPH kita, dan berapa luas lereng curamnya?" #

Kepala Seksi membawa citra satelit dan peta tinggi KPH Contoh. Pertanyaannya dua: seberapa hijau tiap KPH, dan berapa hektare lahan yang lerengnya curam? Anda bisa menjawabnya dengan banyak klik di QGIS. Dengan skrip, satu kali jalan sudah cukup, dan Anda bisa mengulanginya untuk citra tahun depan. Tujuh skrip di bab ini membawa Anda dari memeriksa isi citra sampai tabel ringkasan per KPH.

Konsep: raster dalam tiga kalimat #

Raster adalah gambar yang tersusun dari kotak-kotak kecil bernama piksel, seperti keramik di lantai. Setiap kotak menyimpan satu angka, misalnya kecerahan cahaya atau tinggi tanah. Satu lembar angka disebut band, dan citra satelit biasanya terdiri dari beberapa band yang ditumpuk.

Lima istilah yang perlu Anda kenal sekarang:

• Piksel: satu kotak raster. Pada data latihan ini, satu piksel mewakili 2 m x 2 m di lapangan.
• Band: satu lembar angka. Citra latihan memakai biru, hijau, merah, dan inframerah dekat.
• DEM (Digital Elevation Model): raster yang angkanya adalah tinggi tanah.
• CRS: sistem koordinat yang dipakai layer. Data latihan memakai UTM zona 49 Selatan (EPSG:32749), satuannya meter.
• NoData: nilai penanda "kotak ini kosong". Kotak seperti ini tidak ikut dihitung.

Peta jalan bab ini:

Skrip 16.1 dan 16.4 berdiri sendiri. Skrip 16.7 membutuhkan hasil 16.2, 16.3, 16.5, dan 16.6, jadi jalankan semuanya berurutan.

Dua baris lokasi di awal tiap skrip. Skrip bab ini membaca lokasi data dan folder hasil dari pengaturan sistem (variabel lingkungan, lihat Bab 8). Bila pengaturan itu tidak ada, skrip memakai alamat cadangan. Alamat cadangan hanya contoh, bukan alamat di komputer Anda:

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)

Ada dua cara menyesuaikannya. Cara pertama, atur variabel lingkungan B3_DATA dan B3_OUT ke lokasi di komputer Anda. Cara kedua, ganti alamat cadangan di ujung baris dengan lokasi folder paket-b3 dan folder hasil Anda. Bila memilih cara kedua, ubah di setiap skrip dengan isi yang sama.

Baris os.makedirs(OUT, exist_ok=True) membuat folder hasil bila belum ada (lihat Bab 8). Anda tidak perlu membuatnya lebih dulu. Skrip 16.5 hanya memakai OUT, jadi tidak punya baris DATA. Skrip 16.7 hanya membaca hasil skrip lain, jadi tidak membuat folder.

Skrip 16.1: Memeriksa isi sebuah citra #

Sebelum memakai citra, kenali dulu seperti Anda mengecek label kemasan. Berkasnya b3_01_info.py; di dalamnya tertulis "Skrip 1.1".

Tujuan. Memeriksa bahwa citra terbaca, lalu menampilkan ukuran, jumlah band, ukuran piksel, CRS, tipe data, dan statistik dasarnya.

Bahan. Paket data paket-b3, berkas Citra_B04.tif (band merah). Skrip memakai processing, alat bawaan Konsol Python QGIS. Cara menjalankannya ada di Bab 3.

Istilah baru di skrip ini:

• QgsRasterLayer: Perkakas pembuka raster.
• tipe data: Jenis angka yang disimpan tiap piksel.
• statistik raster: Angka terkecil, terbesar, dan rata-rata.
• kamus hasil: Alat mengembalikan pasangan nama dan nilai.
• raise SystemExit("pesan"): Menghentikan skrip seketika dan menampilkan pesan. Mirip raise di Bab 10.
• os.environ.get(...) dan import os: Baca pengaturan sistem, lihat Bab 8.

Skrip utuh:

# Skrip 1.1: memeriksa isi sebuah raster (jalankan di Python Console QGIS). Penulis: Badar Mubarok Yogaswara
# Ubah DATA ke folder paket-b3 di komputer Anda.
import os
import processing
from qgis.core import QgsRasterLayer

DATA = os.environ.get("B3_DATA", r"C:/KPH_Contoh/paket-b3")

lay = QgsRasterLayer(DATA + "/Citra_B04.tif", "Citra_B04")
if not lay.isValid():
    raise SystemExit("Raster tidak terbaca. Periksa DATA: " + DATA)
print("Terbaca:", lay.isValid())
print("Ukuran (kolom x baris):", lay.width(), "x", lay.height())
print("Jumlah band:", lay.bandCount())
print("Ukuran piksel (m):", lay.rasterUnitsPerPixelX(), "x", lay.rasterUnitsPerPixelY())
print("CRS:", lay.crs().authid())
print("Tipe data band 1:", lay.dataProvider().dataType(1).name)

hasil = processing.run("native:rasterlayerstatistics", {"INPUT": lay, "BAND": 1})
print("Jumlah piksel:", hasil["COUNT"])
print("Min:", hasil["MIN"], " Maks:", hasil["MAX"], " Rata-rata:", round(hasil["MEAN"], 1))

Penjelasan baris demi baris.

Hasil. Layar menampilkan (uji pada QGIS 4.0.2, dengan DATA diganti ke folder salinan paket data):

Terbaca: True
Ukuran (kolom x baris): 200 x 200
Jumlah band: 1
Ukuran piksel (m): 2.0 x 2.0
CRS: EPSG:32749
Tipe data band 1: UInt16
Jumlah piksel: 40000
Min: 227.0  Maks: 2061.0  Rata-rata: 449.0

Dari hasil ini Anda tahu: citra berukuran 200 x 200 piksel, jadi 40000 kotak. Tiap piksel 2 meter, sehingga citra mencakup 400 m x 400 m, yaitu 16 hektare. Tipe UInt16 berarti bilangan bulat tanpa tanda minus, 16 bit.

Bagian yang boleh diubah.

• Baris 7, lokasi data. Atur B3_DATA atau ganti alamat cadangan ke komputer Anda.
• Baris 9, nama berkas. Ganti Citra_B04.tif untuk memeriksa citra lain. Dari uji: pada Citra_B08.tif (inframerah dekat) hasilnya min 160.0, maks 4878.0, rata-rata 3494.4.
• Baris 17 dan 19, nomor band. Pada citra empat band (Skrip 16.2), ganti angka 1 untuk memeriksa band lain.

Jebakan.

• Jalur salah menghentikan skrip. Bila DATA menunjuk folder yang salah atau kosong, baris 10-11 menghentikan skrip dengan pesan Raster tidak terbaca. Periksa DATA: diikuti isi DATA. Tanpa pengaman itu, skrip terus berjalan dan mencetak angka janggal: Terbaca: False, ukuran 0, dan nilai terkecil serta terbesar berupa angka raksasa. Bila Anda melihat pesan pengaman, periksa baris 7 dan nama berkas di baris 9.
• Nomor band berlebih. Dari uji, mengganti BAND menjadi 2 pada citra satu band berhenti dengan QgsProcessingException: Invalid band number for BAND (2): Valid values for input raster are 1 to 1.
• Satuan tidak dicek. Teks "(m)" di baris 15 hanya tulisan. Skrip tidak memeriksa apakah CRS memakai meter. Pada CRS berderajat, angkanya derajat.

Latihan variasi.

1. Periksa Citra_B08.tif, bukan Citra_B04.tif.
2. Tambahkan satu baris yang mencetak luas citra dalam hektare.

Contoh jawaban. Nomor 1: ganti nama berkas di baris 9 dan nama layer di sebelahnya. Layar menampilkan min 160.0, maks 4878.0, dan rata-rata 3494.4. Nomor 2, tambahkan di akhir skrip:

print("Luas citra (ha):", hasil["COUNT"] * lay.rasterUnitsPerPixelX() * lay.rasterUnitsPerPixelY() / 10000)

Untuk Citra_B04.tif, layar menampilkan Luas citra (ha): 16.0.

Skrip 16.2: Menumpuk empat band menjadi satu citra #

Empat berkas satu-band ibarat empat lembar transparan. Skrip ini menumpuknya menjadi satu berkas empat lembar, supaya mudah dirakit menjadi gambar berwarna atau dihitung bersama. Berkasnya b3_02_komposit.py; di dalamnya tertulis "Skrip 2.1".

Tujuan. Menggabungkan citra biru, hijau, merah, dan inframerah dekat menjadi satu citra empat band bernama Citra_4band.tif.

Bahan. Paket paket-b3: Citra_B02.tif, Citra_B03.tif, Citra_B04.tif, dan Citra_B08.tif. Folder hasil OUT dibuat otomatis bila belum ada. Hasilnya dipakai Skrip 16.6.

Istilah baru di skrip ini:

• os.path.exists dan os.remove: Memeriksa apakah berkas ada, lalu menghapusnya. os.path.exists dibahas di Bab 8.
• komposit: Citra multi-band.
• urutan band: Nomor band ditentukan urutan berkas.
• GDAL: Pustaka pengolah raster di dalam QGIS.
• daftar: List, Bab 5.

Skrip utuh:

# Skrip 2.1: menumpuk empat band menjadi satu raster multiband. Penulis: Badar Mubarok Yogaswara
import os
import processing
from qgis.core import QgsRasterLayer

DATA = os.environ.get("B3_DATA", r"C:/KPH_Contoh/paket-b3")
OUT = os.environ.get("B3_OUT", r"C:/kerja_b3")          # folder hasil (dibuat otomatis bila belum ada)
os.makedirs(OUT, exist_ok=True)

# Urutan penting: band 1 = biru (B02), 2 = hijau (B03), 3 = merah (B04), 4 = inframerah dekat (B08)
urut = [DATA + "/Citra_B02.tif", DATA + "/Citra_B03.tif", DATA + "/Citra_B04.tif", DATA + "/Citra_B08.tif"]
HASIL = OUT + "/Citra_4band.tif"
if os.path.exists(HASIL):
    os.remove(HASIL)          # gdal:merge tidak menimpa hasil lama
hasil = processing.run("gdal:merge", {
    "INPUT": urut,
    "SEPARATE": True,        # satu berkas = satu band
    "DATA_TYPE": 2,          # 2 = UInt16, sama dengan data asal
    "OUTPUT": OUT + "/Citra_4band.tif"})
stack = QgsRasterLayer(hasil["OUTPUT"], "Citra_4band")
if not stack.isValid():
    raise SystemExit("Citra_4band.tif tidak terbentuk. Periksa DATA dan OUT.")
print("Jumlah band:", stack.bandCount(), "| tipe:", stack.dataProvider().dataType(1).name)
print("Warna asli  : merah=band 3, hijau=band 2, biru=band 1")
print("Warna palsu : merah=band 4, hijau=band 3, biru=band 2")

Penjelasan baris demi baris.

Hasil. Layar menampilkan (uji pada QGIS 4.0.2):

Jumlah band: 4 | tipe: UInt16
Warna asli  : merah=band 3, hijau=band 2, biru=band 1
Warna palsu : merah=band 4, hijau=band 3, biru=band 2

Bukti urutan band benar: rata-rata band 1 sampai 4 pada hasil, yaitu 293.9, 607.0, 449.0, dan 3494.4, sama dengan rata-rata berkas B02, B03, B04, dan B08. Keempatnya diukur dengan alat statistik pada Skrip 16.1.

Bagian yang boleh diubah.

• Baris 6-7, lokasi. Atur B3_DATA dan B3_OUT, atau ganti alamat cadangannya.
• Baris 11, urutan dan isi daftar. Urutan menentukan nomor band. Dari uji: membuang Citra_B08.tif menghasilkan Jumlah band: 3.
• Baris 18, DATA_TYPE. Angka adalah nomor pilihan pada alat: 0 Byte, 1 Int16, 2 UInt16, 3 UInt32, 4 Int32, 5 Float32. Urutan ini berlaku untuk gdal:merge. [CEK] pada versi QGIS Anda.
• Baris 12 dan 19, nama berkas hasil. Bila diganti, ubah di kedua baris, dan Skrip 16.6 harus ikut diubah.

Jebakan.

• Urutan band salah tanpa pesan. Dari uji: menukar B08 dan B04 di baris 11 tetap menghasilkan 4 band tanpa galat. Akibatnya terlihat di Skrip 16.6, yaitu NDVI terbalik tanda.
• Berkas hasil lama dihapus otomatis. Baris 13-14 menghapus Citra_4band.tif yang lama setiap skrip dijalankan. Tanpa baris itu, gdal:merge tidak menimpa berkas lama. Contohnya hasil percobaan tiga band tetap dipakai walau Anda memberi empat berkas, dan Skrip 16.6 kemudian gagal. Akibat baris ini: bila Anda menyimpan berkas lain dengan nama yang sama, berkas itu ikut terhapus. Pakai nama lain di baris 12 dan 19 bila ingin menyimpannya.
• Pesan pengaman. Bila baris 22 menghentikan skrip, periksa DATA dan OUT, lalu pastikan keempat berkas Citra_B02, B03, B04, dan B08 ada.

Latihan variasi.

1. Gabungkan hanya tiga band (biru, hijau, merah).
2. Pastikan hasilnya benar dengan mencetak rata-rata band 3.

Contoh jawaban. Nomor 1: hapus , DATA + "/Citra_B08.tif" dari baris 11. Berkas hasil lama dihapus otomatis oleh skrip. Layar menampilkan Jumlah band: 3 | tipe: UInt16. Nomor 2: tambahkan di akhir skrip:

st = processing.run("native:rasterlayerstatistics", {"INPUT": stack, "BAND": 3})
print("Rata-rata band 3:", round(st["MEAN"], 1))

Layar menampilkan Rata-rata band 3: 449.0, sama dengan Citra_B04.tif.

Skrip 16.3: Dari satu peta tinggi menjadi empat turunan #

Satu DEM bisa melahirkan banyak peta: bayangan bukit, kemiringan, arah hadap lereng, dan garis kontur. Skrip ini membuat keempatnya sekaligus. Berkasnya b3_03_dem.py; di dalamnya tertulis "Skrip 3.1".

Tujuan. Dari DEM_KPH.tif, membuat Hillshade.tif, Lereng_derajat.tif, Arah_lereng.tif, dan Kontur_5m.gpkg, lalu mencetak statistik elevasi dan lereng.

Bahan. Paket paket-b3: DEM_KPH.tif. Folder hasil OUT dibuat otomatis bila belum ada. Lereng_derajat.tif dipakai Skrip 16.5.

Istilah baru di skrip ini:

• hillshade: Bayangan bukit, seolah matahari menyinari relief.
• kemiringan: slope, derajat.
• arah lereng: aspect, ke mana lereng menghadap.
• kontur: Garis tinggi sama.
• azimut: Arah datang cahaya, dalam derajat.

Skrip utuh:

# Skrip 3.1: hillshade, kemiringan, arah lereng, dan kontur dari DEM. Penulis: Badar Mubarok Yogaswara
import os
import processing

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)
dem = DATA + "/DEM_KPH.tif"

st = processing.run("native:rasterlayerstatistics", {"INPUT": dem, "BAND": 1})
print("Elevasi min/maks (m):", round(st["MIN"], 2), "/", round(st["MAX"], 2))

processing.run("native:hillshade", {"INPUT": dem, "Z_FACTOR": 1, "AZIMUTH": 300, "V_ANGLE": 40,
                                    "OUTPUT": OUT + "/Hillshade.tif"})
processing.run("native:slope", {"INPUT": dem, "Z_FACTOR": 1, "OUTPUT": OUT + "/Lereng_derajat.tif"})
processing.run("native:aspect", {"INPUT": dem, "Z_FACTOR": 1, "OUTPUT": OUT + "/Arah_lereng.tif"})
processing.run("gdal:contour", {"INPUT": dem, "BAND": 1, "INTERVAL": 5, "FIELD_NAME": "ELEV",
                                "OUTPUT": OUT + "/Kontur_5m.gpkg"})

sl = processing.run("native:rasterlayerstatistics", {"INPUT": OUT + "/Lereng_derajat.tif", "BAND": 1})
print("Lereng min/maks/rata-rata (derajat):", round(sl["MIN"], 2), "/", round(sl["MAX"], 2), "/", round(sl["MEAN"], 2))

Penjelasan baris demi baris.

Hasil. Layar menampilkan (uji pada QGIS 4.0.2):

Elevasi min/maks (m): 40.03 / 74.25
Lereng min/maks/rata-rata (derajat): 0.01 / 21.89 / 6.39

Skrip ini hanya mencetak dua baris. Hasil utamanya berupa empat berkas di folder OUT. Dari uji, Kontur_5m.gpkg berisi 17 garis kontur, dengan kolom fid, ID, dan ELEV. Lereng terbesar 21,89 derajat, jadi lahan ini landai sampai agak curam.

Bagian yang boleh diubah.

• Baris 13, AZIMUTH dan V_ANGLE. Arah dan tinggi cahaya. Ubah untuk menggeser bayangan. Penulis memakai 300 dan 40. Nilai lazim di banyak peta adalah 315 dan 45. [CEK]
• Baris 13, 15, dan 16, Z_FACTOR. Pengali tinggi. Angka 1 benar bila satuan tinggi dan satuan horizontal sama (sama-sama meter).
• Baris 17, INTERVAL. Jarak antar kontur. Dari uji: 10 menghasilkan 6 garis kontur, bukan 17.
• Baris 17, FIELD_NAME. Nama kolom tinggi pada garis kontur.

Jebakan.

• Folder hasil dibuat oleh baris 7. Tanpa os.makedirs, alat pertama yang menulis hasil gagal dengan pesan GDAL ERROR 4 ... No such file or directory. Skrip lalu berhenti dengan Could not load source layer ... not found, yang hanya akibat berkas tidak pernah dibuat. Dengan baris 7, folder hasil terbentuk sendiri.
• Satuan tinggi dan horizontal berbeda. Bila DEM Anda memakai satuan tinggi yang berbeda dari satuan jarak, Z_FACTOR harus disesuaikan. Skrip selalu memakai 1 dan tidak membahasnya. [CEK]
• Nilai kembalian tidak disimpan. Baris 13-18 tidak menyimpan hasil processing.run di kotak. Tidak masalah karena hasilnya sudah ke berkas.
• Hasil lama diganti. Dari uji, menjalankan ulang memperbarui Hillshade.tif (waktu simpan berubah). Berbeda dengan gdal:merge di Skrip 16.2.

Latihan variasi.

1. Buat kontur tiap 10 meter, bukan 5.
2. Ganti arah cahaya ke 315 dengan tinggi 45 derajat.

Contoh jawaban. Nomor 1: ganti 5 menjadi 10 di baris 17 (dan nama Kontur_5m di baris 18 menjadi Kontur_10m agar tidak menimpa). Hasilnya berkas dengan 6 garis kontur. Nomor 2: ganti 300 menjadi 315 dan 40 menjadi 45 di baris 13. Layar tidak berubah, tetapi Hillshade.tif berbeda. Bandingkan keduanya dengan menampilkannya di QGIS.

Skrip 16.4: Memotong raster dengan poligon #

Memotong raster seperti menggunting kain mengikuti pola. Bagian di luar pola dibuang, atau ditandai kosong. Berkasnya b3_04_clip.py; di dalamnya tertulis "Skrip 4.1".

Tujuan. Memotong DEM dengan poligon Area_Kajian, menghitung piksel berisi dan luasnya, lalu memotong DEM dengan satu KPH saja.

Bahan. Paket paket-b3: DEM_KPH.tif, Area_Kajian.gpkg, dan Batas_KPH.gpkg (layer Batas_KPH dengan kolom NAMA_KPH). Folder hasil OUT dibuat otomatis bila belum ada.

Istilah baru di skrip ini:

• mask: Poligon pemotong.
• crop to cutline: Kotak raster dipangkas mengikuti poligon.
• nilai unik: Daftar nilai yang ada beserta jumlahnya.
• setSubsetString: Filter layer, sama dengan Skrip 15.2.

Skrip utuh:

# 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"])

Penjelasan baris demi baris.

Hasil. Layar menampilkan (uji pada QGIS 4.0.2):

Piksel total: 24800 | NoData: 6100
Piksel berisi nilai: 18700 = luas (ha): 7.48
KPH Alpha: piksel total 10000 | NoData 0

Pada sebagian komputer, baris Warning 1: Field NAMA_KPH of width 255 truncated to 254. ikut tampil. Itu peringatan GDAL yang tidak berbahaya. Bacalah hasilnya: poligon kajian meninggalkan 6100 kotak kosong. KPH Alpha berbentuk persegi 200 m x 200 m, jadi 100 x 100 piksel penuh, tanpa NoData.

Bagian yang boleh diubah.

• Baris 6-7, lokasi. Atur B3_DATA dan B3_OUT, atau ganti alamat cadangannya.
• Baris 25, nama KPH. Ganti 'KPH Alpha' dengan 'KPH Beta' atau 'KPH Gamma'. Ejaannya harus sama persis dengan isi kolom.
• Baris 14 dan 27, NODATA. Pilih angka yang tidak mungkin muncul sebagai tinggi sungguhan.
• Baris 15-16, CROP_TO_CUTLINE dan KEEP_RESOLUTION. Ubah ke False untuk melihat akibatnya pada ukuran hasil.
• Baris 21, angka 4. Luas satu piksel dalam meter persegi (2 m x 2 m).

Jebakan.

• Angka 4 tertulis tetap. Faktor luas mengandaikan piksel 2 m. Bila piksel Anda 10 m, luas salah tanpa pesan. Ambil ukuran piksel dari Skrip 16.1.
• Filter ikut terbawa. Baris 26-28 memakai layer yang difilter. Bila Anda menghapus baris 25, semua KPH ikut terpotong. Dari uji: tanpa filter, hasilnya 40000 piksel total (seluruh DEM), karena ketiga KPH menutup seluruh DEM.
• Data tidak ditemukan. Folder hasil dibuat otomatis, tetapi DATA tidak. Bila DATA salah, alat berhenti dengan galat QgsProcessingException: Unable to execute algorithm, disusul Could not load source layer for INPUT: .../DEM_KPH.tif not found (diuji). Folder hasil tetap terbentuk, tetapi kosong. Periksa baris 6.
• Nilai -9999. Bila DEM asli punya nilai itu, ia ikut dianggap kosong. [kemungkinan]
• Filter hanya berlaku di sesi ini. setSubsetString mengubah objek batas, bukan berkas.

Latihan variasi.

1. Potong DEM dengan KPH Beta. Berapa piksel total dan NoData-nya?
2. Cetak luas KPH Alpha dalam hektare dari hasil potongan.

Contoh jawaban. Nomor 1: ganti 'KPH Alpha' di baris 25 dengan 'KPH Beta'. Layar menampilkan piksel total 10000 | NoData 0. Nomor 2: tambahkan di akhir skrip:

print("Luas KPH Alpha (ha):", (u2["TOTAL_PIXEL_COUNT"] - u2["NODATA_PIXEL_COUNT"]) * 4 / 10000)

Layar menampilkan Luas KPH Alpha (ha): 4.0.

Skrip 16.5: Mengelompokkan lereng menjadi tiga kelas #

Peta kemiringan berisi puluhan ribu angka berbeda. Atasan lebih mudah membaca tiga kelas: datar, landai, curam. Berkasnya b3_05_reklas.py; di dalamnya tertulis "Skrip 5.1".

Tujuan. Mengubah kemiringan lereng menjadi tiga kelas, lalu menghitung jumlah piksel dan luas tiap kelas.

Bahan. Lereng_derajat.tif di folder OUT, hasil Skrip 16.3. Jalankan Skrip 16.3 lebih dulu. Skrip ini tidak memakai DATA.

Istilah baru di skrip ini:

• reklasifikasi: Mengganti rentang nilai menjadi nomor kelas.
• tabel kelas: Batas bawah, batas atas, nomor baru.
• tabel hasil: Layer berisi baris, bukan peta.
• Byte: Bilangan bulat 0 sampai 255.

Skrip utuh:

# Skrip 5.1: reklasifikasi kemiringan lereng menjadi tiga kelas dan hitung luasnya. Penulis: Badar Mubarok Yogaswara
# Jalankan Skrip 3.1 lebih dulu (butuh Lereng_derajat.tif).
import os
import processing
from qgis.core import QgsVectorLayer

OUT = os.environ.get("B3_OUT", r"C:/kerja_b3")
os.makedirs(OUT, exist_ok=True)

# Tabel: batas bawah, batas atas, kelas baru. Aturan batas: bawah <= nilai < atas (RANGE_BOUNDARIES = 1)
tabel = [0, 5, 1,       # kelas 1 = datar   (0 sampai kurang dari 5 derajat)
         5, 15, 2,      # kelas 2 = landai  (5 sampai kurang dari 15 derajat)
         15, 90, 3]     # kelas 3 = curam   (15 derajat ke atas)
r = processing.run("native:reclassifybytable", {
    "INPUT_RASTER": OUT + "/Lereng_derajat.tif", "RASTER_BAND": 1,
    "TABLE": tabel, "RANGE_BOUNDARIES": 1, "NODATA_FOR_MISSING": True,
    "NO_DATA": 255, "DATA_TYPE": 0,          # 0 = Byte
    "OUTPUT": OUT + "/Lereng_kelas.tif"})

u = processing.run("native:rasterlayeruniquevaluesreport", {"INPUT": r["OUTPUT"], "BAND": 1,
                                                            "OUTPUT_TABLE": "TEMPORARY_OUTPUT"})
nama = {1: "Datar", 2: "Landai", 3: "Curam"}
print("Kelas | jumlah piksel | luas (ha)")
for f in u["OUTPUT_TABLE"].getFeatures():
    k = int(f["value"])
    print(k, nama.get(k, "NoData"), f["count"], round(f["count"] * 4 / 10000, 2))

Penjelasan baris demi baris.

Hasil. Layar menampilkan (uji pada QGIS 4.0.2):

Kelas | jumlah piksel | luas (ha)
1 Datar 19018 7.61
2 Landai 19902 7.96
3 Curam 1080 0.43

Cek: 19018 + 19902 + 1080 = 40000, seluruh piksel citra. Luasnya 7,61 + 7,96 + 0,43 = 16 hektare.

Bagian yang boleh diubah.

• Baris 8, lokasi. Atur B3_OUT atau ganti alamat cadangan di baris 7. Harus sama dengan OUT di Skrip 16.3.
• Baris 11-13, batas kelas. Ubah angkanya untuk mengubah luas tiap kelas. Dari uji: dengan batas 8 dan 20 derajat, hasilnya Datar 23930 piksel (9,57 ha), Landai 15981 (6,39 ha), dan Curam 89 (0,04 ha).
• Baris 16, RANGE_BOUNDARIES. Pilihannya: 0 untuk bawah < nilai <= atas, 1 untuk bawah <= nilai < atas, 2 untuk kedua batas ikut, 3 untuk kedua batas tidak ikut. Dari uji, pada data ini keempat pilihan menghasilkan hitungan yang sama, karena nyaris tidak ada piksel yang bernilai tepat di batas.
• Baris 17, NO_DATA dan DATA_TYPE. Angka 255 cocok dengan tipe Byte. Jangan memakai kelas bernomor 255.
• Baris 22, nama kelas. Ganti kata sesuai kebutuhan, dan tambah bila kelas bertambah.

Jebakan.

• Lupa Skrip 16.3. Folder OUT dibuat otomatis, tetapi masih kosong. Berkas lereng tidak ada, dan baris 14 berhenti dengan: _core.QgsProcessingException: Unable to execute algorithm lalu Could not load source layer for INPUT_RASTER: .../Lereng_derajat.tif not found
• Batas atas 90. Nilai lereng di atas 90 derajat jadi NoData. Lereng tidak mungkin melewati 90, jadi aman.
• Angka 4 tertulis tetap. Luas piksel 4 meter persegi, seperti di Skrip 16.4.
• Tabel sementara hilang. TEMPORARY_OUTPUT hanya ada selama sesi.

Latihan variasi.

1. Ubah batas kelas menjadi 8 dan 20 derajat.
2. Tambahkan satu kelas keempat, "Sangat curam", untuk lereng 20 derajat ke atas.

Contoh jawaban. Nomor 1: ganti isi tabel menjadi [0, 8, 1, 8, 20, 2, 20, 90, 3]. Layar menampilkan 9,57 ha, 6,39 ha, dan 0,04 ha. Nomor 2: ganti tabel menjadi [0, 5, 1, 5, 15, 2, 15, 20, 3, 20, 90, 4] dan nama menjadi {1: "Datar", 2: "Landai", 3: "Curam", 4: "Sangat curam"}. Layar menampilkan, untuk dua baris terakhir, 3 Curam 991 0.4 dan 4 Sangat curam 89 0.04. Jumlahnya 1080 piksel, sama dengan kelas curam semula.

Skrip 16.6: Menghitung NDVI dan merata-ratakan per KPH #

Tanaman sehat memantulkan banyak cahaya inframerah dekat dan menyerap cahaya merah. NDVI membandingkan keduanya dalam satu angka. Berkasnya b3_06_ndvi.py; di dalamnya tertulis "Skrip 6.1".

Tujuan. Menghitung NDVI dari band inframerah dekat dan merah, mengelompokkannya menjadi empat kelas, menghitung luas tiap kelas, lalu mencari rata-rata NDVI tiap KPH.

Bahan. Citra_4band.tif di folder OUT, hasil Skrip 16.2, dan Batas_KPH.gpkg di DATA. Jalankan Skrip 16.2 lebih dulu.

Istilah baru di skrip ini:

• NDVI: Indeks vegetasi, hasil membandingkan pantulan inframerah dekat dan merah.
• kalkulator raster: QgsRasterCalculator.
• entri: Satu band yang diberi nama pendek.
• statistik zonal: Statistik raster per poligon.
• extent: Batas luar raster.

Skrip utuh:

# Skrip 6.1: menghitung NDVI, mengelompokkannya, dan merata-ratakan per KPH. Penulis: Badar Mubarok Yogaswara
# Jalankan Skrip 2.1 lebih dulu (butuh Citra_4band.tif).
import os
import processing
from qgis.core import QgsRasterLayer, QgsProject
from qgis.analysis import QgsRasterCalculator, QgsRasterCalculatorEntry

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)

stack = QgsRasterLayer(OUT + "/Citra_4band.tif", "Citra_4band")
nir = QgsRasterCalculatorEntry()
nir.ref, nir.raster, nir.bandNumber = "nir@1", stack, 4      # band 4 = inframerah dekat
red = QgsRasterCalculatorEntry()
red.ref, red.raster, red.bandNumber = "red@1", stack, 3      # band 3 = merah

rumus = '("nir@1" - "red@1") / ("nir@1" + "red@1")'
calc = QgsRasterCalculator(rumus, OUT + "/NDVI.tif", "GTiff", stack.extent(), stack.width(), stack.height(),
                           [nir, red], QgsProject.instance().transformContext())
print("Kode hasil (0 = sukses):", calc.processCalculation())

st = processing.run("native:rasterlayerstatistics", {"INPUT": OUT + "/NDVI.tif", "BAND": 1})
print("NDVI min/maks/rata-rata:", round(st["MIN"], 3), "/", round(st["MAX"], 3), "/", round(st["MEAN"], 3))

# Kelompokkan: sampai 0 = air; 0-0,3 = lahan terbuka; 0,3-0,6 = vegetasi jarang; di atas 0,6 = vegetasi rapat
processing.run("native:reclassifybytable", {
    "INPUT_RASTER": OUT + "/NDVI.tif", "RASTER_BAND": 1,
    "TABLE": [-1, 0, 1,  0, 0.3, 2,  0.3, 0.6, 3,  0.6, 1.01, 4],
    "RANGE_BOUNDARIES": 0, "NODATA_FOR_MISSING": True, "NO_DATA": 255, "DATA_TYPE": 0,
    "OUTPUT": OUT + "/NDVI_kelas.tif"})
u = processing.run("native:rasterlayeruniquevaluesreport", {"INPUT": OUT + "/NDVI_kelas.tif", "BAND": 1,
                                                            "OUTPUT_TABLE": "TEMPORARY_OUTPUT"})
nama = {1: "Air", 2: "Lahan terbuka", 3: "Vegetasi jarang", 4: "Vegetasi rapat"}
for f in u["OUTPUT_TABLE"].getFeatures():
    k = int(f["value"])
    print(k, nama.get(k, "NoData"), f["count"], "piksel =", round(f["count"] * 4 / 10000, 2), "ha")

# Rata-rata NDVI per KPH
z = processing.run("native:zonalstatisticsfb", {
    "INPUT": DATA + "/Batas_KPH.gpkg|layername=Batas_KPH", "INPUT_RASTER": OUT + "/NDVI.tif", "RASTER_BAND": 1,
    "COLUMN_PREFIX": "ndvi_", "STATISTICS": [2], "OUTPUT": "TEMPORARY_OUTPUT"})
for f in z["OUTPUT"].getFeatures():
    print(f["NAMA_KPH"], "NDVI rata-rata:", round(f["ndvi_mean"], 3))

Penjelasan baris demi baris.

Hasil. Layar menampilkan (uji pada QGIS 4.0.2):

Kode hasil (0 = sukses): 0
NDVI min/maks/rata-rata: -0.328 / 0.903 / 0.735
1 Air 1203 piksel = 0.48 ha
2 Lahan terbuka 1767 piksel = 0.71 ha
3 Vegetasi jarang 4994 piksel = 2.0 ha
4 Vegetasi rapat 32036 piksel = 12.81 ha
KPH Alpha NDVI rata-rata: 0.814
KPH Beta NDVI rata-rata: 0.614
KPH Gamma NDVI rata-rata: 0.755

Cek: 1203 + 1767 + 4994 + 32036 = 40000 piksel, dan 0,48 + 0,71 + 2,0 + 12,81 = 16 hektare. Rata-rata tertinggi di KPH Alpha, terendah di KPH Beta. Nama kelas "Air" hanyalah label penulis untuk NDVI nol ke bawah.

Bagian yang boleh diubah.

• Baris 14 dan 16, nomor band. Band inframerah dekat dan merah. Harus cocok dengan urutan di Skrip 16.2.
• Baris 18, rumus. Ganti dengan indeks lain, selama memakai nama entri yang sama.
• Baris 29, batas kelas. Dari uji: mengganti 0.6 menjadi 0.7 menurunkan vegetasi rapat menjadi 30791 piksel (12,32 ha).
• Baris 41-42, awalan dan statistik. Ganti ndvi_ bila perlu, dan ingat kolom hasilnya mengikuti awalan itu (ndvi_mean).

Jebakan.

• Band tertukar membalik hasil tanpa pesan. Dari uji: Skrip 16.2 dengan urutan B08 dan B04 tertukar tetap memberi kode hasil 0. Namun hampir seluruh KPH menjadi "Air" dan rata-rata NDVI bertanda minus (lihat blok di bawah daftar ini). NDVI negatif untuk hutan adalah tanda ada yang keliru.
• Lupa Skrip 16.2. Berkas empat band tidak ada, dan layar menampilkan galat bahwa sumber raster tidak dapat dimuat.
• Pembagi nol. Bila NIR dan merah sama-sama 0 pada satu piksel, rumus membagi nol. Dari uji pada citra buatan (blok 10 x 10 piksel dengan NIR dan merah 0): tidak ada galat, kode hasil tetap 0, dan 100 piksel itu menjadi NoData, jadi tidak ikut terhitung di tabel kelas (jumlahnya 39900, bukan 40000).
• Kode 0 = sukses bergantung versi QGIS. [CEK]
• Pesan "Warning 1". Peringatan GDAL tentang lebar kolom NAMA_KPH bukan galat.

Layar bila band tertukar:

Kode hasil (0 = sukses): 0
NDVI min/maks/rata-rata: -0.903 / 0.328 / -0.735
1 Air 38797 piksel = 15.52 ha
2 Lahan terbuka 1191 piksel = 0.48 ha
3 Vegetasi jarang 12 piksel = 0.0 ha
KPH Alpha NDVI rata-rata: -0.814
KPH Beta NDVI rata-rata: -0.614
KPH Gamma NDVI rata-rata: -0.755

Latihan variasi.

1. Jadikan batas vegetasi rapat 0,7, bukan 0,6.
2. Cetak NDVI rata-rata tiap KPH dengan tiga desimal, tetapi hanya KPH yang di atas 0,7.

Contoh jawaban. Nomor 1: pada tabel di baris 29, ganti 0.3, 0.6, 3, 0.6, 1.01, 4 menjadi 0.3, 0.7, 3, 0.7, 1.01, 4. Layar menampilkan 3 Vegetasi jarang 6239 piksel = 2.5 ha dan 4 Vegetasi rapat 30791 piksel = 12.32 ha. Nomor 2: ubah baris 44 menjadi if f["ndvi_mean"] > 0.7: lalu cetak di bawahnya. Layar menampilkan hanya KPH Alpha (0.814) dan KPH Gamma (0.755).

Skrip 16.7: Satu tabel ringkasan untuk Kepala Seksi #

Semua hasil tadi disatukan dalam satu tabel: tiap KPH dengan luas, NDVI rata-rata, luas lereng curam, dan luas vegetasi rapat. Berkas b3_07_proyek.py; di dalamnya tertulis "Skrip 7.1".

Tujuan. Mencetak tabel ringkasan per KPH dengan memotong raster hasil sebelumnya untuk tiap KPH.

Bahan. Di folder OUT: NDVI.tif dan NDVI_kelas.tif (Skrip 16.6) serta Lereng_kelas.tif (Skrip 16.5). Di DATA: Batas_KPH.gpkg. Jalankan Skrip 16.2, 16.3, 16.5, dan 16.6 lebih dulu.

Istilah baru di skrip ini:

• return: Mengembalikan hasil sebuah fungsi (Bab 7).
• tuple: Daftar tetap, Bab 5.
• %s: Menyisipkan teks, Bab 9.
• docstring: Catatan di bawah def.

Skrip utuh:

# Skrip 7.1: ringkasan per KPH (luas, NDVI rata-rata, luas lereng curam, luas vegetasi rapat). Penulis: Badar Mubarok Yogaswara
# Jalankan Skrip 2.1, 3.1, 5.1, dan 6.1 lebih dulu (butuh NDVI.tif, NDVI_kelas.tif, Lereng_kelas.tif di OUT).
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")
BATAS = DATA + "/Batas_KPH.gpkg|layername=Batas_KPH"


def luas_kelas(raster, mask, kelas):
    """Luas (ha) piksel bernilai 'kelas' di dalam mask."""
    pot = processing.run("gdal:cliprasterbymasklayer", {
        "INPUT": raster, "MASK": mask, "NODATA": 255, "CROP_TO_CUTLINE": True,
        "KEEP_RESOLUTION": True, "OUTPUT": "TEMPORARY_OUTPUT"})
    u = processing.run("native:rasterlayeruniquevaluesreport", {
        "INPUT": pot["OUTPUT"], "BAND": 1, "OUTPUT_TABLE": "TEMPORARY_OUTPUT"})
    for f in u["OUTPUT_TABLE"].getFeatures():
        if int(f["value"]) == kelas:
            return f["count"] * 4 / 10000
    return 0.0


print("KPH | luas (ha) | NDVI rata-rata | curam (ha) | vegetasi rapat (ha)")
for nama in ("KPH Alpha", "KPH Beta", "KPH Gamma"):
    lyr = QgsVectorLayer(BATAS, "batas", "ogr")
    lyr.setSubsetString("\"NAMA_KPH\" = '%s'" % nama)
    luas = sum(f.geometry().area() for f in lyr.getFeatures()) / 10000
    z = processing.run("native:zonalstatisticsfb", {
        "INPUT": lyr, "INPUT_RASTER": OUT + "/NDVI.tif", "RASTER_BAND": 1,
        "COLUMN_PREFIX": "n_", "STATISTICS": [2], "OUTPUT": "TEMPORARY_OUTPUT"})
    ndvi = [f["n_mean"] for f in z["OUTPUT"].getFeatures()][0]
    curam = luas_kelas(OUT + "/Lereng_kelas.tif", lyr, 3)
    rapat = luas_kelas(OUT + "/NDVI_kelas.tif", lyr, 4)
    print(nama, "|", round(luas, 2), "|", round(ndvi, 3), "|", round(curam, 2), "|", round(rapat, 2))

Penjelasan baris demi baris.

Hasil. Layar menampilkan (uji pada QGIS 4.0.2):

KPH | luas (ha) | NDVI rata-rata | curam (ha) | vegetasi rapat (ha)
KPH Alpha | 4.0 | 0.814 | 0.0 | 3.75
KPH Beta | 4.0 | 0.614 | 0.0 | 2.0
KPH Gamma | 8.0 | 0.755 | 0.43 | 7.06

Pada sebagian komputer, baris Warning 1: Field NAMA_KPH ... ikut tampil beberapa kali. Itu peringatan GDAL, bukan galat. Cek silang dengan skrip sebelumnya: lereng curam 0,43 hektare sama dengan hasil Skrip 16.5, dan vegetasi rapat 3,75 + 2,0 + 7,06 = 12,81 hektare sama dengan hasil Skrip 16.6. Seluruh lereng curam berada di KPH Gamma.

Bagian yang boleh diubah.

• Baris 7-8, lokasi. Atur B3_DATA dan B3_OUT, atau ganti alamat cadangannya. OUT harus sama dengan yang dipakai Skrip 16.3 sampai 16.6.
• Baris 26, daftar nama KPH. Harus sama persis dengan isi kolom NAMA_KPH. Tambah nama bila ada KPH lain.
• Baris 34-35, nomor kelas. Ganti 3 dan 4 sesuai kelas yang ingin dilihat. Dari uji: mengganti 3 menjadi 2 (landai) memberi 1,67, 0,29, dan 6,0 hektare, yang jumlahnya 7,96, sama dengan Skrip 16.5.
• Baris 21, faktor 4. Luas piksel dalam meter persegi.
• Baris 15, NODATA: 255. Harus sama dengan NO_DATA pada skrip yang membuat rasternya. Skrip 16.4 memakai -9999 karena itu DEM dengan nilai tinggi, sedangkan kelas bertipe Byte memakai 255.

Jebakan.

• Nama KPH tertulis tetap. Bila nama tidak ada di data, layer batas kosong. Dari uji dengan "KPH Delta", skrip berhenti di baris 33 dengan IndexError: list index out of range, karena daftar hasil kosong.
• Skrip pendahulu belum dijalankan, atau OUT berbeda. Dari uji: galat Could not load source layer for INPUT_RASTER: .../NDVI.tif not found.
• Faktor 4 tertulis tetap. Sama seperti Skrip 16.4 dan 16.5.
• Hasil luas_kelas untuk kelas yang tidak ada adalah 0.0. Itu benar, bukan galat. Seperti KPH Alpha dan Beta untuk kelas curam.
• Peringatan lebar kolom. Bukan masalah.

Latihan variasi.

1. Tampilkan luas lereng landai (kelas 2), bukan curam.
2. Ganti daftar KPH menjadi hanya "KPH Gamma".

Contoh jawaban. Nomor 1: ganti 3 menjadi 2 di baris 34, dan kata "curam" di baris 25 menjadi "landai". Layar menampilkan KPH Alpha | 4.0 | 0.814 | 1.67 | 3.75, KPH Beta | 4.0 | 0.614 | 0.29 | 2.0, dan KPH Gamma | 8.0 | 0.755 | 6.0 | 7.06. Nomor 2: tulis ("KPH Gamma",) di baris 26, dengan koma setelah nama. Layar menampilkan satu baris data: KPH Gamma | 8.0 | 0.755 | 0.43 | 7.06. Tanpa koma itu, Python membaca tanda kurung sebagai pengelompokan biasa dan memecah tulisan menjadi huruf. [kemungkinan]

Cek paham #

1. Mengapa urutan berkas di baris 11 Skrip 16.2 penting?
2. Apa yang terjadi pada Skrip 16.1 bila jalur DATA salah, dan baris mana yang menanganinya?
3. Dari mana angka 4 pada rumus luas berasal, dan kapan ia menjadi salah?
4. Skrip mana yang harus dijalankan lebih dulu sebelum Skrip 16.7?

Jawaban.

1. Urutan berkas menentukan nomor band: 1 biru, 2 hijau, 3 merah, 4 inframerah dekat. Skrip 16.6 memakai band 3 dan 4. Bila tertukar, NDVI terbalik tanda.
2. Skrip berhenti dengan pesan Raster tidak terbaca. Periksa DATA: .... Itu hasil pengaman di baris 10-11 (raise SystemExit).
3. Dari ukuran piksel 2 m x 2 m, yaitu 4 meter persegi per piksel. Bila piksel data Anda berukuran lain, angka itu harus diganti.
4. Skrip 16.2, 16.3, 16.5, dan 16.6, karena Skrip 16.7 memakai NDVI.tif, NDVI_kelas.tif, dan Lereng_kelas.tif.

Kesalahan umum #

• Membiarkan alamat cadangan DATA dan OUT. Keduanya hanya alamat contoh, bukan alamat di komputer Anda. Perbaikan: atur B3_DATA dan B3_OUT, atau ganti alamat cadangan di tiap skrip dengan isi yang sama.
• Mengabaikan pesan pengaman. Skrip 16.1 dan 16.2 berhenti dengan pesan bila citra tidak terbaca. Perbaikan: baca pesannya, lalu periksa DATA dan OUT.
• Memakai OUT yang berbeda antarskrip. Skrip 16.5, 16.6, dan 16.7 mencari berkas di OUT. Perbaikan: atur B3_OUT sekali, supaya semua skrip memakai folder yang sama.
• Menukar urutan band. Hasil salah tanpa pesan. Perbaikan: cocokkan rata-rata tiap band dengan berkas asalnya.
• Memakai angka 4 untuk data lain. Luas salah tanpa pesan. Perbaikan: ambil ukuran piksel dari Skrip 16.1 dan hitung ulang.
• Menjalankan skrip tidak berurutan. Skrip 16.5, 16.6, dan 16.7 membaca berkas buatan skrip lain. Perbaikan: ikuti peta ketergantungan di awal bab.

Ringkasan dan latihan #

• Raster adalah kisi piksel. Satu lembar angka disebut band.
• QgsRasterLayer membuka raster, dan processing.run menjalankan alat. Hasilnya berupa kamus.
• gdal:merge menumpuk band, native:slope dan kawan-kawannya menurunkan DEM, gdal:cliprasterbymasklayer memotong.
• Reklasifikasi mengubah rentang nilai menjadi kelas, dan laporan nilai unik menghitung piksel tiap kelas.
• NDVI = (inframerah dekat - merah) / (inframerah dekat + merah), dihitung dengan QgsRasterCalculator.
• Pengaman di awal skrip: os.makedirs membuat folder hasil, dan raise SystemExit menghentikan skrip bila citra tidak terbentuk atau tidak terbaca.
• Tiga hal yang diganti untuk data baru: lokasi data dan folder hasil, nama berkas, dan nomor band atau batas kelas.

Latihan. Pada Skrip 16.1, periksa Citra_B03.tif (band hijau). Lalu pada Skrip 16.6, ubah batas vegetasi rapat menjadi 0,5 dan catat luasnya. Berapa hektare vegetasi rapat sekarang dibandingkan 12,81?

Tabel perbandingan: pekerjaan raster di menu dan di skrip #