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

BAB 3: Skrip Batch: Satu Resep untuk Banyak Berkas

#Studi kasus: "Lima puluh ubin foto drone"

Tim drone menyerahkan puluhan ubin foto udara. Kepala Seksi minta NDVI untuk setiap ubin sebelum akhir pekan. Mengklik raster calculator satu per satu akan menghabiskan sore Anda. Lebih buruk lagi, QGIS sering membeku saat satu ubin diproses. Di buku ini dipakai tiga ubin agar mudah dicoba, tetapi skripnya sama untuk lima puluh.

#Konsep: batch dalam tiga kalimat

Batch berarti menjalankan pekerjaan yang sama berulang kali pada banyak masukan. Anda menulis resepnya sekali, lalu perulangan for memasaknya untuk tiap berkas. Untuk pekerjaan berat, QgsTask menyerahkan resep itu ke "juru masak" di latar belakang, sementara Anda, sang "pelayan", tetap bisa melayani QGIS.

Istilah baru bab ini:

  • Batch: pekerjaan sama, banyak masukan.
  • Pola berkas (glob): cara mencari banyak berkas sekaligus, misalnya *.tif.
  • Latar belakang: jalur kerja terpisah yang tidak membekukan layar QGIS.
  • QgsTask: tugas latar belakang milik QGIS.
  • Fungsi selesai: fungsi yang dipanggil setelah tugas berakhir.
Ilustrasi 3.1: Satu resep, banyak berkas
Skema: tiga ubin masuk ke satu perulangan for yang menghitung NDVI dan menyimpan tiga hasil

#Bagian A: QGIS

#Bagian A: Batch di QGIS

#A.1 Tanpa menulis kode: batch bawaan

QGIS punya batch bawaan. Klik kanan sebuah alat di Processing Toolbox, lalu pilih Execute as batch process. Anda mengisi tabel: satu baris untuk satu pekerjaan. Cara ini cocok untuk beberapa berkas. Untuk puluhan berkas, mengisi tabel itu melelahkan, dan nama keluaran sulit diatur. Di sinilah skrip unggul.

#A.2 Skrip batch sederhana

Skrip 3.1 mengumpulkan semua .tif di folder, menghitung NDVI tiap ubin, lalu mencetak rata-ratanya. Pada ubin ini band 1 adalah NIR dan band 2 adalah Merah.

PYTHON (PYQGIS) 3.1
# [SKRIP 3.1: Batch sederhana: NDVI untuk semua ubin ortofoto dalam satu folder]
# Penulis: Badar Mubarok Yogaswara
import glob
import os
import processing
from qgis.core import QgsRasterLayer

DATA = "D:/Latihan/paket-i4"
HASIL = "D:/Latihan/hasil_i4"
folder_masuk = DATA + "/Ortofoto_Mentah"
folder_keluar = HASIL + "/NDVI_ubin"
os.makedirs(folder_keluar, exist_ok=True)       # buat folder bila belum ada

# 1. Kumpulkan semua berkas .tif, diurutkan agar hasilnya selalu berurutan sama
daftar = sorted(glob.glob(folder_masuk + "/*.tif"))
print("Ditemukan", len(daftar), "ubin")

# 2. Ulangi pekerjaan yang sama untuk tiap berkas. Pada ubin ini band 1 = NIR dan band 2 = Merah.
berhasil = []
for jalur in daftar:
    nama = os.path.basename(jalur)                         # contoh: ubin_01.tif
    keluar = folder_keluar + "/NDVI_" + nama
    ubin = QgsRasterLayer(jalur, "ubin", "gdal")           # nama "ubin" dipakai di ekspresi
    processing.run("native:rastercalc", {
        "EXPRESSION": '("ubin@1" - "ubin@2") / ("ubin@1" + "ubin@2")',
        "LAYERS": [ubin], "OUTPUT": keluar})
    hasil = QgsRasterLayer(keluar, "hasil", "gdal")
    rata = processing.run("native:rasterlayerstatistics", {"INPUT": hasil, "BAND": 1})["MEAN"]
    print("  %s -> %s  | NDVI rata-rata %.3f" % (nama, os.path.basename(keluar), rata))
    berhasil.append(keluar)

print("Selesai:", len(berhasil), "dari", len(daftar), "berkas")

Hasil uji di QGIS 4.0.2:

HASIL UJI SKRIP 3.1
Ditemukan 3 ubin
  ubin_01.tif -> NDVI_ubin_01.tif  | NDVI rata-rata 0.471
  ubin_02.tif -> NDVI_ubin_02.tif  | NDVI rata-rata 0.471
  ubin_03.tif -> NDVI_ubin_03.tif  | NDVI rata-rata 0.467
Selesai: 3 dari 3 berkas

#A.3 Satu putaran untuk tiap fitur

Perulangan tidak harus atas berkas. Skrip 3.2 berputar atas fitur: untuk tiap KPH, skrip memilih satu poligon lalu memotong CHM dengannya. Hasilnya satu raster per KPH.

Ada satu trik penting. Nama KPH disisipkan ke ekspresi pemilih dengan %s. Tanda kutip tunggal di sekitar nama ditulis \', supaya Python tidak menganggapnya penutup teks.

PYTHON (PYQGIS) 3.2
# [SKRIP 3.2: Batch per fitur: potong raster tinggi tajuk (CHM) untuk tiap KPH]
# Penulis: Badar Mubarok Yogaswara
import os
import processing
from qgis.core import QgsVectorLayer, QgsRasterLayer

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

batas = QgsVectorLayer(DATA + "/Batas_KPH.gpkg", "Batas_KPH", "ogr")
dsm = QgsRasterLayer(DATA + "/DSM_Drone.tif", "DSM", "gdal")
dtm = QgsRasterLayer(DATA + "/DTM_Drone.tif", "DTM", "gdal")

# 1. CHM untuk seluruh area (satu kali saja, di luar pengulangan)
chm_jalur = processing.run("native:rastercalc", {
    "EXPRESSION": '"DSM@1" - "DTM@1"', "LAYERS": [dsm, dtm], "OUTPUT": "TEMPORARY_OUTPUT"})["OUTPUT"]

# 2. Satu putaran untuk tiap KPH: pilih satu poligon, lalu potong CHM dengannya
for fitur in batas.getFeatures():
    nama = fitur["NAMA_KPH"]                                 # contoh: KPH Alpha
    aman = nama.replace(" ", "_")                            # untuk nama berkas: KPH_Alpha
    satu_kph = processing.run("native:extractbyexpression", {
        "INPUT": batas, "EXPRESSION": '"NAMA_KPH" = \'%s\'' % nama,
        "OUTPUT": "TEMPORARY_OUTPUT"})["OUTPUT"]
    keluar = folder_keluar + "/CHM_" + aman + ".tif"
    processing.run("gdal:cliprasterbymasklayer", {
        "INPUT": chm_jalur, "MASK": satu_kph, "CROP_TO_CUTLINE": True,
        "NODATA": -9999, "OUTPUT": keluar})
    hasil = QgsRasterLayer(keluar, aman, "gdal")
    st = processing.run("native:rasterlayerstatistics", {"INPUT": hasil, "BAND": 1})
    print("%-10s %4d x %4d piksel | rata-rata %.2f m | tertinggi %.2f m" % (
        aman, hasil.width(), hasil.height(), st["MEAN"], st["MAX"]))
print("Berkas di folder:", sorted(os.listdir(folder_keluar)))

Hasil uji di QGIS 4.0.2:

HASIL UJI SKRIP 3.2
KPH_Alpha   400 x  400 piksel | rata-rata 2.60 m | tertinggi 13.94 m
KPH_Beta    400 x  400 piksel | rata-rata 2.69 m | tertinggi 13.96 m
KPH_Gamma   800 x  400 piksel | rata-rata 2.56 m | tertinggi 13.96 m
Berkas di folder: ['CHM_KPH_Alpha.tif', 'CHM_KPH_Alpha.tif.aux.xml', 'CHM_KPH_Beta.tif', 'CHM_KPH_Beta.tif.aux.xml', 'CHM_KPH_Gamma.tif']

Ukuran hasil bisa dicek dengan akal sehat. Alpha dan Beta berukuran 200 m x 200 m dengan piksel 0,5 m, yaitu 400 x 400 piksel. Gamma dua kali lebih lebar, yaitu 800 x 400 piksel.

#A.4 Batch di latar belakang dengan QgsTask

Perulangan di Skrip 3.1 membuat layar QGIS "menunggu" sampai semua ubin selesai. Pada lima puluh ubin, itu terasa seperti macet. Skrip 3.3 memindahkan pekerjaan yang sama ke latar belakang.

Ilustrasi 3.2: Program utama dan latar belakang
Skema: program utama mengirim tugas ke latar belakang, latar belakang mengerjakan alat berat, lalu hasil dikirim kembali ke fungsi selesai

Susunannya selalu dua fungsi:

  1. Fungsi tugas (hitung_ndvi) berjalan di latar belakang. Argumen pertamanya selalu task. Fungsi ini boleh menjalankan alat dan membaca berkas. Fungsi ini tidak boleh menyentuh proyek atau antarmuka QGIS.
  2. Fungsi selesai (selesai) berjalan di program utama setelah tugas berakhir. Di sinilah hasil ditambahkan ke proyek atau ditampilkan.
PYTHON (PYQGIS) 3.3
# [SKRIP 3.3: Batch di latar belakang dengan QgsTask supaya QGIS tidak macet]
# Penulis: Badar Mubarok Yogaswara
import glob
import os
import processing
from qgis.core import QgsApplication, QgsTask, QgsMessageLog, Qgis, QgsRasterLayer, QgsProject

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


def hitung_ndvi(task, folder_in, folder_out):
    """Berjalan di LATAR BELAKANG. Jangan menyentuh proyek atau antarmuka di sini."""
    daftar = sorted(glob.glob(folder_in + "/*.tif"))
    keluar = []
    for i, jalur in enumerate(daftar):
        if task.isCanceled():                       # pengguna menekan Cancel
            return None
        tujuan = folder_out + "/NDVI_" + os.path.basename(jalur)
        ubin = QgsRasterLayer(jalur, "ubin", "gdal")
        processing.run("native:rastercalc", {
            "EXPRESSION": '("ubin@1" - "ubin@2") / ("ubin@1" + "ubin@2")',
            "LAYERS": [ubin], "OUTPUT": tujuan})
        keluar.append(tujuan)
        task.setProgress(100 * (i + 1) / len(daftar))    # bilah kemajuan di pojok kanan bawah QGIS
    return keluar


def selesai(exception, hasil=None):
    """Berjalan di PROGRAM UTAMA setelah tugas berakhir. Aman untuk menyentuh proyek."""
    if exception is not None:
        QgsMessageLog.logMessage("Tugas gagal: %s" % exception, "Batch NDVI", Qgis.Critical)
    elif hasil is None:
        QgsMessageLog.logMessage("Tugas dibatalkan", "Batch NDVI", Qgis.Warning)
    else:
        for jalur in hasil:
            nama = os.path.basename(jalur)[:-4]
            QgsProject.instance().addMapLayer(QgsRasterLayer(jalur, nama, "gdal"))
        QgsMessageLog.logMessage("Selesai: %d berkas" % len(hasil), "Batch NDVI", Qgis.Success)
        print("Tugas selesai:", len(hasil), "berkas NDVI ditambahkan ke proyek")


tugas = QgsTask.fromFunction("Hitung NDVI semua ubin", hitung_ndvi, on_finished=selesai,
                             folder_in=folder_masuk, folder_out=folder_keluar)
QgsApplication.taskManager().addTask(tugas)
print("Tugas dikirim ke latar belakang. QGIS tetap bisa dipakai.")

Hasil uji di QGIS 4.0.2:

HASIL UJI SKRIP 3.3
Tugas dikirim ke latar belakang. QGIS tetap bisa dipakai.
Tugas selesai: 3 berkas NDVI ditambahkan ke proyek

#Bagian B: ArcGIS Pro

#Bagian B: Batch di ArcGIS Pro

Pro tidak punya padanan langsung QgsTask. Pekerjaan banyak berkas ditulis sebagai skrip batch biasa. Alat arcpy.ListRasters mendaftar raster di arcpy.env.workspace. Perulangan for memprosesnya satu per satu. Untuk pekerjaan panjang, kemajuan bisa ditampilkan lewat arcpy.SetProgressor bila skrip dijadikan alat skrip. [CEK] Esri juga menyediakan opsi paralel (Parallel Processing Factor) yang perlu Anda coba sendiri. [CEK]

PYTHON (ARCPY, PYTHON 3)
# [SKRIP 3.1-Pro: Batch: NDVI semua ubin dan potongan CHM per KPH dengan ArcPy]
# Penulis: Badar Mubarok Yogaswara
# [CEK] Cocokkan dengan versi Anda. Butuh ekstensi Spatial Analyst.
import os
import arcpy
from arcpy.sa import Raster, Float, ExtractByMask

DATA = r"D:\Latihan\paket-i4"
HASIL = r"D:\Latihan\hasil_i4"
arcpy.env.overwriteOutput = True
arcpy.CheckOutExtension("Spatial")
os.makedirs(HASIL + r"\NDVI_ubin", exist_ok=True)

# A. Semua ubin dalam satu folder (band 1 = NIR, band 2 = Merah)
arcpy.env.workspace = DATA + r"\Ortofoto_Mentah"
daftar = arcpy.ListRasters("*", "TIF")
for i, nama in enumerate(daftar, start=1):
    b1 = Float(Raster(os.path.join(DATA, "Ortofoto_Mentah", nama, "Band_1")))
    b2 = Float(Raster(os.path.join(DATA, "Ortofoto_Mentah", nama, "Band_2")))
    ((b1 - b2) / (b1 + b2)).save(os.path.join(HASIL, "NDVI_ubin", "NDVI_" + nama))
    print("%d/%d selesai: %s" % (i, len(daftar), nama))

# B. Satu putaran per KPH: pilih satu poligon, potong CHM dengannya
chm = Raster(DATA + r"\DSM_Drone.tif") - Raster(DATA + r"\DTM_Drone.tif")
batas = DATA + r"\shp\Batas_KPH.shp"
with arcpy.da.SearchCursor(batas, ["NAMA_KPH"]) as k:
    nama_kph = [b[0] for b in k]
for nama in nama_kph:
    arcpy.analysis.Select(batas, r"memory\satu_kph", "NAMA_KPH = '%s'" % nama)
    keluar = os.path.join(HASIL, "CHM_" + nama.replace(" ", "_") + ".tif")
    ExtractByMask(chm, r"memory\satu_kph").save(keluar)
    print("CHM tersimpan untuk", nama)
arcpy.CheckInExtension("Spatial")

#Bagian C: ArcMap 10.8

#Bagian C: Batch di ArcMap 10.8

ArcMap memakai Python 2.7, jadi os.makedirs tidak punya exist_ok dan harus diperiksa lebih dulu dengan os.path.exists. ArcMap tidak punya padanan QgsTask. Esri menyediakan pengaturan pemrosesan latar belakang di opsi geoprocessing, yang berlaku untuk alat yang dijalankan dari antarmuka, bukan untuk skrip mandiri. [CEK]

PYTHON (ARCPY, PYTHON 2.7)
# -*- coding: utf-8 -*-
# [SKRIP 3.1-ArcMap: Batch: NDVI semua ubin dan potongan CHM per KPH (Python 2.7)]
# Penulis: Badar Mubarok Yogaswara
# [CEK] Cocokkan dengan versi Anda. Butuh ekstensi Spatial Analyst.
import os
import arcpy
from arcpy.sa import Raster, Float, ExtractByMask

DATA = r"D:\Latihan\paket-i4"
HASIL = r"D:\Latihan\hasil_i4"
arcpy.env.overwriteOutput = True
arcpy.CheckOutExtension("Spatial")
if not os.path.exists(HASIL + r"\NDVI_ubin"):          # Python 2.7 tidak punya exist_ok
    os.makedirs(HASIL + r"\NDVI_ubin")

# A. Semua ubin dalam satu folder (band 1 = NIR, band 2 = Merah)
arcpy.env.workspace = DATA + r"\Ortofoto_Mentah"
daftar = arcpy.ListRasters("*", "TIF")
for i, nama in enumerate(daftar, 1):
    b1 = Float(Raster(os.path.join(DATA, "Ortofoto_Mentah", nama, "Band_1")))
    b2 = Float(Raster(os.path.join(DATA, "Ortofoto_Mentah", nama, "Band_2")))
    ((b1 - b2) / (b1 + b2)).save(os.path.join(HASIL, "NDVI_ubin", "NDVI_" + nama))
    print("%d/%d selesai: %s" % (i, len(daftar), nama))

# B. Satu putaran per KPH: pilih satu poligon, potong CHM dengannya
chm = Raster(DATA + r"\DSM_Drone.tif") - Raster(DATA + r"\DTM_Drone.tif")
batas = DATA + r"\shp\Batas_KPH.shp"
with arcpy.da.SearchCursor(batas, ["NAMA_KPH"]) as k:
    nama_kph = [b[0] for b in k]
for nama in nama_kph:
    arcpy.Select_analysis(batas, r"in_memory\satu_kph", "\"NAMA_KPH\" = '%s'" % nama)
    keluar = os.path.join(HASIL, "CHM_" + nama.replace(" ", "_") + ".tif")
    ExtractByMask(chm, r"in_memory\satu_kph").save(keluar)
    print("CHM tersimpan untuk %s" % nama)
arcpy.CheckInExtension("Spatial")

#Cek paham

  1. Mengapa daftar berkas dibungkus sorted(...)?
  2. Di fungsi mana layer hasil boleh ditambahkan ke proyek pada skrip QgsTask?
  3. Apa gunanya task.isCanceled()?

Jawaban:

  1. Agar urutan pemrosesan dan hasilnya selalu sama setiap kali dijalankan.
  2. Di fungsi selesai, yang berjalan di program utama. Fungsi tugas berjalan di latar belakang dan tidak boleh menyentuh proyek.
  3. Untuk memeriksa apakah pengguna menekan Cancel. Bila ya, fungsi tugas berhenti dengan rapi.

#Kesalahan umum

  • Mengira urutan band. Rumus NDVI memakai band NIR dan Merah. Pada ubin ini band 1 = NIR dan band 2 = Merah. Periksa dulu pada data Anda.
  • Satu berkas rusak menghentikan seluruh batch. Ganti folder pada Skrip 3.1 menjadi Ortofoto_Campur. Tiga ubin pertama selesai, lalu ubin rusak menghentikan skrip dengan pesan Error creating output file. Pesan itu tidak menyebut berkas rusak. Bab 4 menjawab masalah ini.
  • Menyentuh proyek dari latar belakang. Pindahkan ke fungsi selesai.
  • Folder keluaran belum ada. Buat dengan os.makedirs.

#Ringkasan dan latihan

Ringkasan: perulangan for menjalankan resep yang sama untuk banyak berkas atau fitur. glob mencari berkas, dan sorted menertibkan urutan. QgsTask memindahkan pekerjaan berat ke latar belakang: fungsi tugas bekerja, fungsi selesai menampilkan hasil.

Latihan:

  1. Ubah Skrip 3.1 agar nama hasil berbunyi ndvi_01.tif, ndvi_02.tif, dan seterusnya (petunjuk: enumerate).
  2. Ubah Skrip 3.3 agar fungsi selesai mencetak rata-rata NDVI tiap hasil. Gunakan native:rasterlayerstatistics.

#Tabel perbandingan: batch

HalQGISArcGIS ProArcMap 10.8
Batch tanpa kodeExecute as batch processBatch pada dialog alat [CEK]Batch pada menu klik kanan alat [CEK]
Daftar rasterglob.glob("*.tif")arcpy.ListRastersarcpy.ListRasters
Perulangan atas fiturgetFeatures()arcpy.da.SearchCursorarcpy.da.SearchCursor
Tugas latar belakangQgsTask (diuji)Tidak ada padanan langsung; skrip batch dan alat skrip [CEK]Tidak ada padanan langsung; opsi latar belakang [CEK]
Bilah kemajuantask.setProgressarcpy.SetProgressor [CEK]Tidak ada padanan skrip [CEK]