BAB 30: Empat Pekerjaan Klasik di ArcGIS: Skrip Masterclass dalam Bahasa ArcPy

Studi kasus: "Empat tugas rutin, satu kantor cabang" #

Kantor cabang memakai ArcGIS, dan Kepala Seksi punya daftar tugas rutin. Pertama, mengubah dua band citra satelit menjadi peta kehijauan (NDVI). Kedua, menghitung seberapa tepat peta klasifikasi dibanding titik survei lapangan. Ketiga, menghitung NDVI untuk banyak foto drone sekaligus. Keempat, mencetak satu PDF laporan kebakaran untuk setiap KPH. Anda sudah mengenal empat pekerjaan ini di versi QGIS. Sekarang Anda menuliskannya dalam bahasa ArcPy, di ArcGIS Pro dan di ArcMap.

Konsep: pekerjaan sama, bahasa beda #

Bayangkan satu resep nasi goreng yang ditulis dalam dua buku masak. Bahannya sama dan urutannya sama. Hanya nama alat dapurnya yang berbeda. Begitu pula skrip di bab ini. Setiap pekerjaan hadir dalam dua versi: ArcGIS Pro (Python 3) lebih dulu, lalu ArcMap 10.8 (Python 2.7).

Lima istilah baru yang perlu Anda kenal sekarang:

• Ekstensi (extension): izin tambahan di ArcGIS untuk kelompok alat tertentu. Alat raster di bab ini memerlukan ekstensi Spatial Analyst. Anda "meminjam" izin itu di awal skrip dan "mengembalikannya" di akhir.
• Objek raster: nilai di Python yang mewakili satu raster. Objek ini bisa dijumlah, dikurangi, dan dibagi seperti angka biasa. Cara kerja hitungan ini dinamakan map algebra.
• arcpy.env.workspace: folder kerja bawaan. Nama data tanpa jalur dicari di folder ini.
• Geodatabase (berkas .gdb): satu wadah untuk banyak layer. Anda memakainya di Skrip 30.2.
• Map Series (di Pro) dan Data Driven Pages (di ArcMap): fitur yang membuat satu halaman peta untuk setiap fitur di sebuah layer indeks. Anda memakainya di Skrip 30.4.

Semua keluaran di bab ini perkiraan. Penulis belum pernah menjalankan skrip ArcPy ini di ArcGIS. Angka di bagian "Hasil" dihitung dengan alat lain (GDAL dan Python) dari data sintetis di paket data, jadi tulisan persisnya di komputer Anda bisa sedikit berbeda. Semuanya bertanda [CEK].

Skrip 30.1: Memotong citra dan menghitung NDVI #

Berkas bab1_skrip_1_1_pro.py dan bab1_skrip_1_1_arcmap.py. Di dalam berkasnya tertulis "Skrip 1.1".

Tujuan. Memotong dua band citra dengan batas KPH, lalu menghitung NDVI dari kedua band itu dan menyimpannya.

Bahan. Tiga data di folder C:\data_kph: B04.tif (band merah), B08.tif (band inframerah dekat, NIR), dan Batas_KPH.shp. Di paket data buku ini, B04.tif dan B08.tif ada di paket-data\Sentinel_Lokal, sedangkan Batas_KPH.shp ada di paket-i4\shp (bersama berkas pendampingnya .dbf, .shx, dan .prj). Salin ketiganya ke satu folder. Ekstensi Spatial Analyst harus tersedia. Jalankan dari jendela Python, atau dari terminal Python milik ArcGIS (lihat Bab 29).

Versi ArcGIS Pro (Python 3) #

# [SKRIP 1.1-Pro: Memotong citra dan menghitung NDVI dengan ArcPy]
# Penulis: Badar Mubarok Yogaswara
# [CEK] Cocokkan dengan versi Anda. Butuh ekstensi Spatial Analyst atau Image Analyst.
import arcpy
from arcpy.sa import Raster, Float, ExtractByMask

arcpy.CheckOutExtension("Spatial")
arcpy.env.overwriteOutput = True
arcpy.env.workspace = r"C:\data_kph"      # ganti dengan folder Anda

# 1. Masukan: band Red, band NIR, dan batas area (CRS harus sama, lihat 1.2)
red = Raster("B04.tif")
nir = Raster("B08.tif")
batas = "Batas_KPH.shp"

# 2. Potong kedua band sesuai batas
red_clip = ExtractByMask(red, batas)
nir_clip = ExtractByMask(nir, batas)

# 3. NDVI = (NIR - RED) / (NIR + RED)
ndvi = (Float(nir_clip) - Float(red_clip)) / (Float(nir_clip) + Float(red_clip))
ndvi.save("NDVI_KPH.tif")

arcpy.CheckInExtension("Spatial")
print("Selesai: NDVI_KPH.tif tersimpan.")

Penjelasan baris demi baris. Nomor baris mengikuti berkas di atas.

Dokumentasi Esri menyebut ExtractByMask mengembalikan objek raster, dan objek itu disimpan dengan .save(). Karena itu baris 22 diperlukan. [CEK]

Versi ArcMap 10.8 (Python 2.7) #

# [SKRIP 1.1-ArcMap: Memotong citra dan menghitung NDVI dengan ArcPy]
# Penulis: Badar Mubarok Yogaswara
# [CEK] Cocokkan dengan versi Anda. Butuh ekstensi Spatial Analyst.
from __future__ import print_function   # agar print() berlaku di Python 2.7
import arcpy
from arcpy.sa import Raster, Float, ExtractByMask

arcpy.CheckOutExtension("Spatial")
arcpy.env.overwriteOutput = True
arcpy.env.workspace = r"C:\data_kph"      # ganti dengan folder Anda

red = Raster("B04.tif")
nir = Raster("B08.tif")
batas = "Batas_KPH.shp"

red_clip = ExtractByMask(red, batas)
nir_clip = ExtractByMask(nir, batas)

ndvi = (Float(nir_clip) - Float(red_clip)) / (Float(nir_clip) + Float(red_clip))
ndvi.save("NDVI_KPH.tif")

arcpy.CheckInExtension("Spatial")
print("Selesai: NDVI_KPH.tif tersimpan.")

Hasil. Layar menampilkan satu baris pesan. Berkas NDVI_KPH.tif muncul di C:\data_kph. Berikut perkiraannya, dihitung dengan GDAL dari data contoh. [CEK]

Selesai: NDVI_KPH.tif tersimpan.

Batas KPH di data contoh menutupi seluruh citra (400 m kali 400 m, 200 kali 200 piksel pada 2 m). Rata-rata NDVI seluruh piksel kira-kira 0,505. Kisaran nilai NDVI selalu antara -1 dan 1, dan semakin tinggi berarti semakin hijau.

Bagian yang boleh diubah. Nomor baris mengikuti versi Pro.

• Baris 9, folder kerja. Ganti C:\data_kph dengan folder data Anda. Pertahankan huruf r.
• Baris 12-14, nama berkas. Ganti B04.tif, B08.tif, dan Batas_KPH.shp dengan data Anda. Pastikan band merah dan NIR tidak tertukar.
• Baris 21, rumus. Rumus ini boleh diganti rumus indeks lain yang memakai dua band.
• Baris 22, nama keluaran. Ganti NDVI_KPH.tif untuk nama lain.

Jebakan.

• Sistem koordinat harus sama. Komentar di baris 11 mengingatkan hal ini. Bila citra dan batas berbeda sistem koordinat, hasil potongan bisa kosong atau salah tempat. Cek dulu sifat kedua data. [CEK]
• Band tertukar. Bila B04 dan B08 tertukar, NDVI berubah tanda tanpa pemberitahuan. Pastikan B04 merah dan B08 NIR.
• Ekstensi tidak dikembalikan saat galat. Baris 24 hanya tercapai bila baris-baris sebelumnya berhasil. Bila skrip berhenti di tengah, izin Spatial Analyst tetap dipinjam sampai Anda menutup sesi. Skrip I4 di Bab 31 memperbaikinya dengan try dan finally.
• Nama keluaran tanpa jalur. ndvi.save("NDVI_KPH.tif") mengandalkan folder kerja. Periksa bahwa berkas muncul di C:\data_kph. [CEK]
• Keluaran menimpa tanpa tanya. Karena overwriteOutput = True, berkas lama bernama sama hilang.

Latihan variasi. Tampilkan ukuran raster merah sebelum NDVI dihitung. Contoh jawaban, satu baris tambahan sebelum baris 21:

print("Ukuran:", red.width, "x", red.height, "piksel")

Dari data contoh, layar menampilkan Ukuran: 200 x 200 piksel, lalu pesan selesai. [CEK]

Skrip 30.2: Akurasi titik survei lawan peta klasifikasi #

Berkas bab2_skrip_2_1_pro.py dan bab2_skrip_2_1_arcmap.py. Di dalam berkasnya tertulis "Skrip 2.1".

Tujuan. Menarik kelas peta klasifikasi ke setiap titik survei, lalu menghitung berapa persen titik yang kondisinya cocok dengan kelas peta.

Bahan. Geodatabase C:\data_kph\Survei.gdb berisi dua layer: Titik_Survei (10 titik, kolom Kondisi) dan Poligon_Sawah (kolom Kelas_Citra). Di paket data buku ini, keduanya berbentuk GeoPackage di paket-data: Survei.gpkg (layer Titik_Survei) dan Poligon_Sawah.gpkg. Anda perlu memasukkan keduanya ke geodatabase buatan Anda dengan alat impor bawaan ArcGIS. Skrip tidak melakukan itu. [CEK] Skrip ini sepadan dengan skrip PyQGIS di Bab 24.

Versi ArcGIS Pro (Python 3) #

# [SKRIP 2.1-Pro: Validasi Titik Ground-Truth vs Poligon Klasifikasi dengan ArcPy]
# Penulis: Badar Mubarok Yogaswara
# [CEK] Cocokkan dengan versi Anda. Spatial Join tersedia pada semua level lisensi.
import arcpy

arcpy.env.workspace = r"C:\data_kph\Survei.gdb"      # ganti dengan lokasi Anda
arcpy.env.overwriteOutput = True

survei = "Titik_Survei"
klasifikasi = "Poligon_Sawah"
hasil = "Validasi_Hasil"

print("Memulai validasi lapangan otomatis...")

# 1. Spatial Join: tarik atribut poligon ke titik survei yang berada di dalamnya
arcpy.analysis.SpatialJoin(survei, klasifikasi, hasil,
                           "JOIN_ONE_TO_ONE", "KEEP_ALL", "", "INTERSECT")

# 2. Hitung akurasi: bandingkan Kondisi dengan Kelas_Citra
total_titik = 0
titik_cocok = 0
with arcpy.da.SearchCursor(hasil, ["Kondisi", "Kelas_Citra"]) as kursor:
    for kondisi, kelas in kursor:
        total_titik += 1
        if str(kondisi).lower() == str(kelas).lower():
            titik_cocok += 1

akurasi = (titik_cocok / total_titik) * 100 if total_titik > 0 else 0
print(f"Laporan Validasi: {titik_cocok} dari {total_titik} titik cocok.")
print(f"Akurasi Keseluruhan (Overall Accuracy): {akurasi:.2f}%")

Penjelasan baris demi baris. Nomor baris mengikuti berkas di atas.

Dokumentasi Esri mencatat urutan isian SpatialJoin: target_features, join_features, out_feature_class, lalu isian pilihan join_operation, join_type, field_mapping, match_option. Urutan itu cocok dengan baris 16-17. Setelah gabungan, alat menambahkan kolom Join_Count. [CEK]

Versi ArcMap 10.8 (Python 2.7) #

# [SKRIP 2.1-ArcMap: Validasi Titik Ground-Truth vs Poligon Klasifikasi dengan ArcPy]
# Penulis: Badar Mubarok Yogaswara
# [CEK] Cocokkan dengan versi Anda. Spatial Join tersedia pada semua level lisensi.
from __future__ import print_function   # agar print() berlaku di Python 2.7
import arcpy

arcpy.env.workspace = r"C:\data_kph\Survei.gdb"      # ganti dengan lokasi Anda
arcpy.env.overwriteOutput = True

survei = "Titik_Survei"
klasifikasi = "Poligon_Sawah"
hasil = "Validasi_Hasil"

print("Memulai validasi lapangan otomatis...")

# 1. Spatial Join: tarik atribut poligon ke titik survei yang berada di dalamnya
arcpy.SpatialJoin_analysis(survei, klasifikasi, hasil,
                           "JOIN_ONE_TO_ONE", "KEEP_ALL", "", "INTERSECT")

# 2. Hitung akurasi: bandingkan Kondisi dengan Kelas_Citra
total_titik = 0
titik_cocok = 0
with arcpy.da.SearchCursor(hasil, ["Kondisi", "Kelas_Citra"]) as kursor:
    for kondisi, kelas in kursor:
        total_titik += 1
        if str(kondisi).lower() == str(kelas).lower():
            titik_cocok += 1

akurasi = (float(titik_cocok) / total_titik) * 100 if total_titik > 0 else 0
print("Laporan Validasi: {} dari {} titik cocok.".format(titik_cocok, total_titik))
print("Akurasi Keseluruhan (Overall Accuracy): {:.2f}%".format(akurasi))

Hasil. Layar menampilkan tiga baris. Dari data contoh, 9 dari 10 titik cocok (satu titik sengaja dibuat tidak cocok). Perkiraan ini dihitung dengan GDAL dan Python dari Survei.gpkg dan Poligon_Sawah.gpkg. [CEK]

Memulai validasi lapangan otomatis...
Laporan Validasi: 9 dari 10 titik cocok.
Akurasi Keseluruhan (Overall Accuracy): 90.00%

Layer Validasi_Hasil juga muncul di dalam Survei.gdb. Logika Python skrip ini sudah diuji dengan modul ArcPy tiruan buatan sendiri yang membaca paket data, dan hasilnya cocok dengan perkiraan di atas. Perilaku alat ArcGIS yang sebenarnya belum teruji. [CEK]

Bagian yang boleh diubah. Nomor baris mengikuti versi Pro.

• Baris 6, geodatabase. Ganti dengan geodatabase milik Anda.
• Baris 9-11, nama layer. Ganti Titik_Survei, Poligon_Sawah, dan Validasi_Hasil sesuai data Anda.
• Baris 22, nama kolom. Ganti Kondisi dan Kelas_Citra dengan dua kolom yang ingin Anda bandingkan.
• Baris 17, pengaturan gabungan. Ganti "INTERSECT" dengan cara pencocokan lain yang didukung alat. Daftarnya ada di dokumentasi alat. [CEK]

Jebakan.

• Titik di luar poligon. Karena memakai KEEP_ALL, titik tanpa pasangan tetap ada di hasil, dengan Kelas_Citra kosong. str(None).lower() menghasilkan teks none, yang tidak sama dengan kondisi mana pun. Titik itu dihitung tidak cocok. Hasil str(None).lower() sudah diuji di Python 3.
• Nama kolom kembar. Bila kedua layer punya kolom bernama sama, nama di hasil gabungan bisa berubah, dan baris 22 gagal menemukan kolomnya. [CEK]
• Pembagian di ArcMap tanpa float. Di Python 3, (9 / 10) * 100 menghasilkan 90.0. Pembagian bilangan bulat gaya Python 2.7 (diperagakan dengan (9 // 10) * 100 di Python 3) menghasilkan 0.
• Geodatabase belum ada. Skrip tidak membuat Survei.gdb dan tidak mengimpor data. Bila geodatabase kosong, kursor gagal.
• Hasil gabungan menimpa layer lama. overwriteOutput = True membuat Validasi_Hasil yang lama diganti.

Latihan variasi. Tambahkan satu baris laporan: jumlah titik yang tidak cocok. Contoh jawaban, tepat sebelum baris cetak akurasi (baris 30):

print(f"Tidak cocok: {total_titik - titik_cocok} titik.")

Dari data contoh, layar menampilkan Tidak cocok: 1 titik. Untuk ArcMap, ganti dengan print("Tidak cocok: {} titik.".format(total_titik - titik_cocok)). [CEK]

Skrip 30.3: NDVI untuk semua ubin foto drone #

Berkas bab3_skrip_3_1_pro.py dan bab3_skrip_3_1_arcmap.py. Di dalam berkasnya tertulis "Skrip 3.1".

Tujuan. Menghitung NDVI untuk setiap ubin foto drone dalam satu folder, lalu menyimpan hasilnya di folder lain, dengan laporan kemajuan.

Bahan. Folder D:\Data_Drone\Ortofoto_Mentah berisi ubin .tif dua band (band 1 inframerah dekat, band 2 merah). Di paket data, tiga ubin ada di paket-i4\Ortofoto_Mentah. Siapkan juga folder kosong D:\Data_Drone\NDVI_Hasil: skrip tidak membuat folder ini. Ekstensi Spatial Analyst harus tersedia. Padanan PyQGIS-nya ada di Bab 24.

Versi ArcGIS Pro (Python 3) #

# [SKRIP 3.1-Pro: Batch Processing Citra Drone dengan ArcPy]
# Penulis: Badar Mubarok Yogaswara
# [CEK] Cocokkan dengan versi Anda. Butuh ekstensi Spatial Analyst atau Image Analyst.
import os
import arcpy
from arcpy.sa import Raster, Float

arcpy.CheckOutExtension("Spatial")
arcpy.env.overwriteOutput = True

folder_input = r"D:\Data_Drone\Ortofoto_Mentah"
folder_output = r"D:\Data_Drone\NDVI_Hasil"
arcpy.env.workspace = folder_input

daftar = arcpy.ListRasters("*", "TIF")
total = len(daftar)
print(f"Memulai batch {total} ubin...")

for i, nama in enumerate(daftar, start=1):
    # Urutan band mengikuti Skrip 3.1 QGIS: pastikan Band_1 = NIR dan Band_2 = Red
    band1 = Float(Raster(os.path.join(folder_input, nama, "Band_1")))
    band2 = Float(Raster(os.path.join(folder_input, nama, "Band_2")))
    ndvi = (band1 - band2) / (band1 + band2)
    ndvi.save(os.path.join(folder_output, "NDVI_" + nama))
    print(f"{i}/{total} selesai: {nama}")

arcpy.CheckInExtension("Spatial")
print("Batch selesai.")

Penjelasan baris demi baris. Nomor baris mengikuti berkas di atas. Hanya bagian baru yang dijelaskan. Baris import, ekstensi, dan overwriteOutput sudah Anda kenal di Skrip 30.1.

Dokumentasi Esri menyebut ListRasters memerlukan folder kerja yang sudah diatur lebih dulu, dan filter namanya tidak peka huruf besar-kecil. Karena itu baris 13 harus mendahului baris 15. [CEK]

Versi ArcMap 10.8 (Python 2.7) #

# [SKRIP 3.1-ArcMap: Batch Processing Citra Drone dengan ArcPy]
# Penulis: Badar Mubarok Yogaswara
# [CEK] Cocokkan dengan versi Anda. Butuh ekstensi Spatial Analyst.
from __future__ import print_function   # agar print() berlaku di Python 2.7
import os
import arcpy
from arcpy.sa import Raster, Float

arcpy.CheckOutExtension("Spatial")
arcpy.env.overwriteOutput = True

folder_input = r"D:\Data_Drone\Ortofoto_Mentah"
folder_output = r"D:\Data_Drone\NDVI_Hasil"
arcpy.env.workspace = folder_input

daftar = arcpy.ListRasters("*", "TIF")
total = len(daftar)
print("Memulai batch {} ubin...".format(total))

for i, nama in enumerate(daftar, 1):
    # Urutan band mengikuti Skrip 3.1 QGIS: pastikan Band_1 = NIR dan Band_2 = Red
    band1 = Float(Raster(os.path.join(folder_input, nama, "Band_1")))
    band2 = Float(Raster(os.path.join(folder_input, nama, "Band_2")))
    ndvi = (band1 - band2) / (band1 + band2)
    ndvi.save(os.path.join(folder_output, "NDVI_" + nama))
    print("{}/{} selesai: {}".format(i, total, nama))

arcpy.CheckInExtension("Spatial")
print("Batch selesai.")

Hasil. Layar menampilkan lima baris, dan tiga berkas NDVI_ubin_01.tif sampai NDVI_ubin_03.tif muncul di D:\Data_Drone\NDVI_Hasil. Berikut perkiraannya. Nama ubin dan urutannya mengikuti folder contoh. [CEK]

Memulai batch 3 ubin...
1/3 selesai: ubin_01.tif
2/3 selesai: ubin_02.tif
3/3 selesai: ubin_03.tif
Batch selesai.

Rata-rata NDVI tiap ubin, dihitung dengan GDAL dari data contoh, adalah kira-kira 0,471 untuk ubin_01, 0,471 untuk ubin_02, dan 0,467 untuk ubin_03. [CEK] Logika Python skrip sudah diuji dengan modul ArcPy tiruan buatan sendiri yang membaca paket data. Tiga berkas keluaran terbentuk dan cetakan kemajuannya sama dengan perkiraan di atas. Perilaku alat ArcGIS yang sebenarnya belum teruji.

Bagian yang boleh diubah. Nomor baris mengikuti versi Pro.

• Baris 11-12, dua folder. Ganti dengan folder Anda. Pertahankan huruf r.
• Baris 15, pola nama dan format. Ganti "TIF" untuk format lain, atau ganti "*" dengan pola nama, misalnya "ubin_0*".
• Baris 21-22, nomor band. Tukar Band_1 dan Band_2 bila urutan band data Anda berbeda.
• Baris 24, awalan nama. Ganti "NDVI_" untuk awalan lain.

Jebakan.

• Folder keluaran tidak dibuat. Bila D:\Data_Drone\NDVI_Hasil belum ada, simpan di baris 24 kemungkinan besar gagal. Bab 31 menunjukkan cara membuat folder lewat os.makedirs. [CEK]
• Band tertukar. Bila band 1 sebenarnya merah, NDVI berubah tanda dan tidak ada peringatan. Cek urutan band dulu.
• Hanya .tif. Berkas berakhiran lain tidak ikut diproses karena filter "TIF". [CEK]
• Berjalan di depan. Skrip ini menahan ArcGIS sampai selesai. Untuk puluhan ubin, jalankan dari terminal Python di luar aplikasi [CEK]. Versi QGIS-nya memakai tugas latar belakang (Bab 24).
• Ekstensi tidak dikembalikan saat galat. Sama dengan Skrip 30.1.

Latihan variasi. Proses hanya dua ubin pertama untuk uji coba. Contoh jawaban, ganti baris 15 dengan:

daftar = arcpy.ListRasters("*", "TIF")[:2]

Tanda [:2] mengambil dua isi pertama daftar. Nilai total otomatis menjadi 2. Layar menampilkan Memulai batch 2 ubin..., 1/2 selesai: ubin_01.tif, 2/2 selesai: ubin_02.tif, lalu Batch selesai. Jawaban ini sudah diuji dengan modul ArcPy tiruan; perilaku ArcGIS yang sebenarnya belum teruji. [CEK]

Skrip 30.4: Satu PDF laporan untuk setiap KPH #

Berkas bab6_skrip_6_1_pro.py, bab6_skrip_6_1_arcmap.py, dan teks pendamping bab6_teks_dinamis.txt. Di dalam berkasnya tertulis "Skrip 6.1".

Tujuan. Mengekspor satu PDF per KPH dari tata letak peta yang halamannya dibuat satu per satu untuk tiap KPH.

Bahan. Skrip ini tidak bisa dijalankan dari terminal biasa. Ia memakai proyek (atau dokumen peta) yang sedang terbuka, jadi jalankan dari jendela Python di dalam aplikasi. Di Pro, proyek Anda harus punya tata letak bernama persis Laporan Kebakaran Hutan dengan Map Series aktif, memakai layer indeks yang punya kolom NAMA_KPH (misalnya Batas_KPH, 3 fitur). Di ArcMap, dokumen peta Anda harus punya Data Driven Pages aktif dengan layer indeks yang sama. Siapkan folder D:\Laporan. Skrip tidak membuatnya. Padanan PyQGIS-nya (Atlas) ada di Bab 24.

Versi ArcGIS Pro (Python 3) #

# [SKRIP 6.1-Pro: Ekspor PDF per KPH dari Map Series]
# Penulis: Badar Mubarok Yogaswara
# [CEK] Cocokkan dengan versi Anda.
import os
import arcpy

aprx = arcpy.mp.ArcGISProject("CURRENT")
lyt = aprx.listLayouts("Laporan Kebakaran Hutan")[0]
ms = lyt.mapSeries
folder_out = r"D:\Laporan"

for no in range(1, ms.pageCount + 1):
    ms.currentPageNumber = no
    nama = ms.pageRow.NAMA_KPH
    lyt.exportToPDF(os.path.join(folder_out, f"{nama}_Laporan_Bulan_Ini.pdf"))
    print(f"Selesai: {nama}")

Penjelasan baris demi baris. Nomor baris mengikuti berkas di atas.

Dokumentasi Esri mencatat pageCount bersifat baca saja, currentPageNumber boleh dibaca dan diisi, dan pageRow memberi atribut fitur indeks halaman yang aktif. Dokumentasi yang sama mencatat Layout.exportToPDF(out_pdf, ...) dan ArcGISProject.listLayouts({wildcard}). Semua itu cocok dengan skrip. [CEK]

Versi ArcMap 10.8 (Python 2.7) #

# [SKRIP 6.1-ArcMap: Ekspor PDF per KPH dari Data Driven Pages]
# Penulis: Badar Mubarok Yogaswara
# [CEK] Cocokkan dengan versi Anda.
from __future__ import print_function   # agar print() berlaku di Python 2.7
import os
import arcpy

mxd = arcpy.mapping.MapDocument("CURRENT")
ddp = mxd.dataDrivenPages
folder_out = r"D:\Laporan"

for no in range(1, ddp.pageCount + 1):
    ddp.currentPageID = no
    nama = ddp.pageRow.getValue("NAMA_KPH")
    arcpy.mapping.ExportToPDF(mxd, os.path.join(folder_out, nama + "_Laporan_Bulan_Ini.pdf"))
    print("Selesai: {}".format(nama))

del mxd

Dokumentasi ArcMap menyebut DataDrivenPages diakses lewat MapDocument, dengan sifat pageCount, currentPageID, dan pageRow. Fungsi ExportToPDF(map_document, out_pdf, ...) juga tercatat di sana. Cara membaca nilai lewat getValue belum penulis cocokkan dengan dokumentasi. [CEK]

Teks dinamis untuk tata letak ArcMap #

Berkas bab6_teks_dinamis.txt bukan skrip Python. Isinya empat baris teks untuk elemen teks di tata letak ArcMap. Tanda <dyn .../> diganti isi kolom milik halaman yang tampil.

LAPORAN KEBAKARAN HUTAN KPH <dyn type="page" property="attribute" field="NAMA_KPH"/>
-----------------------------------
Total Luas KPH: <dyn type="page" property="attribute" field="Luas_Ha"/> Ha
Jumlah Titik Api Aktif: <dyn type="page" property="attribute" field="Join_Count"/> Titik.

Tiap tanda <dyn type="page" property="attribute" field="..."/> mengambil isi satu kolom dari layer indeks. Tiga kolom yang dipakai: NAMA_KPH, Luas_Ha, dan Join_Count. Kolom Luas_Ha dan Join_Count harus sudah ada di layer indeks. Skrip di bab ini tidak membuatnya. Di data contoh, luas ketiga KPH adalah 4, 4, dan 8 hektare, dan jumlah titik panas masing-masing 2, 5, dan 5. [CEK] Di Pro, teks dinamis serupa tersedia, tetapi berkas pendampingnya tidak ada di folder skrip. [CEK]

Hasil. Layar menampilkan satu baris untuk tiap halaman, dan PDF muncul di D:\Laporan. Dari data contoh, ada tiga halaman. Berikut perkiraannya. Urutan halaman mengikuti urutan fitur di layer indeks. [CEK]

Selesai: KPH Alpha
Selesai: KPH Beta
Selesai: KPH Gamma

Tiga berkas PDF terbentuk: KPH Alpha_Laporan_Bulan_Ini.pdf, KPH Beta_Laporan_Bulan_Ini.pdf, dan KPH Gamma_Laporan_Bulan_Ini.pdf. Nama berkas memuat spasi karena nama KPH memuat spasi. Di versi ArcMap, cetakannya sama persis, sebab print("Selesai: {}".format(nama)) menghasilkan teks yang sama.

Bagian yang boleh diubah. Nomor baris mengikuti versi Pro.

• Baris 8, nama tata letak. Ganti dengan nama tata letak Anda, persis sama termasuk huruf besar dan spasi.
• Baris 10, folder keluaran. Ganti D:\Laporan dan pastikan folder itu ada.
• Baris 14, nama kolom. Ganti NAMA_KPH dengan kolom lain di layer indeks.
• Baris 15, nama berkas. Ganti _Laporan_Bulan_Ini.pdf untuk akhiran lain.

Jebakan.

• Nama tata letak harus persis. Bila tidak ada tata letak dengan nama itu, listLayouts(...) menghasilkan daftar kosong dan [0] gagal dengan IndexError: list index out of range. Pesan itu sudah diuji pada daftar kosong di Python 3. [CEK]
• Map Series harus aktif. Dokumentasi Esri menyebut lyt.mapSeries bernilai None bila Map Series tidak ada atau tidak aktif. Maka baris 12 gagal dengan AttributeError: 'NoneType' object has no attribute 'pageCount'. Pesan Python itu sudah diuji pada nilai None.
• Folder D:\Laporan harus ada. Skrip tidak membuatnya. [CEK]
• Dari terminal biasa gagal. Kata "CURRENT" hanya berarti bila ada proyek atau dokumen peta yang terbuka.
• Nama berkas bisa tertimpa. Dua KPH dengan nama sama menghasilkan satu PDF.

Latihan variasi. Beri nomor urut di depan nama berkas, supaya PDF terurut rapi. Contoh jawaban, ganti baris 15 dengan:

    lyt.exportToPDF(os.path.join(folder_out, f"{no:02d}_{nama}.pdf"))

Bagian :02d memaksa nomor memakai dua angka dengan nol di depan. Nama berkasnya menjadi 01_KPH Alpha.pdf, 02_KPH Beta.pdf, dan 03_KPH Gamma.pdf. Format :02d sudah diuji di Python 3, tetapi ekspor PDF-nya belum teruji. [CEK]

Cek paham #

1. Mengapa skrip Pro dan skrip ArcMap pada Skrip 30.2 menghitung persen akurasi dengan cara sedikit berbeda?
2. Apa yang terjadi bila folder D:\Data_Drone\NDVI_Hasil belum ada saat Skrip 30.3 dijalankan?
3. Mengapa Skrip 30.4 tidak bisa dijalankan dari terminal biasa?

Jawaban.

1. Di Python 2.7 (ArcMap), pembagian dua bilangan bulat menghasilkan bilangan bulat, sehingga 9 dibagi 10 menjadi 0. Versi ArcMap memakai float(titik_cocok) agar hasilnya desimal.
2. Skrip tidak membuat folder itu. Simpan hasil di baris 24 kemungkinan besar gagal. Buat foldernya lebih dulu. [CEK]
3. Karena memakai "CURRENT", yaitu proyek atau dokumen peta yang sedang terbuka. Di terminal biasa tidak ada yang terbuka.

Kesalahan umum #

• Lupa meminjam ekstensi. Alat raster seperti ExtractByMask memerlukan Spatial Analyst. Perbaikan: cek baris arcpy.CheckOutExtension("Spatial") ada di awal.
• Menyalin gaya Python 3 ke ArcMap. Contohnya f-string. Perbaikan: pakai .format atau % di ArcMap.
• Memakai jalur tanpa huruf r. Perbaikan: tulis r"C:\data_kph".
• Menjalankan skrip ber-CURRENT di luar aplikasi. Perbaikan: jalankan dari jendela Python.

Ringkasan dan latihan #

• Empat pekerjaan klasik dituliskan dalam ArcPy: NDVI dari dua band, akurasi titik survei, batch NDVI, dan PDF per KPH.
• Pola raster: pinjam ekstensi, buka band, hitung dengan Float, simpan, kembalikan ekstensi.
• Pola vektor: SpatialJoin, lalu kursor SearchCursor per baris.
• Beda ArcMap yang utama: from __future__ import print_function, float(...) untuk pembagian, .format, nama alat bergaya _analysis, dan arcpy.mapping.
• Semua skrip di bab ini belum dijalankan di ArcGIS dan bertanda [CEK].

Latihan. Ubah Skrip 30.1 supaya memotong dan menghitung NDVI untuk data milik Anda sendiri. Bagian mana yang Anda ganti, dan bagian mana yang tetap?

Tabel perbandingan: PyQGIS dan ArcPy pada empat pekerjaan #