BAB 18: Dari Peta ke Jawaban: Buffer, Overlay, Join, dan Model dengan Skrip

Studi kasus: "Berapa hektare yang benar-benar bisa dipanen?" #

Kepala Seksi bertanya: "Dari 16 hektare petak kita, berapa hektare yang berada di luar jarak aman sungai dan permukiman? Rinci per KPH dan per jenis tanaman. Berapa titik panas yang jatuh di areal itu?" Jawabannya tidak terbaca di layar. Anda harus menghitungnya.

Alat penghitungnya sudah ada di QGIS. Skrip menyatukan alat-alat itu menjadi satu resep yang bisa diulang, misalnya bila Kepala Seksi minta jarak aman diganti 50 meter.

Konsep: geoprocessing dalam tiga kalimat #

Geoprocessing adalah mengolah layer lewat sebuah alat untuk menghasilkan layer baru. Bayangkan mesin penggiling: bahan masuk, mesin bekerja menurut pengaturan, hasil giling keluar, dan bahan aslinya tidak berubah. Di skrip, satu perintah, processing.run, memanggil mesin itu.

Lima istilah yang perlu Anda kenal sekarang:

• Alat (algorithm): satu mesin dengan tugas tertentu, misalnya native:buffer. Nama alat selalu ditulis dalam tanda kutip.
• Parameter: pengaturan sebuah alat, ditulis dalam kamus (Bab 5), misalnya jarak buffer.
• Layer sementara: hasil yang hanya ada di memori. Di skrip ditulis "memory:". Layer ini hilang saat QGIS ditutup.
• Geometri valid: bentuk objek yang tidak rusak, misalnya poligon yang tidak menyilang dirinya sendiri.
• processing.run: perintah memanggil alat. Hasilnya berupa kamus. Layer hasil diambil dengan ["OUTPUT"].

Sebelum mulai: bahan, jalur, dan urutan #

Bahan. Siapkan dua folder data.

• paket-i1: Petak.gpkg (16 petak, kolom KODE_PETAK, JENIS, TAHUN_TANAM), Sungai.gpkg (2 garis: Sungai Tengah berorde 1 dan Anak Sungai Barat berorde 2), Permukiman.gpkg (2 poligon), Jalan.gpkg (2 garis), dan Perhutanan_Sosial.gpkg (1 poligon).
• paket-data: Batas_KPH.gpkg (3 KPH: Alpha, Beta, Gamma) dan Hotspot_Nasional.gpkg (12 titik panas).

Semua layer memakai sistem koordinat EPSG:32749 (UTM zona 49 Selatan, satuan meter). Kode petak di seri ini berbentuk P07, tanpa tanda hubung.

Jalur. Skrip memakai tiga lokasi di komputer penulis: PAKET = r"C:\latihan\paket-i1", DATA = r"C:\latihan\paket-data", dan (di Skrip 18.6 dan 18.7) HASIL = r"C:\latihan\hasil". Ganti dengan lokasi folder Anda.

Nama berkas disambung dengan garis miring terbalik, misalnya PAKET + r"\Petak.gpkg". Itu benar di Windows. Di macOS dan Ubuntu pemisah folder adalah garis miring biasa, jadi Anda perlu mengganti r"\Petak.gpkg" menjadi "/Petak.gpkg" di semua baris pemuatan. Bagian macOS dan Ubuntu ini belum saya uji. [CEK]

Urutan. Skrip 18.1 sampai 18.6 berdiri sendiri, karena masing-masing memuat datanya sendiri. Skrip 18.7 membaca berkas model yang ditulis Skrip 18.6, jadi jalankan 18.6 lebih dulu.

Menjalankan. Buka Konsol Python QGIS atau Editor Skrip (Bab 3), tempel isi skrip, lalu jalankan.

Pengujian. Semua hasil di bab ini dijalankan di QGIS 4.0.2 tanpa jendela utama. Saya menjalankan salinan skrip dengan tiga lokasi di atas diganti ke lokasi salinan data. Kode yang ditampilkan di bab ini tetap isi berkas asli.

Skrip 18.1: Memeriksa data sebelum dianalisis #

Berkas 1_1_cek_data.py. Di dalam berkasnya tidak tertulis nomor skrip, hanya nama berkas. Skrip ini punya dua tugas: memeriksa data, lalu memperkenalkan "blok pembuka" yang muncul di enam skrip bab ini.

Tujuan. Memuat empat layer, memeriksa sistem koordinat dan jumlah objeknya, memeriksa kevalidan geometri petak, lalu mengukur panjang sungai yang masuk KPH Gamma.

Bahan. paket-i1 (Petak, Sungai, Permukiman) dan paket-data (Batas_KPH). Skrip berdiri sendiri.

Sebelum menghitung apa pun, periksa dulu bahannya. Itu seperti memeriksa timbangan sebelum menimbang. Bila data rusak, semua angka sesudahnya menyesatkan.

Skrip utuh:

# -*- coding: utf-8 -*-
# Penulis: Badar Mubarok Yogaswara. Dijalankan di Python Console QGIS. Data sintetis, bukan kondisi lapangan nyata.
from qgis.core import QgsVectorLayer, QgsProject, QgsProperty
import processing

PAKET = r"C:\latihan\paket-i1"      # ubah sesuai letak folder di komputer Anda
DATA = r"C:\latihan\paket-data"


def muat(berkas, lapisan, nama):
    """Memuat satu layer GeoPackage dan menambahkannya ke proyek."""
    lyr = QgsVectorLayer("%s|layername=%s" % (berkas, lapisan), nama, "ogr")
    if not lyr.isValid():
        raise RuntimeError("Layer tidak terbaca: " + berkas)
    QgsProject.instance().addMapLayer(lyr)
    return lyr


def luas_ha(lyr):
    """Jumlah luas semua poligon dalam hektare (CRS harus UTM, satuan meter)."""
    return sum(f.geometry().area() for f in lyr.getFeatures()) / 10000.0


def buf(lyr, jarak):
    """Buffer dengan dissolve, hasil sementara."""
    return processing.run("native:buffer", {"INPUT": lyr, "DISTANCE": jarak, "DISSOLVE": True, "OUTPUT": "memory:"})["OUTPUT"]


petak = muat(PAKET + r"\Petak.gpkg", "Petak", "Petak")
sungai = muat(PAKET + r"\Sungai.gpkg", "Sungai", "Sungai")
kampung = muat(PAKET + r"\Permukiman.gpkg", "Permukiman", "Permukiman")
kph = muat(DATA + r"\Batas_KPH.gpkg", "Batas_KPH", "Batas_KPH")

# 1. Periksa CRS dan jumlah objek tiap layer
for lyr in (petak, sungai, kampung, kph):
    print(lyr.name(), lyr.crs().authid(), lyr.featureCount(), "objek")

# 2. Periksa kevalidan geometri (laporan, bukan perbaikan)
hasil = processing.run("native:checkvalidity", {"INPUT_LAYER": petak, "METHOD": 2, "VALID_OUTPUT": "memory:", "INVALID_OUTPUT": "memory:", "ERROR_OUTPUT": "memory:"})
print("Petak valid:", hasil["VALID_COUNT"], "tidak valid:", hasil["INVALID_COUNT"])

# 3. Potong sungai dengan KPH Gamma saja (Clip), hasil sementara
gamma = processing.run("native:extractbyexpression", {"INPUT": kph, "EXPRESSION": "\"NAMA_KPH\" = 'KPH Gamma'", "OUTPUT": "memory:"})["OUTPUT"]
potong = processing.run("native:clip", {"INPUT": sungai, "OVERLAY": gamma, "OUTPUT": "memory:"})["OUTPUT"]
panjang = sum(f.geometry().length() for f in potong.getFeatures())
print("Sungai di KPH Gamma: %d ruas, %.1f m" % (potong.featureCount(), panjang))

Blok pembuka: baris 1-28 #

Baris 1-28 sama persis di enam skrip: 1_1, 2_1, 3_1, 4_1, 5_1, dan 7_1. Penulis menyalinnya ke tiap berkas supaya setiap skrip bisa dijalankan sendiri. Penjelasannya hanya di sini. Di skrip lain, Anda cukup membaca "baris 1-28: blok pembuka".

Gambar di atas hasil membaca isi skrip. Fungsi muat dipakai semua. Fungsi luas_ha dan buf hanya dipakai sebagian, dan QgsProperty hanya di Skrip 18.2. Kelebihan itu tidak berbahaya, hanya mengisi baris.

Isi khusus Skrip 18.1: baris 29-46 #

Hasil. Layar menampilkan (berdasarkan uji pada QGIS 4.0.2):

Petak EPSG:32749 16 objek
Sungai EPSG:32749 2 objek
Permukiman EPSG:32749 2 objek
Batas_KPH EPSG:32749 3 objek
Petak valid: 16 tidak valid: 0
Sungai di KPH Gamma: 2 ruas, 570.3 m

Dari hasil ini Anda tahu: keempat layer memakai EPSG:32749 (meter), 16 petak semuanya valid, dan sungai yang masuk KPH Gamma sepanjang sekitar 570 meter. Empat layer baru juga muncul di panel Layers.

Bagian yang boleh diubah.

• Baris 6-7, lokasi data. Ganti dengan folder Anda.
• Baris 29-32, nama berkas dan layer. Ganti untuk memeriksa layer lain, misalnya Jalan.gpkg dengan layer Jalan.
• Baris 43, nama KPH. Ganti 'KPH Gamma' dengan 'KPH Beta' atau 'KPH Alpha'.
• Baris 39, METHOD. Angka ini memilih cara pemeriksaan: 0 untuk pengaturan digitasi, 1 untuk QGIS, 2 untuk GEOS.

Jebakan.

• Jalur salah. Bila PAKET menunjuk folder yang tidak ada, skrip berhenti di baris 29 dengan RuntimeError: Layer tidak terbaca: C:\latihan\paket-i1\Petak.gpkg. Pesan ini jelas dan menyebut berkas yang bermasalah. Fungsi muat memang dibuat untuk itu.
• Layer ganda. Setiap skrip dijalankan, empat layer ditambahkan lagi ke proyek. Menjalankan dua kali menghasilkan layer Batas_KPH ganda, Permukiman ganda, dan seterusnya (uji saya: delapan layer setelah dua kali jalan). Hapus layer lama atau mulai dari proyek baru.
• Fungsi yang tidak dipakai. luas_ha dan buf ada di blok pembuka, tetapi skrip ini tidak memanggilnya. Itu bukan galat. Baca skrip dengan kritis: yang ada belum tentu dipakai.
• Ukuran dalam meter. Skrip menganggap satuan peta meter. Bila layer Anda berkoordinat derajat, panjang dan luas menjadi tidak berarti. Periksa baris EPSG: pada hasil.

Latihan variasi. Ubah skrip supaya mengukur sungai di KPH Beta.

Contoh jawaban. Ganti baris 43, lalu ganti juga kata KPH Gamma di dalam tanda kutip pada baris 46 menjadi KPH Beta, karena tulisan itu ditulis tetap:

gamma = processing.run("native:extractbyexpression", {"INPUT": kph, "EXPRESSION": "\"NAMA_KPH\" = 'KPH Beta'", "OUTPUT": "memory:"})["OUTPUT"]

Layar menampilkan Sungai di KPH Beta: 1 ruas, 26.9 m. Bila baris 46 tidak diganti, layar tetap menampilkan Sungai di KPH Gamma: 1 ruas, 26.9 m, padahal angkanya milik KPH Beta. Nama kotak gamma tidak perlu diganti, karena Python tidak peduli artinya. Untuk keterbacaan, sebaiknya Anda menggantinya juga.

Skrip 18.2: Zona aman di sekeliling sungai #

Berkas 2_1_buffer_sempadan.py. Skrip ini hanya berisi hal baru dari baris 29 ke bawah.

Tujuan. Membuat buffer sungai 25 meter dengan dan tanpa penggabungan (dissolve), membuat buffer dengan jarak yang berbeda menurut orde sungai, lalu membuat buffer permukiman 20 meter dan membandingkan luasnya.

Bahan. Sungai.gpkg (kolom ORDE) dan Permukiman.gpkg dari paket-i1.

Buffer adalah zona selebar jarak tertentu di sekeliling sebuah objek. Bayangkan Anda berjalan menyusuri sungai sambil menggambar garis 25 meter di kiri dan kanannya. Dissolve artinya menggabungkan zona-zona yang bertumpuk menjadi satu bidang. Tanpa dissolve, daerah yang tumpang tindih terhitung dua kali.

Skrip utuh:

# -*- coding: utf-8 -*-
# Penulis: Badar Mubarok Yogaswara. Dijalankan di Python Console QGIS. Data sintetis, bukan kondisi lapangan nyata.
from qgis.core import QgsVectorLayer, QgsProject, QgsProperty
import processing

PAKET = r"C:\latihan\paket-i1"      # ubah sesuai letak folder di komputer Anda
DATA = r"C:\latihan\paket-data"


def muat(berkas, lapisan, nama):
    """Memuat satu layer GeoPackage dan menambahkannya ke proyek."""
    lyr = QgsVectorLayer("%s|layername=%s" % (berkas, lapisan), nama, "ogr")
    if not lyr.isValid():
        raise RuntimeError("Layer tidak terbaca: " + berkas)
    QgsProject.instance().addMapLayer(lyr)
    return lyr


def luas_ha(lyr):
    """Jumlah luas semua poligon dalam hektare (CRS harus UTM, satuan meter)."""
    return sum(f.geometry().area() for f in lyr.getFeatures()) / 10000.0


def buf(lyr, jarak):
    """Buffer dengan dissolve, hasil sementara."""
    return processing.run("native:buffer", {"INPUT": lyr, "DISTANCE": jarak, "DISSOLVE": True, "OUTPUT": "memory:"})["OUTPUT"]


sungai = muat(PAKET + r"\Sungai.gpkg", "Sungai", "Sungai")
kampung = muat(PAKET + r"\Permukiman.gpkg", "Permukiman", "Permukiman")

# Buffer sungai 25 m, tanpa dissolve dan dengan dissolve
tanpa = processing.run("native:buffer", {"INPUT": sungai, "DISTANCE": 25, "DISSOLVE": False, "OUTPUT": "memory:"})["OUTPUT"]
gabung = buf(sungai, 25)
print("Tanpa dissolve: %d objek, %.2f ha (tumpang tindih terhitung dua kali)" % (tanpa.featureCount(), luas_ha(tanpa)))
print("Dengan dissolve: %d objek, %.2f ha" % (gabung.featureCount(), luas_ha(gabung)))

# Jarak berbeda menurut orde sungai: 30 m untuk orde 1, 15 m untuk orde 2
rumus = 'CASE WHEN "ORDE" = 1 THEN 30 ELSE 15 END'
variabel = processing.run("native:buffer", {"INPUT": sungai, "DISTANCE": QgsProperty.fromExpression(rumus), "DISSOLVE": True, "OUTPUT": "memory:"})["OUTPUT"]
print("Jarak variabel: %.2f ha" % luas_ha(variabel))

# Buffer permukiman 20 m
print("Buffer permukiman 20 m: %.2f ha" % luas_ha(buf(kampung, 20)))

Penjelasan baris demi baris. Baris 1-28 adalah blok pembuka, sama dengan Skrip 18.1. Hanya muat, luas_ha, buf, dan QgsProperty yang dipakai di sini.

Hasil. Layar menampilkan (berdasarkan uji pada QGIS 4.0.2):

Tanpa dissolve: 2 objek, 3.37 ha (tumpang tindih terhitung dua kali)
Dengan dissolve: 1 objek, 3.12 ha
Jarak variabel: 3.28 ha
Buffer permukiman 20 m: 1.18 ha

Selisih 3,37 dan 3,12, yaitu 0,25 hektare, adalah daerah tumpang tindih yang terhitung dua kali bila tanpa dissolve. Jarak variabel menghasilkan 3,28 hektare, lebih lebar dari 25 meter karena sungai berorde 1 mendapat 30 meter.

Bagian yang boleh diubah.

• Baris 33-34, jarak 25. Ganti di kedua baris sekaligus. Bila hanya satu, perbandingan dissolve menjadi tidak adil.
• Baris 39, rumus. Ganti angka 30 dan 15 untuk mengubah lebar tiap orde. Angka itu hanya contoh latihan, bukan ketentuan resmi sempadan sungai.
• Baris 44, jarak 20. Ganti jarak buffer permukiman.
• Baris 33, "DISSOLVE": False. Ganti True untuk menggabungkan hasil.

Jebakan.

• Kolom salah ketik di dalam rumus. Bila Anda menulis "ORDEE", tidak ada galat. QGIS menganggap kolom itu kosong, sehingga semua sungai jatuh ke jarak 15 meter. Uji saya: Jarak variabel: 1.22 ha, jauh lebih kecil dari 3,28. Periksa angka hasil, bukan hanya ada atau tidaknya galat.
• Satuan harus meter. Fungsi luas_ha membaca luas dalam satuan peta. Pada layer berkoordinat derajat, angkanya tidak berarti.
• Dua tempat untuk satu angka. Jarak 25 ada di baris 33 dan 34. Mengubah satu baris saja membuat tampilan "tanpa" dan "dengan" dissolve tidak setara.
• Tumpang tindih. Tanpa dissolve, luas terlalu besar (3,37 lawan 3,12). Selalu gabungkan sebelum menghitung luas total sempadan.

Latihan variasi. Lebarkan sempadan sungai menjadi 50 meter, baik tanpa maupun dengan dissolve.

Contoh jawaban. Ganti angka 25 pada baris 33 dan 34:

tanpa = processing.run("native:buffer", {"INPUT": sungai, "DISTANCE": 50, "DISSOLVE": False, "OUTPUT": "memory:"})["OUTPUT"]
gabung = buf(sungai, 50)

Layar menampilkan Tanpa dissolve: 2 objek, 7.50 ha ... dan Dengan dissolve: 1 objek, 6.52 ha. Selisihnya 0,98 hektare, jauh lebih besar daripada selisih pada 25 meter (0,25 hektare). Itu wajar: sempadan yang lebih lebar saling bertumpuk lebih banyak.

Skrip 18.3: Dua peta bertumpuk, tiga cara membaca tumpukannya #

Berkas 3_1_overlay.py. Seperti Skrip 18.2, baris 1-28 adalah blok pembuka.

Tujuan. Memperlihatkan tiga operasi tumpang susun: Intersect, Difference, dan Union. Skrip menghitung luas bagian petak yang masuk sempadan, luas areal efektif, serta luas blok perhutanan sosial di dalam dan di luar petak.

Bahan. Petak.gpkg, Sungai.gpkg, Permukiman.gpkg, dan Perhutanan_Sosial.gpkg dari paket-i1.

Bayangkan dua lembar plastik transparan yang ditumpuk: lembar petak dan lembar sempadan. Ada tiga cara menggunting hasilnya, seperti pada gambar.

• Intersect: ambil hanya bagian yang bertumpuk.
• Difference: ambil bagian lembar pertama yang tidak tertutup lembar kedua. Urutan masukan menentukan hasil.
• Union: ambil seluruh bagian kedua lembar, dengan garis potong di tempat mereka bertumpuk.

Skrip utuh:

# -*- coding: utf-8 -*-
# Penulis: Badar Mubarok Yogaswara. Dijalankan di Python Console QGIS. Data sintetis, bukan kondisi lapangan nyata.
from qgis.core import QgsVectorLayer, QgsProject, QgsProperty
import processing

PAKET = r"C:\latihan\paket-i1"      # ubah sesuai letak folder di komputer Anda
DATA = r"C:\latihan\paket-data"


def muat(berkas, lapisan, nama):
    """Memuat satu layer GeoPackage dan menambahkannya ke proyek."""
    lyr = QgsVectorLayer("%s|layername=%s" % (berkas, lapisan), nama, "ogr")
    if not lyr.isValid():
        raise RuntimeError("Layer tidak terbaca: " + berkas)
    QgsProject.instance().addMapLayer(lyr)
    return lyr


def luas_ha(lyr):
    """Jumlah luas semua poligon dalam hektare (CRS harus UTM, satuan meter)."""
    return sum(f.geometry().area() for f in lyr.getFeatures()) / 10000.0


def buf(lyr, jarak):
    """Buffer dengan dissolve, hasil sementara."""
    return processing.run("native:buffer", {"INPUT": lyr, "DISTANCE": jarak, "DISSOLVE": True, "OUTPUT": "memory:"})["OUTPUT"]


petak = muat(PAKET + r"\Petak.gpkg", "Petak", "Petak")
sungai = muat(PAKET + r"\Sungai.gpkg", "Sungai", "Sungai")
kampung = muat(PAKET + r"\Permukiman.gpkg", "Permukiman", "Permukiman")
ps = muat(PAKET + r"\Perhutanan_Sosial.gpkg", "Perhutanan_Sosial", "Perhutanan_Sosial")

sempadan = buf(sungai, 25)
zona_kampung = buf(kampung, 20)

# Intersect: bagian petak yang masuk sempadan sungai
irisan = processing.run("native:intersection", {"INPUT": petak, "OVERLAY": sempadan, "OUTPUT": "memory:"})["OUTPUT"]
print("Intersect: %d potongan, %.2f ha" % (irisan.featureCount(), luas_ha(irisan)))

# Difference: petak dikurangi sempadan, lalu dikurangi zona permukiman
sisa = processing.run("native:difference", {"INPUT": petak, "OVERLAY": sempadan, "OUTPUT": "memory:"})["OUTPUT"]
efektif = processing.run("native:difference", {"INPUT": sisa, "OVERLAY": zona_kampung, "OUTPUT": "memory:"})["OUTPUT"]
print("Difference sungai: %.2f ha; ditambah permukiman: %.2f ha" % (luas_ha(sisa), luas_ha(efektif)))

# Union: petak digabung dengan blok perhutanan sosial
gab = processing.run("native:union", {"INPUT": petak, "OVERLAY": ps, "OUTPUT": "memory:"})["OUTPUT"]
luar = [f for f in gab.getFeatures() if f["KODE_PETAK"] is None]
print("Union: %d objek, %d di luar petak (%.2f ha)" % (gab.featureCount(), len(luar), sum(f.geometry().area() for f in luar) / 10000.0))
irisan_ps = processing.run("native:intersection", {"INPUT": petak, "OVERLAY": ps, "OUTPUT": "memory:"})["OUTPUT"]
print("Petak dalam perhutanan sosial: %d potongan, %.2f ha dari %.2f ha blok" % (irisan_ps.featureCount(), luas_ha(irisan_ps), luas_ha(ps)))

Penjelasan baris demi baris. Baris 1-28 adalah blok pembuka. Fungsi buf dan luas_ha dipakai di sini.

Hasil. Layar menampilkan (berdasarkan uji pada QGIS 4.0.2):

Intersect: 8 potongan, 2.83 ha
Difference sungai: 13.17 ha; ditambah permukiman: 12.39 ha
Union: 23 objek, 1 di luar petak (0.29 ha)
Petak dalam perhutanan sosial: 6 potongan, 3.47 ha dari 3.77 ha blok

Cara membaca. Dari 16 hektare petak, sempadan sungai mengambil 2,83 hektare, dan sisanya 13,17 hektare. Zona permukiman mengambil 0,78 hektare lagi, sehingga areal efektif 12,39 hektare. Cek angka: 16,00 dikurangi 2,83 sama dengan 13,17, dan 13,17 dikurangi 12,39 sama dengan 0,78.

Pada blok perhutanan sosial, 3,47 hektare ada di dalam petak dan 0,29 hektare di luar. Jumlahnya, 3,76, sama dengan 3,77 setelah pembulatan.

Bagian yang boleh diubah.

• Baris 34-35, jarak buffer. Ganti 25 dan 20.
• Baris 29-32, layer. Ganti dengan layer lain, tetapi kolom KODE_PETAK di baris 48 harus ada di layer petak Anda.
• Baris 38, 42, 43, 47, 50, urutan masukan. Tukar INPUT dan OVERLAY untuk melihat akibatnya pada Difference.

Jebakan.

• Urutan masukan Difference. Difference tidak simetris. Saya menukar masukan baris 42 (sempadan sebagai INPUT, petak sebagai OVERLAY). Hasilnya Difference sungai: 0.29 ha, yaitu bagian sempadan yang berada di luar petak. Untuk areal efektif, petaklah yang harus jadi INPUT.
• Penanda "di luar petak" memakai kolom kosong. Baris 48 menganggap KODE_PETAK kosong berarti "bukan dari petak". Bila data petak Anda sendiri punya kode yang kosong, hitungannya salah.
• Jumlah objek berbeda antar alat. Intersect menghasilkan 8 potongan. Union menghasilkan 23 objek. Jumlah objek tidak sama dengan jumlah petak. Jangan menghitung petak dari jumlah potongan.
• Satuan meter. Sama dengan Skrip 18.2: luas_ha hanya benar untuk layer berkoordinat meter.

Latihan variasi. Hitung luas bagian petak yang masuk zona permukiman (zona_kampung) dengan Intersect.

Contoh jawaban. Tambahkan dua baris di akhir skrip:

irisan_kampung = processing.run("native:intersection", {"INPUT": petak, "OVERLAY": zona_kampung, "OUTPUT": "memory:"})["OUTPUT"]
print("Petak dalam zona permukiman: %d potongan, %.2f ha" % (irisan_kampung.featureCount(), luas_ha(irisan_kampung)))

Layar menampilkan Petak dalam zona permukiman: 5 potongan, 1.18 ha. Angka 1,18 sama dengan luas buffer permukiman di Skrip 18.2. Artinya seluruh zona itu berada di dalam petak.

Skrip 18.4: Melebur batas dan merekap luas #

Berkas 4_1_dissolve_rekap.py. Baris 1-28 adalah blok pembuka.

Tujuan. Dari areal efektif (petak dikurangi sempadan sungai dan zona permukiman), membuat rekap luas per jenis tanaman dengan dissolve, memecah bidang gabungan menjadi bagian-bagian, lalu membandingkannya dengan statistik berkelompok.

Bahan. Petak.gpkg (kolom JENIS, TAHUN_TANAM), Sungai.gpkg, dan Permukiman.gpkg dari paket-i1.

Dissolve melebur batas antar objek yang punya isi kolom sama, seperti menghapus tembok antar kamar yang sama fungsinya. Petak berjenis Jati yang berjauhan lalu menjadi satu objek. Objek hasil seperti itu disebut multipart (satu objek terdiri dari banyak bagian terpisah). Lawannya singlepart: satu objek, satu bagian.

Skrip utuh:

# -*- coding: utf-8 -*-
# Penulis: Badar Mubarok Yogaswara. Dijalankan di Python Console QGIS. Data sintetis, bukan kondisi lapangan nyata.
from qgis.core import QgsVectorLayer, QgsProject, QgsProperty
import processing

PAKET = r"C:\latihan\paket-i1"      # ubah sesuai letak folder di komputer Anda
DATA = r"C:\latihan\paket-data"


def muat(berkas, lapisan, nama):
    """Memuat satu layer GeoPackage dan menambahkannya ke proyek."""
    lyr = QgsVectorLayer("%s|layername=%s" % (berkas, lapisan), nama, "ogr")
    if not lyr.isValid():
        raise RuntimeError("Layer tidak terbaca: " + berkas)
    QgsProject.instance().addMapLayer(lyr)
    return lyr


def luas_ha(lyr):
    """Jumlah luas semua poligon dalam hektare (CRS harus UTM, satuan meter)."""
    return sum(f.geometry().area() for f in lyr.getFeatures()) / 10000.0


def buf(lyr, jarak):
    """Buffer dengan dissolve, hasil sementara."""
    return processing.run("native:buffer", {"INPUT": lyr, "DISTANCE": jarak, "DISSOLVE": True, "OUTPUT": "memory:"})["OUTPUT"]


petak = muat(PAKET + r"\Petak.gpkg", "Petak", "Petak")
sungai = muat(PAKET + r"\Sungai.gpkg", "Sungai", "Sungai")
kampung = muat(PAKET + r"\Permukiman.gpkg", "Permukiman", "Permukiman")

sisa = processing.run("native:difference", {"INPUT": petak, "OVERLAY": buf(sungai, 25), "OUTPUT": "memory:"})["OUTPUT"]
efektif = processing.run("native:difference", {"INPUT": sisa, "OVERLAY": buf(kampung, 20), "OUTPUT": "memory:"})["OUTPUT"]

# Dissolve per jenis tanaman. Nama kolom diperiksa dulu: salah ketik tidak boleh lolos diam-diam
KOLOM = "JENIS"
assert petak.fields().indexOf(KOLOM) >= 0, "Kolom tidak ada di layer Petak: " + KOLOM
per_jenis = processing.run("native:dissolve", {"INPUT": efektif, "FIELD": [KOLOM], "OUTPUT": "memory:"})["OUTPUT"]
awal = processing.run("native:dissolve", {"INPUT": petak, "FIELD": [KOLOM], "OUTPUT": "memory:"})["OUTPUT"]
luas_awal = {f["JENIS"]: f.geometry().area() / 10000.0 for f in awal.getFeatures()}
for f in sorted(per_jenis.getFeatures(), key=lambda f: f["JENIS"]):
    print("%-7s awal %.2f ha, efektif %.2f ha (multipart: %s)" % (f["JENIS"], luas_awal[f["JENIS"]], f.geometry().area() / 10000.0, f.geometry().isMultipart()))

# Dissolve tanpa field: satu poligon untuk seluruh areal efektif
semua = processing.run("native:dissolve", {"INPUT": efektif, "OUTPUT": "memory:"})["OUTPUT"]
print("Dissolve semua: %d objek, %.2f ha" % (semua.featureCount(), luas_ha(semua)))

# Pecah multipart menjadi satu bagian per objek
tunggal = processing.run("native:multiparttosingleparts", {"INPUT": per_jenis, "OUTPUT": "memory:"})["OUTPUT"]
print("Setelah multipart ke singlepart: %d objek" % tunggal.featureCount())

# Atribut lain setelah dissolve: nilai diambil dari satu objek saja (hati-hati)
for f in sorted(per_jenis.getFeatures(), key=lambda f: f["JENIS"]):
    print(f["JENIS"], "TAHUN_TANAM tersisa =", f["TAHUN_TANAM"])

# Rekap luas tanpa dissolve: tambah field luas, lalu Statistics by categories
hitung = processing.run("native:fieldcalculator", {"INPUT": efektif, "FIELD_NAME": "LUAS_HA", "FIELD_TYPE": 0, "FIELD_LENGTH": 10, "FIELD_PRECISION": 4, "FORMULA": "area($geometry) / 10000", "OUTPUT": "memory:"})["OUTPUT"]
statistik = processing.run("qgis:statisticsbycategories", {"INPUT": hitung, "VALUES_FIELD_NAME": "LUAS_HA", "CATEGORIES_FIELD_NAME": ["JENIS"], "OUTPUT": "memory:"})["OUTPUT"]
for f in sorted(statistik.getFeatures(), key=lambda f: f["JENIS"]):
    print("Statistik %-7s %d petak, jumlah luas %.2f ha" % (f["JENIS"], f["count"], f["sum"]))

Penjelasan baris demi baris. Baris 1-28 adalah blok pembuka.

Hasil. Layar menampilkan (berdasarkan uji pada QGIS 4.0.2):

Jati    awal 6.00 ha, efektif 4.65 ha (multipart: True)
Mahoni  awal 4.00 ha, efektif 3.63 ha (multipart: True)
Pinus   awal 3.00 ha, efektif 2.08 ha (multipart: True)
Sengon  awal 3.00 ha, efektif 2.03 ha (multipart: True)
Dissolve semua: 1 objek, 12.39 ha
Setelah multipart ke singlepart: 21 objek
Jati TAHUN_TANAM tersisa = 2008
Mahoni TAHUN_TANAM tersisa = 2015
Pinus TAHUN_TANAM tersisa = 2008
Sengon TAHUN_TANAM tersisa = 2015
Statistik Jati    6 petak, jumlah luas 4.65 ha
Statistik Mahoni  4 petak, jumlah luas 3.63 ha
Statistik Pinus   3 petak, jumlah luas 2.08 ha
Statistik Sengon  3 petak, jumlah luas 2.03 ha

Pada uji tanpa jendela utama, QGIS menambah satu baris peringatan DeprecationWarning sebelum empat baris statistik. Peringatan itu berasal dari alat bawaan QGIS, bukan dari skrip, dan tidak memengaruhi hasil.

Dua cara menghasilkan angka yang sama. Rekap dissolve (4,65, 3,63, 2,08, 2,03) sama dengan jumlah luas dari statistik. Bedanya, dissolve menghasilkan objek geometri baru, sedangkan statistik hanya menghasilkan tabel angka.

Bagian yang boleh diubah.

• Baris 37, KOLOM = "JENIS". Ganti dengan kolom lain, misalnya "TAHUN_TANAM". Nama kolom harus persis sama dengan di tabel atribut. Pengaman baris 38 ikut memeriksanya.
• Baris 33-34, jarak buffer. Ganti 25 dan 20.
• Baris 58, FIELD_PRECISION. Ganti 4 untuk mengubah jumlah angka di belakang koma.
• Baris 59, CATEGORIES_FIELD_NAME. Ganti kolom pengelompok.

Jebakan.

• Pengaman nama kolom (baris 37-38). Alat dissolve tidak protes bila nama kolom salah. Semua objek melebur menjadi satu, dan skrip tetap mencetak angka yang keliru. Baris 38 mencegahnya. Saya mengganti KOLOM menjadi "JENSI". Skrip berhenti dengan AssertionError: Kolom tidak ada di layer Petak: JENSI, sebelum ada angka tercetak. Bila Anda mengganti kolom, ubah KOLOM saja, jangan menulis nama kolom langsung di baris 39-40.
• Pengaman hanya menjaga satu kotak. Baris 38 memeriksa KOLOM pada layer petak. Nama JENIS di baris 41-43, 54, dan 59-61 ditulis langsung dan tidak diperiksa. Setelah mengganti kolom, periksa angka hasilnya.
• Atribut lain setelah dissolve tidak representatif. Setelah dissolve, kolom selain kolom pengelompok mengambil nilai dari satu objek saja. TAHUN_TANAM Jati tertulis 2008, padahal petak Jati ditanam di tahun berbeda (2008, 2012, 2015, 2019). Jangan memakai kolom itu sebagai fakta.
• Multipart sulit dipakai. Hasil dissolve per jenis berupa 4 objek multipart. Bila perlu menghitung bagian, pecah dulu (baris 50).
• Peringatan tanpa dampak. Pesan DeprecationWarning pada alat statistik bukan galat Anda.

Latihan variasi. Rekap luas areal efektif per tahun tanam, bukan per jenis.

Contoh jawaban. Biarkan baris 39-43 seperti semula, karena per_jenis yang dibuat di sana masih dipakai baris 50 dan 54 (bila baris itu dihapus, skrip berhenti dengan NameError: name 'per_jenis' is not defined). Sisipkan tiga baris ini tepat sesudah baris 43:

per_tahun = processing.run("native:dissolve", {"INPUT": efektif, "FIELD": ["TAHUN_TANAM"], "OUTPUT": "memory:"})["OUTPUT"]
for f in sorted(per_tahun.getFeatures(), key=lambda f: f["TAHUN_TANAM"]):
    print(f["TAHUN_TANAM"], "efektif %.2f ha" % (f.geometry().area() / 10000.0))

Empat baris jenis tetap tampil lebih dulu, lalu layar menampilkan 2008 efektif 2.89 ha, 2012 efektif 3.59 ha, 2015 efektif 2.81 ha, dan 2019 efektif 3.09 ha. Jumlahnya 12,38, sama dengan 12,39 setelah pembulatan.

Skrip 18.5: Menghubungkan data dari letaknya #

Berkas 5_1_spatial_join.py. Skrip ini yang terpanjang di bab ini, dan paling padat. Tabel dan hasilnya dikelompokkan per tugas. Baris 1-28 adalah blok pembuka.

Tujuan. Memberi nama KPH pada tiap petak berdasarkan letaknya, menghitung titik panas per petak dan per KPH, memilih petak yang kena titik panas, mencari jalan terdekat dari tiap titik, lalu menghitung titik per jenis tanaman.

Bahan. Petak.gpkg dan Jalan.gpkg dari paket-i1, serta Batas_KPH.gpkg dan Hotspot_Nasional.gpkg dari paket-data.

Spatial join menyalin atribut dari satu layer ke layer lain menurut hubungan letak, bukan menurut kolom yang sama. Bayangkan menempelkan label nama KPH pada tiap petak dengan aturan: "tempel label dari KPH tempat petak ini berada". Aturan hubungan letak itu disebut predikat (misalnya intersect, are within). Ada dua cara menggabung: satu-ke-banyak (satu petak bisa muncul di beberapa baris bila cocok dengan beberapa objek) dan satu-ke-satu (hanya satu baris tiap petak).

Skrip utuh:

# -*- coding: utf-8 -*-
# Penulis: Badar Mubarok Yogaswara. Dijalankan di Python Console QGIS. Data sintetis, bukan kondisi lapangan nyata.
from qgis.core import QgsVectorLayer, QgsProject, QgsProperty
import processing

PAKET = r"C:\latihan\paket-i1"      # ubah sesuai letak folder di komputer Anda
DATA = r"C:\latihan\paket-data"


def muat(berkas, lapisan, nama):
    """Memuat satu layer GeoPackage dan menambahkannya ke proyek."""
    lyr = QgsVectorLayer("%s|layername=%s" % (berkas, lapisan), nama, "ogr")
    if not lyr.isValid():
        raise RuntimeError("Layer tidak terbaca: " + berkas)
    QgsProject.instance().addMapLayer(lyr)
    return lyr


def luas_ha(lyr):
    """Jumlah luas semua poligon dalam hektare (CRS harus UTM, satuan meter)."""
    return sum(f.geometry().area() for f in lyr.getFeatures()) / 10000.0


def buf(lyr, jarak):
    """Buffer dengan dissolve, hasil sementara."""
    return processing.run("native:buffer", {"INPUT": lyr, "DISTANCE": jarak, "DISSOLVE": True, "OUTPUT": "memory:"})["OUTPUT"]


from collections import Counter

petak = muat(PAKET + r"\Petak.gpkg", "Petak", "Petak")
jalan = muat(PAKET + r"\Jalan.gpkg", "Jalan", "Jalan")
kph = muat(DATA + r"\Batas_KPH.gpkg", "Batas_KPH", "Batas_KPH")
titik = muat(DATA + r"\Hotspot_Nasional.gpkg", "Hotspot_Nasional", "Hotspot_Nasional")

# 1. Spatial join petak ke KPH. Predikat 0 = intersect, 5 = are within. Metode 0 = satu-ke-banyak
for nama, pred in (("intersect", 0), ("are within", 5)):
    j = processing.run("native:joinattributesbylocation", {"INPUT": petak, "PREDICATE": [pred], "JOIN": kph, "JOIN_FIELDS": ["NAMA_KPH"], "METHOD": 0, "OUTPUT": "memory:"})["OUTPUT"]
    print("Join dengan %-10s: %d baris (petak ada 16)" % (nama, j.featureCount()))

# 2. Join yang benar: 'are within', ambil satu (metode 1)
petak_kph = processing.run("native:joinattributesbylocation", {"INPUT": petak, "PREDICATE": [5], "JOIN": kph, "JOIN_FIELDS": ["NAMA_KPH"], "METHOD": 1, "OUTPUT": "memory:"})["OUTPUT"]
print("Petak per KPH:", dict(Counter(f["NAMA_KPH"] for f in petak_kph.getFeatures())))

# 3. Hitung titik panas per petak
hitung = processing.run("native:countpointsinpolygon", {"POLYGONS": petak_kph, "POINTS": titik, "FIELD": "JML_TITIK", "OUTPUT": "memory:"})["OUTPUT"]
per_petak = {f["KODE_PETAK"]: int(f["JML_TITIK"]) for f in hitung.getFeatures() if f["JML_TITIK"] > 0}
print("Titik per petak:", per_petak, "jumlah =", sum(per_petak.values()))
per_kph = Counter()
for f in hitung.getFeatures():
    per_kph[f["NAMA_KPH"]] += int(f["JML_TITIK"])
print("Titik per KPH:", dict(per_kph))

# 4. Ringkasan berdasarkan lokasi (join summary): hitung titik per petak
ringkas = processing.run("native:joinbylocationsummary", {"INPUT": petak, "PREDICATE": [0], "JOIN": titik, "JOIN_FIELDS": ["ID_Titik"], "SUMMARIES": [0], "DISCARD_NONMATCHING": False, "OUTPUT": "memory:"})["OUTPUT"]
kolom = ringkas.fields()[ringkas.fields().count() - 1].name()
print("Kolom hasil ringkasan:", kolom)
print("Total count:", int(sum((f[kolom] or 0) for f in ringkas.getFeatures())))

# 5. Select by location: petak yang beririsan dengan titik panas
processing.run("native:selectbylocation", {"INPUT": petak, "PREDICATE": [0], "INTERSECT": titik, "METHOD": 0})
print("Petak terpilih:", sorted(f["KODE_PETAK"] for f in petak.selectedFeatures()))
petak.removeSelection()

# 6. Jalan terdekat dari tiap titik panas
dekat = processing.run("native:joinbynearest", {"INPUT": titik, "INPUT_2": jalan, "FIELDS_TO_COPY": ["NAMA"], "NEIGHBORS": 1, "OUTPUT": "memory:"})["OUTPUT"]
baris = sorted((f["distance"], f["ID_Titik"], f["NAMA"]) for f in dekat.getFeatures())
print("Terdekat ke jalan: %s %.1f m (%s); terjauh: %s %.1f m" % (baris[0][1], baris[0][0], baris[0][2], baris[-1][1], baris[-1][0]))

# 7. Titik panas ke petak (are within, satu-ke-satu), lalu hitung titik per jenis
titik_petak = processing.run("native:joinattributesbylocation", {"INPUT": titik, "PREDICATE": [5], "JOIN": petak, "JOIN_FIELDS": ["KODE_PETAK", "JENIS"], "METHOD": 1, "OUTPUT": "memory:"})["OUTPUT"]
print("Titik di petak P07:", sorted(f["ID_Titik"] for f in titik_petak.getFeatures() if f["KODE_PETAK"] == "P07"))
per_jenis = {}
for f in titik_petak.getFeatures():
    per_jenis.setdefault(f["JENIS"], []).append(f["ID_Titik"])
for jenis, daftar in sorted(per_jenis.items()):
    print("Titik pada %-7s: %d %s" % (jenis, len(daftar), sorted(daftar)))

Penjelasan baris demi baris. Baris 1-28 adalah blok pembuka. Dua perkakas Python baru muncul di sini. Counter adalah penghitung otomatis: ia menghitung berapa kali tiap nilai muncul. setdefault adalah perintah kamus: "ambil isi kunci ini, dan bila belum ada, buat dulu dengan isi awal ini."

Hasil. Layar menampilkan (berdasarkan uji pada QGIS 4.0.2):

Join dengan intersect : 30 baris (petak ada 16)
Join dengan are within: 16 baris (petak ada 16)
Petak per KPH: {'KPH Alpha': 4, 'KPH Beta': 4, 'KPH Gamma': 8}
Titik per petak: {'P05': 1, 'P06': 1, 'P03': 1, 'P07': 3, 'P08': 1, 'P10': 1, 'P11': 1, 'P12': 1, 'P13': 1, 'P14': 1} jumlah = 12
Titik per KPH: {'KPH Alpha': 2, 'KPH Beta': 5, 'KPH Gamma': 5}
Kolom hasil ringkasan: ID_Titik_count
Total count: 12
Petak terpilih: ['P03', 'P05', 'P06', 'P07', 'P08', 'P10', 'P11', 'P12', 'P13', 'P14']
Terdekat ke jalan: HS-002 12.1 m (Jalan Inspeksi); terjauh: HS-001 225.2 m
Titik di petak P07: ['HS-004', 'HS-006', 'HS-011']
Titik pada Jati   : 4 ['HS-001', 'HS-003', 'HS-005', 'HS-010']
Titik pada Mahoni : 2 ['HS-002', 'HS-008']
Titik pada Pinus  : 2 ['HS-009', 'HS-012']
Titik pada Sengon : 4 ['HS-004', 'HS-006', 'HS-007', 'HS-011']

Hal terpenting pada hasil ini: join dengan intersect menghasilkan 30 baris, padahal petak hanya 16. Selisihnya 14 baris. Itu berasal dari petak di tepi batas KPH yang menyentuh KPH tetangganya, sehingga muncul baris tambahan. Dengan are within jumlahnya tepat 16.

Semua jumlah titik juga konsisten: 12 titik, baik dihitung per petak, per KPH, maupun lewat ringkasan.

Jebakan nomor predikat: dua alat, dua daftar #

Angka predikat tidak sama di semua alat. Skrip ini memakai 5 untuk are within pada alat join (baris 38 dan 71), tetapi 0 untuk intersect pada alat seleksi (baris 61). Itu benar.

Namun bila Anda memakai angka 5 di alat seleksi, artinya lain. Saya memeriksa daftarnya dengan membaca definisi alat di QGIS 4.0.2:

Saya juga mengujinya langsung pada data ini. Alat seleksi dengan angka 0 memilih 16 petak, angka 5 memilih 0 petak (karena 5 di alat itu berarti overlap), dan angka 6 memilih 16 petak. Alat join dengan angka 5 mencocokkan 16 petak, sedangkan angka 6 mencocokkan 0 (karena 6 di alat itu berarti cross). Kesimpulannya: selalu periksa daftar predikat tiap alat, jangan menyalin angka dari alat lain.

Bagian yang boleh diubah.

• Baris 38 dan 42, predikat dan metode. Ganti angka predikat sesuai tabel. Ganti METHOD 0 atau 1 sesuai kebutuhan.
• Baris 46, 55, 61, 66, 71, layer pasangan. Ganti layer yang dihubungkan.
• Baris 66, NEIGHBORS. Ganti 1 menjadi 2 untuk mencari dua jalan terdekat. Uji saya: hasilnya 24 baris untuk 12 titik, karena data jalan hanya berisi dua ruas.
• Baris 72, kode petak. Ganti "P07" dengan petak lain.
• Baris 55, SUMMARIES. Ganti 0 dengan nomor ringkasan lain. Nomor 5 adalah sum, dan nomor 6 mean.

Jebakan.

• Satu-ke-banyak menggandakan baris. Join intersect berbaris 30, bukan 16. Baris 37-39 memperlihatkannya sebagai pelajaran. Untuk "satu petak, satu KPH", pakai are within dan METHOD: 1.
• Nama kolom ringkasan ditebak. Baris 56 mengambil kolom paling akhir, bukan menyebut nama pasti. Hasil uji: ID_Titik_count. Bila Anda mengubah JOIN_FIELDS atau SUMMARIES, kolom terakhir bisa berbeda.
• Angka 16 di teks adalah angka mati. Baris 39 menulis "petak ada 16" sebagai teks tetap. Pada data lain, tulisan ini salah.
• Seleksi menandai layer asli. Baris 61 menyeleksi di layer petak. Baris 63 menghapus seleksinya. Bila skrip berhenti sebelum baris 63, petak tetap terpilih di proyek.
• Urutan kamus. Urutan isi kamus Titik per petak tidak urut abjad. Pakai sorted bila butuh urutan.

Latihan variasi. Cetak titik panas yang jatuh di petak P05, bukan P07.

Contoh jawaban. Ganti P07 menjadi P05 pada baris 72:

print("Titik di petak P05:", sorted(f["ID_Titik"] for f in titik_petak.getFeatures() if f["KODE_PETAK"] == "P05"))

Layar menampilkan Titik di petak P05: ['HS-008'].

Skrip 18.6: Merangkai alat menjadi model #

Berkas 6_1_buat_model.py. Skrip ini tidak memuat data. Ia membuat satu berkas model. Skrip ini tidak memakai blok pembuka.

Tujuan. Membuat model Processing bernama areal_efektif lewat kode: petak dikurangi sempadan sungai, lalu dikurangi zona permukiman, lalu dilebur per jenis. Model disimpan ke berkas areal_efektif.model3.

Bahan. Tidak ada data masukan. Skrip menulis satu berkas ke folder HASIL.

Model adalah resep tertulis dari rangkaian alat. Bayangkan resep masakan: bahan (masukan), langkah-langkah (alat), dan hidangan jadi (keluaran).

Sekali ditulis, resep bisa dimasak berkali-kali dengan jumlah bahan berbeda. Di QGIS, resep ini bisa disusun di Model Designer lewat klik. Skrip ini menyusunnya lewat kode.

Skrip utuh:

# -*- coding: utf-8 -*-
# Penulis: Badar Mubarok Yogaswara. Membuat model "Areal Efektif" lewat Python (setara dengan menyusunnya di Model Designer).
# Jalankan di Python Console QGIS. Hasil: berkas areal_efektif.model3 di folder HASIL.
from qgis.core import (QgsProcessingModelAlgorithm, QgsProcessingModelParameter, QgsProcessingModelChildAlgorithm,
                       QgsProcessingModelChildParameterSource, QgsProcessingModelOutput,
                       QgsProcessingParameterFeatureSource, QgsProcessingParameterNumber, QgsProcessing)

HASIL = r"C:\latihan\hasil"      # ubah sesuai komputer Anda; folder dibuat bila belum ada
import os
os.makedirs(HASIL, exist_ok=True)

m = QgsProcessingModelAlgorithm()
m.setName("areal_efektif")
m.setGroup("Latihan I1")

# --- Masukan model (yang diisi pengguna saat menjalankan)
p = QgsProcessingParameterFeatureSource("PETAK", "Petak", [QgsProcessing.TypeVectorPolygon])
m.addModelParameter(p, QgsProcessingModelParameter("PETAK"))
p = QgsProcessingParameterFeatureSource("SUNGAI", "Sungai", [QgsProcessing.TypeVectorLine])
m.addModelParameter(p, QgsProcessingModelParameter("SUNGAI"))
p = QgsProcessingParameterFeatureSource("PERMUKIMAN", "Permukiman", [QgsProcessing.TypeVectorPolygon])
m.addModelParameter(p, QgsProcessingModelParameter("PERMUKIMAN"))
p = QgsProcessingParameterNumber("JARAK_SUNGAI", "Jarak aman sungai (m)", QgsProcessingParameterNumber.Double, 25)
m.addModelParameter(p, QgsProcessingModelParameter("JARAK_SUNGAI"))
p = QgsProcessingParameterNumber("JARAK_PERMUKIMAN", "Jarak aman permukiman (m)", QgsProcessingParameterNumber.Double, 20)
m.addModelParameter(p, QgsProcessingModelParameter("JARAK_PERMUKIMAN"))

MP = QgsProcessingModelChildParameterSource


def langkah(id_, alg, deskripsi, params):
    c = QgsProcessingModelChildAlgorithm(alg)
    c.setChildId(id_)
    c.setDescription(deskripsi)
    for k, v in params.items():
        c.addParameterSources(k, [v])
    return c


# --- Langkah 1 dan 2: buffer sungai dan permukiman (dengan dissolve)
m.addChildAlgorithm(langkah("buf_sungai", "native:buffer", "Buffer sungai", {
    "INPUT": MP.fromModelParameter("SUNGAI"), "DISTANCE": MP.fromModelParameter("JARAK_SUNGAI"),
    "DISSOLVE": MP.fromStaticValue(True)}))
m.addChildAlgorithm(langkah("buf_kampung", "native:buffer", "Buffer permukiman", {
    "INPUT": MP.fromModelParameter("PERMUKIMAN"), "DISTANCE": MP.fromModelParameter("JARAK_PERMUKIMAN"),
    "DISSOLVE": MP.fromStaticValue(True)}))
# --- Langkah 3: petak dikurangi sempadan sungai
m.addChildAlgorithm(langkah("kurang_sungai", "native:difference", "Petak dikurangi sempadan", {
    "INPUT": MP.fromModelParameter("PETAK"), "OVERLAY": MP.fromChildOutput("buf_sungai", "OUTPUT")}))
# --- Langkah 4: hasilnya dikurangi zona permukiman (keluaran model: Petak_Efektif)
c4 = langkah("petak_efektif", "native:difference", "Dikurangi zona permukiman", {
    "INPUT": MP.fromChildOutput("kurang_sungai", "OUTPUT"), "OVERLAY": MP.fromChildOutput("buf_kampung", "OUTPUT")})
c4.setModelOutputs({"Petak_Efektif": QgsProcessingModelOutput("Petak_Efektif", "Petak_Efektif")})
c4.modelOutput("Petak_Efektif").setChildId("petak_efektif")
c4.modelOutput("Petak_Efektif").setChildOutputName("OUTPUT")
m.addChildAlgorithm(c4)
# --- Langkah 5: dissolve per jenis (keluaran model: Efektif_per_Jenis)
c5 = langkah("per_jenis", "native:dissolve", "Dissolve per jenis", {
    "INPUT": MP.fromChildOutput("petak_efektif", "OUTPUT"), "FIELD": MP.fromStaticValue(["JENIS"])})
c5.setModelOutputs({"Efektif_per_Jenis": QgsProcessingModelOutput("Efektif_per_Jenis", "Efektif_per_Jenis")})
c5.modelOutput("Efektif_per_Jenis").setChildId("per_jenis")
c5.modelOutput("Efektif_per_Jenis").setChildOutputName("OUTPUT")
m.addChildAlgorithm(c5)

berkas = os.path.join(HASIL, "areal_efektif.model3")
print("Simpan model:", m.toFile(berkas), berkas)
print("Langkah:", [c.description() for c in m.childAlgorithms().values()])

Penjelasan baris demi baris. Skrip ini punya banyak perkakas baru yang namanya panjang. Anda tidak perlu menghafalnya. Cukup pahami pola: masukan dulu, lalu langkah, lalu simpan.

Hasil. Layar menampilkan (berdasarkan uji pada QGIS 4.0.2):

Simpan model: True ...\hasil\areal_efektif.model3
Langkah: ['Buffer permukiman', 'Buffer sungai', 'Petak dikurangi sempadan', 'Dissolve per jenis', 'Dikurangi zona permukiman']

Jalur di baris pertama mengikuti isi HASIL di komputer Anda. Di sini saya menyingkatnya. Berkas model yang dihasilkan berukuran sekitar 26 kilobyte.

Perhatikan urutan di baris kedua: bukan urutan penulisan di skrip. Pada hasil uji, daftar itu tampil menurut urutan abjad kode langkah (buf_kampung, buf_sungai, kurang_sungai, per_jenis, petak_efektif). Urutan itu tidak memengaruhi cara model bekerja, karena tiap langkah menyebut sumbernya sendiri.

Bagian yang boleh diubah.

• Baris 8, HASIL. Ganti dengan folder Anda.
• Baris 13-14, nama dan grup model. Ganti sesuka Anda.
• Baris 23 dan 25, angka bawaan. Ganti 25 dan 20 untuk mengubah isian awal model.
• Baris 59, ["JENIS"]. Ganti dengan kolom pengelompok lain.

Jebakan.

• Folder dibuat otomatis. Baris 10 membuat folder HASIL bila belum ada. Saya mengujinya dengan folder yang belum ada: folder dibuat, dan berkas model tersimpan. Akibatnya, salah ketik pada HASIL tidak menimbulkan galat. Skrip membuat folder baru di tempat yang salah. Periksa jalur yang dicetak di hasil.
• Nama keluaran model. Nama yang ditulis di skrip ini Petak_Efektif dan Efektif_per_Jenis. Skrip 18.7 memakai nama huruf kecil (petak_efektif). Penjelasannya ada di Skrip 18.7.
• Belum teruji di antarmuka. Berkas model ini berhasil dibaca dan dijalankan oleh Skrip 18.7. Saya belum membukanya di Model Designer, jadi tampilan kotak-kotaknya belum saya lihat. [CEK]
• QgsProcessingParameterNumber.Double. Penulisan ini berjalan baik di QGIS 4.0.2. Pada versi lain nama pilihannya bisa berbeda. [CEK]
• Menimpa model lama. Menjalankan ulang skrip menulis berkas yang sama. Salin dulu bila Anda sudah mengubah model di Model Designer.

Latihan variasi. Ubah jarak bawaan model menjadi 50 meter untuk sungai dan 30 meter untuk permukiman.

Contoh jawaban. Ganti baris 23 dan 25:

p = QgsProcessingParameterNumber("JARAK_SUNGAI", "Jarak aman sungai (m)", QgsProcessingParameterNumber.Double, 50)

p = QgsProcessingParameterNumber("JARAK_PERMUKIMAN", "Jarak aman permukiman (m)", QgsProcessingParameterNumber.Double, 30)

Saya membuka kembali model hasilnya, lalu mencetak nilai bawaannya. Layar menampilkan bawaan: 50 30.

Skrip 18.7: Laporan areal efektif dalam satu skrip #

Berkas 7_1_laporan_terpadu.py. Skrip ini memuat blok pembuka (baris 1-28), lalu menjalankan model dari Skrip 18.6.

Tujuan. Menjalankan model areal efektif dari berkas, memberi nama KPH pada hasilnya, merekap per KPH dan jenis, menyimpan rekap ke GeoPackage, dan menghitung titik panas yang jatuh di areal efektif.

Bahan. Lima layer: Petak, Sungai, Permukiman (paket-i1), Batas_KPH, dan Hotspot_Nasional (paket-data). Jalankan Skrip 18.6 lebih dulu supaya berkas areal_efektif.model3 ada di folder HASIL. Skrip ini menghasilkan tiga berkas GeoPackage di folder itu.

Skrip ini menjawab pertanyaan Kepala Seksi di awal bab. Ia menyatukan semua hal sebelumnya: memuat, menjalankan resep, menghubungkan lokasi, merekap, dan menyeleksi.

Skrip utuh:

# -*- coding: utf-8 -*-
# Penulis: Badar Mubarok Yogaswara. Dijalankan di Python Console QGIS. Data sintetis, bukan kondisi lapangan nyata.
from qgis.core import QgsVectorLayer, QgsProject, QgsProperty
import processing

PAKET = r"C:\latihan\paket-i1"      # ubah sesuai letak folder di komputer Anda
DATA = r"C:\latihan\paket-data"


def muat(berkas, lapisan, nama):
    """Memuat satu layer GeoPackage dan menambahkannya ke proyek."""
    lyr = QgsVectorLayer("%s|layername=%s" % (berkas, lapisan), nama, "ogr")
    if not lyr.isValid():
        raise RuntimeError("Layer tidak terbaca: " + berkas)
    QgsProject.instance().addMapLayer(lyr)
    return lyr


def luas_ha(lyr):
    """Jumlah luas semua poligon dalam hektare (CRS harus UTM, satuan meter)."""
    return sum(f.geometry().area() for f in lyr.getFeatures()) / 10000.0


def buf(lyr, jarak):
    """Buffer dengan dissolve, hasil sementara."""
    return processing.run("native:buffer", {"INPUT": lyr, "DISTANCE": jarak, "DISSOLVE": True, "OUTPUT": "memory:"})["OUTPUT"]


import os
from collections import Counter
from qgis.core import QgsProcessingModelAlgorithm
from qgis.core import QgsProcessingModelAlgorithm

HASIL = r"C:\latihan\hasil"      # folder keluaran; dibuat bila belum ada
os.makedirs(HASIL, exist_ok=True)

petak = muat(PAKET + r"\Petak.gpkg", "Petak", "Petak")
sungai = muat(PAKET + r"\Sungai.gpkg", "Sungai", "Sungai")
kampung = muat(PAKET + r"\Permukiman.gpkg", "Permukiman", "Permukiman")
kph = muat(DATA + r"\Batas_KPH.gpkg", "Batas_KPH", "Batas_KPH")
titik = muat(DATA + r"\Hotspot_Nasional.gpkg", "Hotspot_Nasional", "Hotspot_Nasional")

# 1. Jalankan model dari Bab 6 (berkas .model3 hasil skrip 6_1 atau hasil susunan Anda sendiri)
model = QgsProcessingModelAlgorithm()
assert model.fromFile(HASIL + r"\areal_efektif.model3"), "Model tidak terbaca; jalankan skrip 6_1 lebih dulu"
hasil = processing.run(model, {
    "PETAK": petak, "SUNGAI": sungai, "PERMUKIMAN": kampung,
    "JARAK_SUNGAI": 25, "JARAK_PERMUKIMAN": 20,
    "petak_efektif": HASIL + r"\Petak_Efektif.gpkg",
    "efektif_per_jenis": HASIL + r"\Efektif_per_Jenis.gpkg"})
efektif = QgsVectorLayer(hasil["petak_efektif"], "Petak_Efektif", "ogr")
print("Areal efektif: %.2f ha dari %.2f ha" % (luas_ha(efektif), luas_ha(petak)))

# 2. Spatial join: setiap potongan petak diberi nama KPH
efektif_kph = processing.run("native:joinattributesbylocation", {"INPUT": efektif, "PREDICATE": [5], "JOIN": kph, "JOIN_FIELDS": ["NAMA_KPH"], "METHOD": 1, "OUTPUT": "memory:"})["OUTPUT"]

# 3. Dissolve per KPH dan jenis, simpan ke GeoPackage
rekap = processing.run("native:dissolve", {"INPUT": efektif_kph, "FIELD": ["NAMA_KPH", "JENIS"], "OUTPUT": HASIL + r"\Rekap_KPH_Jenis.gpkg"})["OUTPUT"]
rekap = QgsVectorLayer(rekap, "Rekap_KPH_Jenis", "ogr")
print("%-10s %-7s %s" % ("KPH", "Jenis", "Efektif (ha)"))
for f in sorted(rekap.getFeatures(), key=lambda f: (f["NAMA_KPH"], f["JENIS"])):
    print("%-10s %-7s %.2f" % (f["NAMA_KPH"], f["JENIS"], f.geometry().area() / 10000.0))
per_kph = Counter()
for f in efektif_kph.getFeatures():
    per_kph[f["NAMA_KPH"]] += f.geometry().area() / 10000.0
print("Total per KPH (ha):", {k: round(v, 2) for k, v in sorted(per_kph.items())})

# 4. Select by location: titik panas yang jatuh di areal efektif
processing.run("native:selectbylocation", {"INPUT": titik, "PREDICATE": [0], "INTERSECT": efektif, "METHOD": 0})
print("Titik panas di areal efektif: %d dari %d" % (titik.selectedFeatureCount(), titik.featureCount()))
titik.removeSelection()

Penjelasan baris demi baris. Baris 1-28 adalah blok pembuka.

Hasil. Layar menampilkan (berdasarkan uji pada QGIS 4.0.2):

Areal efektif: 12.39 ha dari 16.00 ha
KPH        Jenis   Efektif (ha)
KPH Alpha  Jati    1.00
KPH Alpha  Mahoni  2.00
KPH Alpha  Pinus   1.00
KPH Beta   Jati    0.80
KPH Beta   Pinus   0.80
KPH Beta   Sengon  1.57
KPH Gamma  Jati    2.85
KPH Gamma  Mahoni  1.63
KPH Gamma  Pinus   0.28
KPH Gamma  Sengon  0.46
Total per KPH (ha): {'KPH Alpha': 4.0, 'KPH Beta': 3.17, 'KPH Gamma': 5.22}
Titik panas di areal efektif: 7 dari 12

Inilah jawaban untuk Kepala Seksi. Dari 16,00 hektare, 12,39 hektare berada di luar jarak aman. Jumlah ketiga KPH, 4,00 ditambah 3,17 ditambah 5,22, sama dengan 12,39. Sebanyak 7 dari 12 titik panas jatuh di areal efektif.

Bagian yang boleh diubah.

• Baris 34, HASIL. Ganti dengan folder Anda. Folder harus sama dengan yang dipakai Skrip 18.6.
• Baris 48, jarak. Angka 25 dan 20 adalah jarak aman. Ganti untuk skenario lain.
• Baris 58, kolom pengelompok. Ganti ["NAMA_KPH", "JENIS"] dengan kolom lain.
• Baris 49-50, nama berkas keluaran. Ganti nama berkas GeoPackage.

Jebakan.

• Pengaman berkas model (baris 45). Bila areal_efektif.model3 tidak ada di HASIL, skrip berhenti di baris 45 dengan AssertionError: Model tidak terbaca; jalankan skrip 6_1 lebih dulu. Pesan menyebut skrip 6_1, yaitu Skrip 18.6 di bab ini. Pada uji saya, folder HASIL kosong. Tanpa assert, kegagalan tidak terlihat di baris itu. Skrip baru berhenti di baris 51 dengan KeyError: 'petak_efektif', yang tidak menyebut berkas model. Bila pesan pengaman muncul, jalankan Skrip 18.6 lebih dulu.
• Nama keluaran model berubah huruf. Skrip 18.6 menamai keluaran Petak_Efektif. Hasil processing.run memakai kunci petak_efektif. Hasil uji saya: kunci yang ada adalah CHILD_INPUTS, CHILD_RESULTS, efektif_per_jenis, dan petak_efektif. Dugaan saya, QGIS membuat nama aman dari nama keluaran, yaitu huruf kecil dan garis bawah. [kemungkinan] Bila Anda mengganti nama keluaran di model, kunci di baris 49-51 juga berubah.
• Potongan yang melintasi batas KPH. Baris 55 memakai are within. Potongan petak yang melintasi batas dua KPH tidak berada seluruhnya di satu KPH. Potongan seperti itu kemungkinan tidak mendapat nama KPH. Pada data ini hal itu tidak terjadi: dari 16 potongan, 0 yang tanpa nama KPH, dan jumlah per KPH sama dengan total areal efektif. Pada data lain, periksa dulu. [kemungkinan]
• Berkas lama ditimpa. Tiga berkas GeoPackage di HASIL ditulis ulang tiap kali.
• Pemanggilan ganda. Baris 31-32 mengimpor hal yang sama dua kali. Bukan galat, hanya kelebihan.
• Komentar menyebut "Bab 6". Baris 43 merujuk ke bab di Buku I1. Di buku ini, modelnya dibuat oleh Skrip 18.6.

Latihan variasi. Jalankan skenario "bagaimana bila" dengan jarak aman sungai 50 meter dan permukiman 30 meter.

Contoh jawaban. Ganti baris 48:

    "JARAK_SUNGAI": 50, "JARAK_PERMUKIMAN": 30,

Layar menampilkan Areal efektif: 9.69 ha dari 16.00 ha, total per KPH {'KPH Alpha': 3.71, 'KPH Beta': 2.34, 'KPH Gamma': 3.64}, dan Titik panas di areal efektif: 2 dari 12. Baris KPH Gamma Pinus 0.00 juga muncul: areal Pinus di KPH Gamma menyusut sampai 0,00 hektare setelah dibulatkan.

Cek paham #

1. Mengapa join intersect pada Skrip 18.5 menghasilkan 30 baris, padahal petak hanya 16?
2. Angka 5 berarti apa pada alat joinattributesbylocation, dan apa pada selectbylocation?
3. Apa fungsi ["OUTPUT"] di ujung processing.run(...)?

Jawaban.

1. Petak di tepi batas KPH menyentuh KPH tetangganya. Dengan intersect dan satu-ke-banyak, petak seperti itu cocok dengan dua KPH, sehingga muncul dua baris.
2. Di joinattributesbylocation, 5 berarti are within. Di selectbylocation, 5 berarti overlap, sedangkan are within bernomor 6.
3. processing.run mengembalikan kamus hasil. ["OUTPUT"] mengambil layer hasil dari kamus itu.

Kesalahan umum #

• Membiarkan jalur contoh. PAKET, DATA, dan HASIL menunjuk folder di komputer penulis. Perbaikan: ganti ketiganya, dan ganti garis miring terbalik bila Anda memakai macOS atau Ubuntu.
• Menyalin nomor predikat dari alat lain. Perbaikan: buka daftar predikat tiap alat, atau gunakan tabel di Skrip 18.5.
• Tidak membaca hasil, hanya mencari galat. Beberapa salah ketik tidak menimbulkan galat, tetapi mengubah hasil (kolom ORDEE pada rumus buffer, atau nama kolom yang tidak dijaga assert). Perbaikan: bandingkan angka hasil dengan perkiraan Anda, dan periksa nama kolom dengan assert seperti di baris 38 Skrip 18.4.
• Menghitung luas pada layer berkoordinat derajat. Perbaikan: pastikan CRS dalam meter sebelum memakai luas_ha.
• Menjalankan Skrip 18.7 sebelum 18.6. Perbaikan: jalankan 18.6 dulu, atau pastikan berkas model ada.

Ringkasan dan latihan #

• Geoprocessing mengolah layer lewat alat. Di skrip, processing.run("nama_alat", {parameter})["OUTPUT"] memanggil satu alat.
• Hasil satu alat boleh langsung menjadi masukan alat berikutnya, sehingga rantai analisis bisa ditulis dalam beberapa baris.
• Buffer membuat zona aman, overlay (Intersect, Difference, Union) menggunting tumpukan, dissolve melebur batas, dan spatial join menghubungkan data menurut letak.
• Nomor predikat berbeda antar alat. Periksa tabelnya.
• Model menyimpan rangkaian alat sebagai resep. Satu model bisa dijalankan berulang dengan jarak berbeda.

Latihan. Buat salinan Skrip 18.7 untuk skenario sempadan sungai 15 meter dan zona permukiman 10 meter. Catat luas areal efektif dan jumlah titik panas yang jatuh di dalamnya. Bagian mana yang Anda ubah?

Tabel perbandingan: alat geoprocessing di bab ini #