Saya telah ditugaskan untuk membuat analisis kesesuaian kondisi gelombang di Teluk Meksiko. Saya memiliki 2 ribu file raster yang masing-masing sekitar 8 MB (2438 kolom, 1749 baris, ukuran sel 1 km). Parameter dalam file raster adalah periode gelombang dan saya ingin mengklasifikasi ulang semua raster sehingga jika 4<= wave period <=9 then make cell = 1
, sel lain = 0. Kemudian jumlah semua raster menjadi raster terakhir dan bagi dengan jumlah total raster untuk memberikan persentase total pengamatan yang sesuai dan hasil ekspor ke dalam beberapa format yang kompatibel dengan ESRI ... mungkin sesuatu yang dapat mendukung pelampung jika perlu. Saya tidak bekerja sama sekali dengan Python atau R, tetapi setelah mencari online sepertinya masuk akal untuk melakukan proses ini dalam salah satu bahasa tersebut. Sejauh ini saya telah menemukan beberapa kode di R, tetapi saya bingung bagaimana cara membuatnya.
library(rgdal)
raster_data <- list.files(path=getwd()) #promt user for dir containing raster files
num_files <- length(raster_data)
for (i in raster_data) { #read in rasters
my_data <- readGDAL(raster_data[i])
Pada titik ini saya bingung apakah saya juga harus mereklasifikasi dan mulai menjumlahkan data dalam loop ini atau tidak. Dugaan saya adalah ya karena kalau tidak saya pikir saya akan kehabisan memori di komputer saya, tetapi tidak yakin. Saya juga tidak yakin tentang cara mereklasifikasi data.
Dalam riset online saya pikir saya akan menggunakan reclass(my_data, c(-Inf,3,0, 4,9,1, 10,Inf,0))
, tetapi apakah itu terlihat benar?
Dan untuk menjumlahkan akan saya gunakan sum(stack(my_data))
dan entah bagaimana output itu. Juga ... jika ini lebih efisien dilakukan atau ditulis dalam Python saya akan terbuka untuk itu juga. Saya benar-benar pemula dalam hal pemrograman.
Jawaban:
Ini adalah cara ringkas untuk melakukannya di R --- di sini tanpa file perantara:
sumber
R
harus disimpan dalam RAM untuk dibuats
. (Kemungkinan RAM dua kali lebih banyak diperlukan karenas
mungkin tidak akan digantiraster_data
.) Saya harap Anda memiliki banyak RAM! Solusi Anda dapat dibuat praktis dengan memecah 2000 kisi menjadi grup yang lebih kecil, melakukan perhitungan untuk setiap grup, dan kemudian menggabungkan perhitungan tersebut.stack
dancalc
bekerja seperti kebanyakanR
fungsi lainnya dengan terlebih dahulu memuat semua data ke dalam RAM.Seperti yang saya perhatikan dalam komentar, umumnya Anda harus menghindari penggunaan R untuk tujuan non-statistik karena masalah kinerja dalam aspek-aspek tertentu (bekerja dengan siklus adalah contoh). Berikut adalah contoh kode untuk Anda di Pyhton (terima kasih untuk artikel ini ) untuk klasifikasi ulang satu file dengan satu band. Anda akan dapat memodifikasinya dengan mudah untuk pemrosesan batch jika Anda tahu cara mendapatkan semua file dari direktori . Perhatikan bahwa raster direpresentasikan sebagai array, jadi Anda dapat menggunakan metode array untuk meningkatkan kinerja bila berlaku. Untuk bekerja dengan array di Python lihat dokumentasi Numpy .
UPD: kode yang saya posting pada awalnya adalah versi terpotong dari filter khusus yang diperlukan per pemrosesan piksel. Tetapi untuk pertanyaan ini, penggunaan Numpy akan meningkatkan kinerja (lihat kode).
sumber
Jangan gunakan readGDAL. Bunyinya menjadi objek Spasial * yang mungkin bukan ide yang baik ..
Gunakan
raster
paket. Itu dapat membaca hal-hal GDAL menjadi objek Raster. Ini adalah hal yang baik.r = raster("/path/to/rasterfile.tif")
akan membacanyar
.Klasifikasi Anda kemudian
t = r > 4 & r <= 9
Pertanyaan besar adalah apakah akan menampilkan ini ke file raster baru dan kemudian melakukan langkah ringkasan di loop lain. Jika Anda memiliki penyimpanan, saya akan menulisnya ke file baru hanya karena jika loop Anda gagal (karena salah satu dari 2000 file tersebut adalah sampah) Anda harus mulai lagi. Jadi gunakan
writeRaster
untuk membuat raster yang di-threshold jika Anda memutuskan untuk melakukannya.Maka loop Anda hanyalah sesuatu seperti
ini.
Manajemen memori R mungkin mengenai Anda di sini - saat Anda melakukannya,
count=count+r
R mungkin akan membuat salinan barucount
. Jadi itu berpotensi 2.000 salinannya - tetapi pengumpulan sampah akan masuk dan membersihkan, dan ini di sini bahwa perulangan R dulu sangat buruk.Dalam hal pengaturan waktu, pada PC 4yo saya, op ambang membutuhkan sekitar 1,5 detik pada raster sebesar angka acak antara nol dan dua puluh. Times 2000 = Anda mengerjakan matematika. Seperti biasa, buat set tes kecil untuk mengembangkan kode, lalu luncurkan pada data Anda yang sebenarnya dan nikmati kopi.
sumber
R
2.15.0), membaca grid 1,6 megapiksel (dalam format floating point ESRI asli) membutuhkan waktu 0,19 detik dan melakukan lima operasi loop - dua perbandingan, konjungsi, tambahan, dan toko - mengambil yang lain 0,09 detik, untuk 0,28 detik per iterasi. Ketika sebaliknya loop dilakukan dengan caching 100 grid sekaligus ke dalam array dan menggunakanrowSums
untuk melakukan total, penggunaan RAM tentu saja naik: tetapi, sama jelasnya, tidak ada pengumpulan sampah. Waktu turun hanya ke 0,26 detik per iterasi.