Saya mencoba melakukan pencocokan histogram menggunakan Python untuk meningkatkan proses mosaik beberapa raster yang tumpang tindih. Saya mendasarkan kode saya pada yang ditemukan di:
http://www.idlcoyote.com/ip_tips/histomatch.html
Sampai saat ini, saya telah berhasil klip area yang tumpang tindih dari dua raster yang berdekatan dan meratakan array.
jadi saya punya dua array 1 dimensi dengan panjang yang sama.
Saya kemudian menulis kode berikut berdasarkan yang ditemukan di situs web di atas. Dalam kode yang ditunjukkan, saya telah mengganti dua dataset sangat kecil untuk gambar gd dan bd.
import matplotlib.pyplot as plt
from scipy.interpolate import interp1d
bins = range(0,100, 10)
gd_hist = [1,2,3,4,5,4,3,2,1]
bd_hist = [2,4,6,8,10,8,6,4,2]
nPixels = len(gd_hist)
# here we are creating the cumulative distribution frequency for the bad image
cdf_bd = []
for k in range(0, len(bins)-1):
b = sum(bd_hist[:k])
cdf_bd.append(float(b)/nPixels)
# here we are creating the cumulative distribution frequency for the good image
cdf_gd = []
for l in range(0, len(bins)-1):
g = sum(gd_hist[:l])
cdf_gd.append(float(g)/nPixels)
# we plot a histogram of the number of
plt.plot(bins[1:], gd_hist, 'g')
plt.plot(bins[1:], bd_hist, 'r--')
plt.show()
# we plot the cumulative distribution frequencies of both images
plt.plot(bins[1:], cdf_gd, 'g')
plt.plot(bins[1:], cdf_bd, 'r--')
plt.show()
z = []
# loop through the bins
for m in range(0, len(bins)-1):
p = [cdf_bd.index(b) for b in cdf_bd if b < cdf_gd[m]]
if len(p) == 0:
z.append(0)
else:
# if p is not empty, find the last value in the list p
lastval = p[len(p)-1]
# find the bin value at index 'lastval'
z.append(bins[lastval])
plt.plot(bins[1:], z, 'g')
plt.show()
# look into the 'bounds_error'
fi = interp1d(bins[1:], z, bounds_error=False, kind='cubic')
plt.plot(bins[1:], gd_hist, 'g')
plt.show
plt.plot(bins[1:], fi(bd_hist), 'r--')
plt.show()
Program saya memplot histogram dan distribusi frekuensi kumulatif dengan sukses ... dan saya pikir saya memiliki bagian untuk mendapatkan fungsi transformasi 'z' benar .... tetapi kemudian ketika saya menggunakan fungsi distribusi 'fi' pada 'bd_hist' untuk mencoba mencocokkannya dengan gd dataset semuanya berbentuk pear.
Saya bukan ahli matematika dan sangat mungkin saya telah mengabaikan sesuatu yang cukup jelas.
cdf_bd = np.cumsum(bd_hist) / float(np.sum(bd_hist))
Jawaban:
Meskipun saya tidak dapat mengomentari penerapan yang disarankan, Anda mungkin ingin memeriksa implementasi pencocokan histogram yang sudah ada yang dilakukan untuk GRASS GIS 7 (di sini tambahan):
https://trac.osgeo.org/grass/browser/grass-addons/grass7/imagery/i.histo.match
Untuk manual dan contoh, lihat
http://grass.osgeo.org/grass70/manuals/addons/i.histo.match.html
Kode diterbitkan di bawah lisensi GPL2 +.
sumber
Sebagai fudge liar; Saya tidak yakin Anda memerlukan PDF jika Anda memiliki data jumlah dalam kategori ...
Bisakah Anda mengonversi jumlah setiap nilai untuk setiap histogram yang berbeda menjadi nilai XY, dan kemudian menggunakan semacam indikator regresi untuk memeriksa kecocokan itu? Yaitu, untuk dua histogram yang sangat identik, analisis korelasi akan memberikan dan R kuadrat dari 1,0.
sumber
beberapa data sampel akan bagus karena dapat bervariasi dari sat ke sat. di sini adalah salah satu skrip sederhana yang saya buat dalam upaya untuk menyamakan histogram:
https://github.com/rupestre-campos/histogram_equalize
Mungkin Anda bisa mendapatkan wawasan.
Ini juga menghitung cdf seperti yang Anda lakukan, tetapi seperti yang saya coba itu akan menjadi gila jika Anda menghitung band-per-band, jadi Anda harus mempertimbangkan seluruh raster.
Sepertinya Anda kehilangan keseimbangan referensi warna dan profil spektral. Juga ada kebutuhan untuk tidak menghitung tidak ada piksel data, harus menghapus kemudian dari jumlah piksel gambar total untuk menghitung pdf yang benar.
Setelah beberapa pengujian saya menyukai hasil visual menggunakan pendekatan raster keseluruhan untuk 3-4 band Landsat8 dan mengkonversi dari 16bit ke 8bit rentang 0-255.
sumber