basalt/scripts/check_gamma_correction.py

98 lines
2.6 KiB
Python
Raw Normal View History

2019-06-13 13:37:17 +02:00
#!/usr/bin/env python3
import sys
import math
import os
2019-07-11 18:49:13 +02:00
import cv2
2019-06-13 13:37:17 +02:00
import numpy as np
from matplotlib import pyplot as plt
dataset_path = sys.argv[1]
print(dataset_path)
timestamps = np.loadtxt(dataset_path + '/mav0/cam0/data.csv', usecols=[0], delimiter=',', dtype=np.int64)
2019-07-11 19:13:03 +02:00
exposures = np.loadtxt(dataset_path + '/mav0/cam0/exposure.csv', usecols=[1], delimiter=',', dtype=np.int64).astype(np.float64) * 1e-6
2019-06-13 13:37:17 +02:00
pixel_avgs = list()
2019-07-11 18:49:13 +02:00
if timestamps.shape[0] != exposures.shape[0]: print("timestamps and exposures do not match")
imgs = []
2019-06-13 13:37:17 +02:00
# check image data.
for timestamp in timestamps:
path = dataset_path + '/mav0/cam0/data/' + str(timestamp)
2019-07-11 18:49:13 +02:00
img = cv2.imread(dataset_path + '/mav0/cam0/data/' + str(timestamp) + '.webp', cv2.IMREAD_GRAYSCALE)[:,:,0]
imgs.append(img)
2019-06-13 13:37:17 +02:00
pixel_avgs.append(np.mean(img))
2019-07-11 18:49:13 +02:00
imgs = np.array(imgs)
print(imgs.shape)
print(imgs.dtype)
2019-07-12 10:56:47 +02:00
num_pixels_by_intensity = np.bincount(imgs.flat)
print('num_pixels_by_intensity', num_pixels_by_intensity)
inv_resp = np.arange(num_pixels_by_intensity.shape[0], dtype=np.float64)
inv_resp[-1] = -1.0 # Use negative numbers to detect saturation
2019-07-11 18:49:13 +02:00
def opt_irradiance():
corrected_imgs = inv_resp[imgs] * exposures[:, np.newaxis, np.newaxis]
times = np.ones_like(corrected_imgs) * (exposures**2)[:, np.newaxis, np.newaxis]
2019-07-11 19:13:03 +02:00
times[corrected_imgs < 0] = 0
corrected_imgs[corrected_imgs < 0] = 0
denom = np.sum(times, axis=0)
2019-07-12 10:56:47 +02:00
idx = (denom != 0)
irr = np.sum(corrected_imgs, axis=0)
irr[idx] /= denom[idx]
2019-07-11 19:13:03 +02:00
irr[denom == 0] = -1.0
2019-07-11 18:49:13 +02:00
return irr
def opt_inv_resp():
generated_imgs = irradiance[np.newaxis, :, :] * exposures[:, np.newaxis, np.newaxis]
2019-07-11 19:13:03 +02:00
num_pixels_by_intensity = np.bincount(imgs.flat, generated_imgs.flat >= 0)
2019-07-12 10:56:47 +02:00
generated_imgs[generated_imgs < 0] = 0
2019-07-11 18:49:13 +02:00
sum_by_intensity = np.bincount(imgs.flat, generated_imgs.flat)
new_inv_resp = inv_resp
idx = np.nonzero(num_pixels_by_intensity > 0)
new_inv_resp[idx] = sum_by_intensity[idx] / num_pixels_by_intensity[idx]
2019-07-12 10:56:47 +02:00
new_inv_resp[-1] = -1.0 # Use negative numbers to detect saturation
return new_inv_resp
2019-07-11 18:49:13 +02:00
def print_error():
generated_imgs = irradiance[np.newaxis, :, :] * exposures[:, np.newaxis, np.newaxis]
generated_imgs -= inv_resp[imgs]
generated_imgs[imgs == 255] = 0
print(np.sum(generated_imgs**2))
2019-07-12 10:56:47 +02:00
for iter in range(5):
2019-07-11 18:49:13 +02:00
irradiance = opt_irradiance()
print_error()
inv_resp = opt_inv_resp()
print_error()
2019-07-11 19:13:03 +02:00
2019-07-11 18:49:13 +02:00
plt.figure()
2019-07-12 10:56:47 +02:00
plt.plot(inv_resp[:-1])
plt.ylabel('Irradiance Value')
plt.xlabel('Image Intensity')
plt.title('Inverse Responce Function')
2019-07-11 18:49:13 +02:00
plt.figure()
plt.imshow(irradiance)
2019-07-12 10:56:47 +02:00
plt.title('Irradiance Image')
2019-06-13 13:37:17 +02:00
plt.show()