Skip to content

01_relinear_jpg

复制本地路径 | 在线编辑

原始文件为 .py 代码，本文是转换后的 Markdown 文件。

```.py
'''
将 JPG 文件线性化，涉及到计算相机响应曲线，OpenCV 也有现成的 API，本文件也进行了计算

本文件的一些计算方式和 OpenCV 不太一样（因为完全参考本文件夹的实验手册，里面有些地方要求和 OpenCV 不一样，比如这里会比较不同的权重函数等等）
在 relinear_as_opencv.py 中复现了和 OpenCV 一样的相机响应曲线
'''

先读取文件、设置曝光时间

import os
os.chdir(os.path.dirname(os.path.abspath(file)))

import cv2
import time
import math
import numpy as np
import matplotlib.pyplot as plt
import scipy.optimize

picnum = 16
images = [cv2.imread(f'../data/door_stack/exposure{i}.jpg', -1) for i in range(picnum)]
exptime = np.array([(2**i)/4096 for i in range(picnum)], dtype=np.float32)

获取图片宽和高

rows, cols = images[0].shape[0:2]

'''
直接用 opencv 的 API 来做
'''

参考: https://learnopencv.com/high-dynamic-range-hdr-imaging-using-opencv-cpp-python/

校准（对齐）各张图片，不过有的时候结果可能会反而不好

cv2.createAlignMTB().process(images, images)

不仅仅只有 Debevec 方法，用的时候再去细看

获取响应函数 -> 制作 HDR

cv_response = cv2.createCalibrateDebevec().process(images, exptime)
result = cv2.createMergeDebevec().process(images, exptime, cv_response)

print(cv_response.shape)

plt.plot(np.log(cv_response[..., 0]))

plt.plot(np.log(cv_response[..., 1]))

plt.plot(np.log(cv_response[..., 2]))

plt.title('camera response by opencv')

plt.show()

ToneMapping，有很多，不列举了，看上面的网页，解释的非常非常好！

tonemapping_1 = cv2.createTonemapReinhard(1.5, 0, 0, 0).process(result)

cv2.imwrite('./opencv_result.png', tonemapping_1*255)

'''
权重函数，最小化要用到
'''
def w_uniform(z, minv=0.05, maxv=0.95):
return 1 if minv<z<maxv else 0

def w_tent(z, minv=0.05, maxv=0.95):
return min(z, 1-z) if minv<z<maxv else 0

def w_gaussian(z, minv=0.05, maxv=0.95):
return np.exp(-4(z-0.5)(z-0.5)/0.25) if minv<z<maxv else 0

def w_photon(z, minv=0.05, maxv=0.95):
return 1 if minv<z<maxv else 0

wfuncs = {
'uniform': w_uniform,
'tent': w_tent,
'gaussian': w_gaussian,
'photon': w_photon
}

'''
自己来求响应，两种方法第一种：最小二乘法矩阵形式
'''
def compute_by_matrix(I, my_lambda, w_func):
N, P = I.shape[0:2]

A = np.zeros((N*P+255, 256+N))
b = np.zeros((N*P+255, 1))

for i in range(N):
    for j in range(P):
        Iij = int(I[i, j])
        w = w_func(Iij / 255)
        A[i*P+j, Iij] = w
        A[i*P+j, 256+i] = -w
        b[i*P+j] = w * math.log2(exptime[j])

for i in range(254):
    A[N*P+i, i+0] = my_lambda * w_func((i+1)/255)
    A[N*P+i, i+1] = -2 * my_lambda * w_func((i+1)/255)
    A[N*P+i, i+2] = my_lambda * w_func((i+1)/255)
A[-1, 128] = 1

x = np.linalg.lstsq(A, b, rcond=None)[0]
return x[0:256, 0]

'''
自己来求响应，两种方法第二种：用 scipy minimize 来迭代求解
'''
def compute_by_scipy(I, my_lambda, w_func):
def compute_func(data):
# func: 前256个数字是 g(i)，后面是 logL
g, logL = data[0:256], data[256:]

    result = 0
    # 第一项
    for ij in range(N):
        for k in range(P):
            pixel = int(I[ij, k])
            w = w_func(pixel / 255)
            result += ( w*(g[pixel] - logL[ij] - math.log2(exptime[k])) )**2
    # 第二项
    for z in range(1, 255):
        Dg = g[z+1] - 2*g[z] + g[z-1]
        w = w_func(z / 255)
        result += my_lambda * (( w * Dg )**2)

    return result

N, P = I.shape[0:2]
# 限制：g(127) = 0
cons= [{'type':'eq', 'fun': lambda m: m[127]}]

# # 限制: g 是单调递增函数
# for i in range(1, 256):
#     # create lambda in for loop 需要注意！！Cite: stackoverflow 2295290
#     cons.append({'type':'ineq', 'fun': lambda m,i=i: m[i] - m[i-1] - 1e-8})

# 迭代次数
opts = {'maxiter':1000, 'disp':False}

# 初始化
this_input = np.zeros(256+N)
res = scipy.optimize.minimize(compute_func, this_input, method='SLSQP', constraints=cons, options=opts)
return res.x[0:256]

在 P 张曝光图像上选择 N 个点

P = picnum
N = 60

I 表示选择的点，注意 I(i,j) 中 i 表示像素位置，j 表示第几张图片

I = np.zeros((N, P))

要求对三个通道像素都进行处理

方式：1.直接平均；2.转成 Y 分量；3.每个分量都取一些点

这里用第三个方式

for channel in range(0, 3):
# 每个通道平均取 N/3 个点
assert(N % 3 == 0)

nowN = N // 3
points = np.random.randint([0, 0], [rows, cols], size=(nowN, 2))
for pidx in range(P):
    onepic = images[pidx][..., channel][points[:, 0], points[:, 1]]
    I[nowN*channel:nowN*(channel+1), pidx] = onepic

# 计算结果，并且比较结果和 opencv 的标定结果

for my_lambda in [10, 20, 30]:

for name, wfunc in wfuncs.items():

# 注意：我们是三个通道合起来处理得出一个响应函数 g；而 OpenCV 是三个通道分别处理

plt.figure()

plt.title(f'{name}-lambda={my_lambda}')

plt.plot(np.log(cv_response[..., 1]), label='Opencv-G')

t1 = time.time()

matrix_response = compute_by_matrix(I, my_lambda=my_lambda, w_func=wfunc)

t2 = time.time()

scipy_response = compute_by_scipy(I, my_lambda=my_lambda, w_func=wfunc)

t3 = time.time()

plt.plot(matrix_response, label=f'My-bymatrix-time={round(t2-t1, 3)}')

plt.plot(scipy_response, label=f'My-byscipy-time={round(t3-t2, 3)}')

plt.legend()

plt.show()

存储结果

good_lambda = {
'uniform': 25,
'tent': 25,
'gaussian': 25,
'photon': 25
}
os.makedirs('../data/response', exist_ok=True)
for name, wfunc in wfuncs.items():
print(name, 'doing...............')
matrix_response = compute_by_matrix(I, my_lambda=good_lambda[name], w_func=wfunc)
np.save(f'../data/response/{name}', np.exp(matrix_response))```

Comments