Saya menggunakan python dan QGIS 2.0. Saya mencoba untuk klip raster di folder dengan satu fitur poligon. Ini pertama kalinya bagi saya menggunakan (katakanlah) "PyQGIS", saya sudah terbiasa dengan arcpy sebelumnya. Lagi pula, saya tidak membuat skrip sederhana saya berfungsi, saran apa pun akan sangat dihargai!
import qgis.core, qgis,utils
QgsApplication.setPrefixPath("C:/OSGeo4W64/apps/qgis", True)
QgsApplication.initQgis()
CLIP= "C:/Users/unim/Documents/Umberto/Universita/PhD/Guglielmin/Permafrost/Alta_Valtellina/Landsat_ita/study_area_foscagno.shp"
INPUT_FOLDER="C:/Users/unimi/Documents/Umberto/Universita/PhD/Guglielmin/Permafrost/Alta_Valtellina/Landsat_ita/LE71930282000259EDC00"
OUTPUT= "C:/Users/unim/Documents/Umberto/Universita/PhD/Guglielmin/Permafrost/Alta_Valtellina/Landsat_ita/foscagno_pyqgis/"
for RASTER in INPUT_FOLDER.tif
do
echo "Processing $RASTER"
gdalwarp -q -cutline CLIP -crop_to_cutline -of GTiff RASTER OUTPUT+ "clip_"+ RASTER
done
QgsApplication.exitQgis()
Di bawah ini adalah peningkatan yang telah saya lakukan sejak sekarang, tidak membuat skrip bekerja, tapi saya pikir saya mungkin semakin dekat ...
import qgis.core, qgis.utils, os, fnmatch
from osgeo import gdal
CLIP= "C:/Users/unimi/Documents/Umberto/Universita/PhD/Guglielmin/Permafrost/Alta_Valtellina/Landsat_ita/study_area_foscagno.shp"
INPUT_FOLDER= "C:/Users/unimi/Documents/Umberto/Universita/PhD/Guglielmin/Permafrost/Alta_Valtellina/Landsat_ita/LE71930282000259EDC00/DNs2Reflectance_LE71930282000259EDC00"
OUTPUT= "C:/Users/unimi/Documents/Umberto/Universita/PhD/Guglielmin/Permafrost/Alta_Valtellina/Landsat_ita/Cloud_mask_AltaValtellina/clip_2_foscagno"
def findRasters (path, filter):
for root, dirs, files in os.walk(path):
for file in fnmatch.filter(files, filter):
yield os.path.join (root, file)
for raster in findRasters (INPUT_FOLDER, '*.tif'):
print (raster)
outRaster = OUTPUT + '/clip_' + raster
cmd = 'gdalwarp -dstnodata 0 -q -cutline CLIP -crop_to_cutline %s %s' % (raster, outRaster)
os.system (cmd)
Saya pikir mungkin ada sesuatu yang salah dalam perintah "gdal", karena fungsi "print" melakukan tugasnya dengan benar, tetapi tidak ada file yang ditulis ke output, saya juga tidak menerima kesalahan. Ngomong-ngomong, sulit untuk menyukai dokumentasi yang mudah tentang pengkodean gdal ...
CLIP
dalamcmd
ekspresi masalahnya. Jika Anda meletakkan variabel dalam string, itu tidak dibaca. Sebagai gantinya, Anda akan menggabungkan string dengan variabel.print(cmd)
di tempatos.system(cmd)
.outRaster
Variabel Anda salah.Jawaban:
Saya setuju dengan Nathan. Anda perlu membuat python seluruh skrip Anda. Jadi gantilah
for
loop Anda dengan sesuatu seperti berikut ini:Catatan 1: Saya kira file raster Anda adalah GeoTIFF (
*.tif
).Catatan 2:
-of GTiff
tidak diperlukan dicmd
, karena ini adalah format keluaran default digdalwarp
.sumber
os.system
, Anda benar.Saya akhirnya berhasil dengan skrip yang sangat sederhana dan bersih ini, yang memanggil GDAL dari Python tanpa mengimpornya (seperti yang disarankan, tetapi menggunakan metode "Call ()" alih-alih "os.system ()". Saya harap ini bisa membantu!
sumber
Versi modifikasi dari solusi umbe1987 untuk pengguna Linux:
sumber