Lewati ke isi
Profil penulisSeri Buku GIS Kehutanan dan Pertanian/ I4
Tampilkan bagian untuk:

BAB 1: Python di Dalam QGIS: Bicara dengan Peta Lewat Perintah

#Studi kasus: "Daftar petak sebelum rapat"

Kepala Seksi menghampiri meja Anda lima menit sebelum rapat. "Tolong daftar petak KPH Contoh, lengkap dengan luasnya dalam hektare." Biasanya Anda membuka tabel atribut, menghitung luas, lalu mengetik ulang di lembar lain. Hari ini Anda mencoba cara baru: satu skrip pendek yang membaca layer lalu menulis daftarnya sendiri. Bulan depan, skrip yang sama bisa dijalankan lagi.

#Konsep: Python di QGIS dalam tiga kalimat

Python adalah bahasa perintah tertulis, dan PyQGIS adalah cara Python berbicara dengan QGIS. Hampir semua yang bisa Anda klik di QGIS punya perintah padanan di PyQGIS. Memakai menu itu seperti memesan makanan lewat pelayan, sedangkan skrip seperti daftar pesanan tertulis yang bisa diserahkan lagi setiap hari.

Istilah baru bab ini:

  • Skrip: berkas teks berisi perintah Python, biasanya berakhiran .py.
  • Konsol Python: jendela di QGIS tempat perintah diketik atau dijalankan.
  • Proyek: proyek QGIS yang sedang terbuka, dilihat dari Python lewat QgsProject.
  • Layer: satu data di proyek, bisa vektor (QgsVectorLayer) atau raster (QgsRasterLayer).
  • Fitur: satu baris di layer vektor, yaitu satu objek peta.
Ilustrasi 1.1: Dua cara memberi perintah Python
Skema: konsol untuk mencoba satu baris, editor skrip untuk perintah yang disimpan, lalu QGIS menjalankan dan memberi hasil

#Bagian A: QGIS

#Bagian A: PyQGIS di QGIS

#A.1 Membuka Konsol Python

  1. Buka QGIS, lalu pilih Plugins ► Python Console. Pintasannya Ctrl+Alt+P. Panel konsol muncul.
  2. Ketik print("Halo KPH Contoh"), lalu tekan Enter. Tulisan Halo KPH Contoh muncul di konsol.
  3. Klik tombol Show Editor di bilah konsol. Editor skrip terbuka di samping konsol.
  4. Tempel sebuah skrip di editor. Jalankan seluruhnya dengan Ctrl+Shift+E. Untuk menjalankan baris yang disorot saja, tekan Ctrl+E.

Hasilnya: Anda punya dua tempat bekerja. Konsol cocok untuk mencoba satu baris. Editor cocok untuk skrip yang akan Anda simpan dan pakai lagi.

#A.2 Lima hal Python yang cukup untuk memulai

Anda tidak perlu menguasai Python untuk mengikuti buku ini. Lima hal di bawah sudah cukup. Istilah barunya: variabel, daftar, kamus, perulangan, dan fungsi.

  • Variabel adalah kotak bernama yang menyimpan nilai.
  • Daftar adalah kumpulan nilai berurutan. Nomor urutnya dimulai dari 0.
  • Kamus adalah pasangan kunci dan nilai, seperti daftar isi berisi nomor halaman.
  • Perulangan (for) mengulang perintah untuk tiap isi daftar.
  • Fungsi adalah resep bernama yang bisa dipakai berulang kali.

Baris yang menjorok ke dalam (indentasi) menandai isi sebuah blok. Empat spasi sudah cukup. Coba skrip berikut di editor. Skrip ini belum memakai QGIS sama sekali.

PYTHON 1.0
# [SKRIP 1.0: Lima hal Python yang cukup untuk memulai]
# Penulis: Badar Mubarok Yogaswara
# Tidak memakai QGIS sama sekali. Boleh dicoba di Python Console baris demi baris.

# 1. Variabel: kotak bernama yang menyimpan nilai
nama = "KPH Alpha"
luas_ha = 4.0
print(nama, "luasnya", luas_ha, "ha")

# 2. Daftar: kumpulan nilai berurutan. Nomor urut dimulai dari 0, bukan 1.
petak = ["KPH Alpha", "KPH Beta", "KPH Gamma"]
print("Petak pertama:", petak[0])
print("Jumlah petak:", len(petak))

# 3. Kamus: pasangan kunci dan nilai, seperti daftar isi dengan nomor halaman
luas = {"KPH Alpha": 4.0, "KPH Beta": 4.0, "KPH Gamma": 8.0}
print("Luas KPH Gamma:", luas["KPH Gamma"])

# 4. Perulangan for dan percabangan if. Baris yang menjorok ke dalam ikut dalam blok.
total = 0
for nama, ha in luas.items():
    total = total + ha
    if ha > 5:
        print(nama, "termasuk petak besar")
print("Total:", total, "ha")


# 5. Fungsi: resep yang diberi nama supaya bisa dipakai berulang kali
def ke_hektar(meter_persegi):
    return meter_persegi / 10000


print("40000 m2 =", ke_hektar(40000), "ha")

Hasil uji di QGIS 4.0.2:

HASIL UJI SKRIP 1.0
KPH Alpha luasnya 4.0 ha
Petak pertama: KPH Alpha
Jumlah petak: 3
Luas KPH Gamma: 8.0
KPH Gamma termasuk petak besar
Total: 16.0 ha
40000 m2 = 4.0 ha

#A.3 Memuat layer dan mengambilnya dengan nama

Bayangkan rak buku. Anda meletakkan buku di rak (memuat layer), lalu mengambilnya lagi dengan menyebut judulnya (mengambil layer dengan nama). Rak itu adalah proyek.

Ilustrasi 1.2: Isi proyek dilihat dari Python
Skema: proyek berisi layer vektor dan raster, layer vektor berisi fitur, dan pencarian nama mengembalikan daftar

Skrip 1.1 memuat satu layer vektor dan satu layer raster. Lalu skrip mengambil keduanya kembali dan mencetak sifatnya.

PYTHON (PYQGIS) 1.1
# [SKRIP 1.1: Memuat layer lewat Python dan mengambilnya kembali dengan nama]
# Penulis: Badar Mubarok Yogaswara
# Jalankan di Plugins > Python Console. Ganti DATA dengan folder paket-i4 di komputer Anda.
from qgis.core import QgsProject, QgsVectorLayer, QgsRasterLayer, QgsWkbTypes

DATA = "D:/Latihan/paket-i4"
HASIL = "D:/Latihan/hasil_i4"

proyek = QgsProject.instance()          # proyek yang sedang terbuka
proyek.removeAllMapLayers()             # kosongkan dulu agar tidak dobel

# 1. Muat satu layer vektor dan satu raster. Argumen: jalur, nama di panel Layers, penyedia data.
batas = QgsVectorLayer(DATA + "/Batas_KPH.gpkg", "Batas_KPH", "ogr")
dsm = QgsRasterLayer(DATA + "/DSM_Drone.tif", "DSM_Drone", "gdal")

# 2. Periksa dulu: layer yang gagal dibuka tidak memicu galat, hanya isValid() bernilai False.
print("Batas_KPH sah?", batas.isValid())
print("DSM_Drone sah?", dsm.isValid())

# 3. Masukkan ke proyek supaya tampil di panel Layers.
proyek.addMapLayer(batas)
proyek.addMapLayer(dsm)

# 4. Ambil layer dengan namanya. mapLayersByName() mengembalikan DAFTAR, ambil isi pertama.
layer = proyek.mapLayersByName("Batas_KPH")[0]
print("Nama:", layer.name())
print("CRS:", layer.crs().authid())
print("Jenis geometri:", QgsWkbTypes.displayString(layer.wkbType()))
print("Jumlah fitur:", layer.featureCount())
print("Nama field:", [f.name() for f in layer.fields()])

# 5. Sifat sebuah raster
r = proyek.mapLayersByName("DSM_Drone")[0]
print("Raster:", r.width(), "x", r.height(), "piksel,", r.bandCount(), "band")
print("Ukuran piksel (m):", r.rasterUnitsPerPixelX())
print("Jumlah layer di proyek:", len(proyek.mapLayers()))

Hasil uji di QGIS 4.0.2:

HASIL UJI SKRIP 1.1
Batas_KPH sah? True
DSM_Drone sah? True
Nama: Batas_KPH
CRS: EPSG:32749
Jenis geometri: Polygon
Jumlah fitur: 3
Nama field: ['fid', 'NAMA_KPH']
Raster: 800 x 800 piksel, 1 band
Ukuran piksel (m): 0.5
Jumlah layer di proyek: 2

#A.4 Membaca fitur dan atribut

Layer vektor seperti tabel. Tiap baris adalah satu fitur. Skrip 1.2 membaca tiap baris, menghitung luas, memilih satu petak dengan ekspresi, lalu menulis hasil ke berkas CSV.

PYTHON (PYQGIS) 1.2
# [SKRIP 1.2: Membaca fitur dan atribut, memilih dengan ekspresi, menulis CSV]
# Penulis: Badar Mubarok Yogaswara
import csv
import os
from qgis.core import QgsProject, QgsVectorLayer

DATA = "D:/Latihan/paket-i4"
HASIL = "D:/Latihan/hasil_i4"
os.makedirs(HASIL, exist_ok=True)

batas = QgsVectorLayer(DATA + "/Batas_KPH.gpkg", "Batas_KPH", "ogr")
QgsProject.instance().addMapLayer(batas)

# 1. Satu fitur = satu baris tabel atribut. getFeatures() membacanya satu per satu.
baris = []
for fitur in batas.getFeatures():
    nama = fitur["NAMA_KPH"]
    luas_ha = fitur.geometry().area() / 10000      # m2 menjadi hektare (CRS harus berupa meter)
    baris.append((nama, round(luas_ha, 2)))
    print(nama, round(luas_ha, 2), "ha")

total = sum(luas for _, luas in baris)
print("Total:", total, "ha")

# 2. Memilih fitur dengan ekspresi, sama seperti Select by Expression
batas.selectByExpression('"NAMA_KPH" = \'KPH Beta\'')
print("Terpilih:", batas.selectedFeatureCount(), "fitur")
for fitur in batas.selectedFeatures():
    print("  ->", fitur["NAMA_KPH"])
batas.removeSelection()

# 3. Menulis hasil ke berkas CSV
jalur = HASIL + "/luas_kph.csv"
with open(jalur, "w", newline="", encoding="utf-8") as f:
    penulis = csv.writer(f)
    penulis.writerow(["NAMA_KPH", "LUAS_HA"])
    penulis.writerows(baris)
print("CSV tersimpan:", os.path.exists(jalur))

Hasil uji di QGIS 4.0.2:

HASIL UJI SKRIP 1.2
KPH Alpha 4.0 ha
KPH Beta 4.0 ha
KPH Gamma 8.0 ha
Total: 16.0 ha
Terpilih: 1 fitur
  -> KPH Beta
CSV tersimpan: True

Angka luas tepat karena layer berkoordinat meter (UTM). Pada layer berderajat (EPSG:4326), geometry().area() menghasilkan derajat persegi, bukan meter persegi. Ingat pelajaran CRS dari seri B1.

#A.5 Menyimpan proyek dan membukanya lagi

Skrip 1.3 menyusun proyek dari nol: mengatur CRS, memuat empat layer, menyimpan sebagai .qgz, mengosongkan proyek, lalu membukanya kembali sebagai bukti.

PYTHON (PYQGIS) 1.3
# [SKRIP 1.3: Menyusun proyek lewat Python lalu menyimpan dan membukanya lagi]
# Penulis: Badar Mubarok Yogaswara
import os
from qgis.core import QgsProject, QgsVectorLayer, QgsRasterLayer, QgsCoordinateReferenceSystem

DATA = "D:/Latihan/paket-i4"
HASIL = "D:/Latihan/hasil_i4"
os.makedirs(HASIL, exist_ok=True)

proyek = QgsProject.instance()
proyek.removeAllMapLayers()
proyek.setCrs(QgsCoordinateReferenceSystem("EPSG:32749"))     # CRS proyek: UTM zona 49S
proyek.setTitle("KPH Contoh - proyek skrip")

daftar = [
    QgsVectorLayer(DATA + "/Batas_KPH.gpkg", "Batas_KPH", "ogr"),
    QgsVectorLayer(DATA + "/Sungai.gpkg", "Sungai", "ogr"),
    QgsVectorLayer(DATA + "/Hotspot_Nasional.gpkg", "Hotspot_Nasional", "ogr"),
    QgsRasterLayer(DATA + "/Citra_KPH.tif", "Citra_KPH", "gdal"),
]
for layer in daftar:
    if layer.isValid():
        proyek.addMapLayer(layer)
    else:
        print("GAGAL dimuat:", layer.name())

berkas = HASIL + "/Proyek_KPH.qgz"
print("Tersimpan?", proyek.write(berkas))          # write() mengembalikan True bila berhasil

# Buka kembali dan buktikan isinya utuh
proyek.clear()
print("Setelah dikosongkan:", len(proyek.mapLayers()), "layer")
print("Terbuka?", proyek.read(berkas))
print("Judul:", proyek.title(), "| CRS:", proyek.crs().authid())
print("Layer:", sorted(l.name() for l in proyek.mapLayers().values()))

Hasil uji di QGIS 4.0.2:

HASIL UJI SKRIP 1.3
Tersimpan? True
Setelah dikosongkan: 0 layer
Terbuka? True
Judul: KPH Contoh - proyek skrip | CRS: EPSG:32749
Layer: ['Batas_KPH', 'Citra_KPH', 'Hotspot_Nasional', 'Sungai']

Hasilnya: proyek yang dibuat lewat skrip sama utuhnya dengan proyek yang dibuat lewat klik.

#Bagian B: ArcGIS Pro

#Bagian B: ArcPy dasar di ArcGIS Pro

ArcPy adalah paket Python milik Esri. ArcGIS Pro memakai Python 3. Jendela perintahnya dibuka lewat tab Analysis ► Python ► Python Window. [CEK]

Padanannya dengan Bagian A: arcpy.ListFeatureClasses() untuk daftar data, arcpy.da.SearchCursor untuk membaca baris, arcpy.Raster untuk raster, dan arcpy.mp.ArcGISProject("CURRENT") untuk proyek yang terbuka. Bacalah skrip berikut dengan membandingkannya dengan Skrip 1.1 sampai 1.3.

PYTHON (ARCPY, PYTHON 3)
# [SKRIP 1.1-Pro: ArcPy dasar: membaca data, atribut, dan proyek]
# Penulis: Badar Mubarok Yogaswara
# [CEK] Cocokkan dengan versi ArcGIS Pro Anda.
import arcpy

DATA = r"D:\Latihan\paket-i4"
arcpy.env.workspace = DATA + r"\shp"          # folder Shapefile; daftar data dibaca dari sini

# 1. Daftar data dan sifatnya
print(arcpy.ListFeatureClasses())              # contoh: ['Batas_KPH.shp', 'Hotspot_Nasional.shp', ...]
batas = "Batas_KPH.shp"
d = arcpy.Describe(batas)
print(d.shapeType, d.spatialReference.name, d.spatialReference.factoryCode)   # Polygon, ..., 32749
print("Jumlah fitur:", int(arcpy.management.GetCount(batas)[0]))
print("Field:", [f.name for f in arcpy.ListFields(batas)])

# 2. Membaca tiap baris dengan kursor (padanan getFeatures() di QGIS)
with arcpy.da.SearchCursor(batas, ["NAMA_KPH", "SHAPE@AREA"]) as kursor:
    for nama, luas in kursor:
        print(nama, round(luas / 10000, 2), "ha")          # luas dalam satuan peta (m2) dibagi 10.000

# 3. Raster
r = arcpy.Raster(DATA + r"\Citra_KPH.tif")
print("Raster:", r.width, "x", r.height, "piksel,", r.bandCount, "band; piksel", r.meanCellWidth, "m")

# 4. Proyek: tambah data lalu ambil layer dengan nama
aprx = arcpy.mp.ArcGISProject("CURRENT")      # proyek yang sedang terbuka
peta = aprx.listMaps()[0]
peta.addDataFromPath(DATA + r"\shp\Batas_KPH.shp")
layer = peta.listLayers("Batas_KPH")[0]       # listLayers() juga mengembalikan DAFTAR
print("Layer ditemukan:", layer.name)
aprx.save()

Beda yang perlu diingat. listLayers("Batas_KPH") juga mengembalikan daftar, sama seperti mapLayersByName. Kursor SearchCursor membaca nilai lewat urutan nama field yang Anda sebut, bukan lewat nama di tiap baris.

#Bagian C: ArcMap 10.8

#Bagian C: ArcPy dasar di ArcMap 10.8

ArcMap memakai Python 2.7. Jendela perintahnya dibuka lewat menu Geoprocessing ► Python. [CEK] Pengelolaan dokumen peta memakai modul arcpy.mapping, bukan arcpy.mp.

Dua beda penting dengan Python 3. Pertama, pembagian dua bilangan bulat menghasilkan bilangan bulat, jadi tulis 10000.0 agar hasilnya desimal. Kedua, tidak ada f-string, jadi teks disusun dengan tanda %.

PYTHON (ARCPY, PYTHON 2.7)
# -*- coding: utf-8 -*-
# [SKRIP 1.1-ArcMap: ArcPy dasar: membaca data, atribut, dan dokumen peta (Python 2.7)]
# Penulis: Badar Mubarok Yogaswara
# [CEK] Cocokkan dengan versi ArcMap Anda.
import arcpy

DATA = r"D:\Latihan\paket-i4"
arcpy.env.workspace = DATA + r"\shp"

# 1. Daftar data dan sifatnya
print(arcpy.ListFeatureClasses())
batas = "Batas_KPH.shp"
d = arcpy.Describe(batas)
print("%s %s %s" % (d.shapeType, d.spatialReference.name, d.spatialReference.factoryCode))
print("Jumlah fitur: %d" % int(arcpy.GetCount_management(batas).getOutput(0)))
print("Field: %s" % [f.name for f in arcpy.ListFields(batas)])

# 2. Membaca tiap baris dengan kursor
with arcpy.da.SearchCursor(batas, ["NAMA_KPH", "SHAPE@AREA"]) as kursor:
    for nama, luas in kursor:
        print("%s %.2f ha" % (nama, luas / 10000.0))        # 10000.0: hasil bagi desimal di Python 2

# 3. Raster
r = arcpy.Raster(DATA + r"\Citra_KPH.tif")
print("Raster: %d x %d piksel, %d band, piksel %s m" % (r.width, r.height, r.bandCount, r.meanCellWidth))

# 4. Dokumen peta: tambah layer lalu ambil dengan nama
mxd = arcpy.mapping.MapDocument("CURRENT")
df = arcpy.mapping.ListDataFrames(mxd)[0]
arcpy.mapping.AddLayer(df, arcpy.mapping.Layer(DATA + r"\shp\Batas_KPH.shp"))
layer = arcpy.mapping.ListLayers(mxd, "Batas_KPH", df)[0]    # juga mengembalikan DAFTAR
print("Layer ditemukan: %s" % layer.name)
mxd.save()

#Cek paham

  1. Mengapa mapLayersByName("Batas_KPH") harus diikuti [0]?
  2. Apa beda "ogr" dan "gdal" saat membuat layer?
  3. Berkas data hilang, tetapi skrip tidak menampilkan galat. Bagaimana cara tahu bahwa layer gagal dibuka?

Jawaban:

  1. Karena hasilnya berupa daftar. Bisa ada beberapa layer dengan nama sama. [0] mengambil yang pertama.
  2. "ogr" untuk data vektor, "gdal" untuk raster.
  3. Periksa layer.isValid(). Nilai False berarti layer tidak terbuka.

#Kesalahan umum

  • Memakai garis miring balik di jalur. Tulis D:/Latihan/paket-i4, bukan D:\Latihan\paket-i4. Bila ingin garis miring balik, awali jalur dengan r seperti r"D:\Latihan".
  • Menjalankan skrip dua kali lalu muncul layer ganda. Kosongkan proyek lebih dulu dengan removeAllMapLayers(), seperti di Skrip 1.1.
  • Salah huruf pada nama layer. "batas_kph" berbeda dari "Batas_KPH". Hasilnya daftar kosong, lalu muncul galat IndexError saat memanggil [0].
  • Indentasi tidak rapi. Baris dalam blok for harus sama menjoroknya.

#Ringkasan dan latihan

Ringkasan: skrip adalah catatan perintah yang bisa diulang. QgsProject.instance() adalah proyek yang terbuka. Layer dimuat dengan QgsVectorLayer atau QgsRasterLayer, diperiksa dengan isValid(), lalu diambil dengan mapLayersByName()[0]. Tiap baris vektor dibaca dengan getFeatures().

Latihan:

  1. Ubah Skrip 1.2 agar hanya mencetak petak yang luasnya lebih dari 5 ha. Hasil yang benar: hanya KPH Gamma.
  2. Muat Sungai.gpkg, lalu cetak panjang tiap sungai dengan fitur.geometry().length(). Hasil yang benar: Sungai Tengah 420,3 m dan Anak Sungai Timur 289,5 m, total 709,8 m.

#Tabel perbandingan: PyQGIS dan ArcPy dasar

HalQGIS (Python 3)ArcGIS Pro (Python 3)ArcMap 10.8 (Python 2.7)
Jendela perintahPlugins ► Python ConsoleAnalysis ► Python ► Python Window [CEK]Geoprocessing ► Python [CEK]
Memuat data ke petaQgsVectorLayer lalu addMapLayerMap.addDataFromPath [CEK]arcpy.mapping.AddLayer [CEK]
Ambil layer dengan namamapLayersByName("x")[0]Map.listLayers("x")[0] [CEK]ListLayers(mxd, "x", df)[0] [CEK]
Baca tiap barisgetFeatures()arcpy.da.SearchCursor [CEK]arcpy.da.SearchCursor [CEK]
Simpan proyekQgsProject.write()aprx.save() [CEK]mxd.save() [CEK]
Cek layer sahisValid()arcpy.Exists() [CEK]arcpy.Exists() [CEK]