import openpyxl
import numpy as np
import matplotlib.pyplot as plt
import cv2
from sklearn.linear_model import LinearRegression
import sys
def read_excel(filename: str, range_str: str, sheet_no: int = 0) -> np.array:
"""Excelから数値データ取得
Args:
filename (str): ファイル名
range_str (str): セル範囲文字列 'A1:B2'
Returns:
np.array: セル値配列
"""
wb = openpyxl.load_workbook(filename)
sheet = wb[wb.sheetnames[sheet_no]]
return np.array(
[
[cell.value for cell in sheet[range_str][i]]
for i in range(len(sheet[range_str]))
]
)
def extend_rgb(rgb_arr: np.array) -> np.array:
"""RGB値の高次項を拡張する
Args:
rgb_arr (np.array): RGB値の配列 1x3, Nx3
Returns:
np.array: 拡張した配列 1x11, Nx11
"""
conv_2d = 0
if rgb_arr.ndim == 1:
rgb_arr = rgb_arr.reshape([1, len(rgb_arr)])
conv_2d = 1
extended = np.array(
[
[
1,
rgb[0],
rgb[1],
rgb[2],
rgb[0] ** 2,
rgb[1] ** 2,
rgb[2] ** 2,
rgb[0] * rgb[1],
rgb[1] * rgb[2],
rgb[2] * rgb[0],
rgb[0] * rgb[1] * rgb[2],
]
for rgb in rgb_arr
]
)
if conv_2d == 1:
return extended[0]
else:
return extended
def lab_func(t: float) -> float:
"""XYZ->LAB変換に使用する暗部補正付き関数
Args:
t (float): 値
Returns:
float: 変換結果
"""
if t > 0.008856:
return t ** (1.0 / 3.0)
else:
return 7.787 * t + 16.0 / 116.0
def xyz_to_lab(xyz_arr: np.array, white_xyz: np.array) -> np.array:
"""XYZからCIELABへ変換
Args:
xyz_arr (np.array): XYZ三刺激値の配列 1x3, Nx3
white_xyz (np.array): 白色点の三刺激値
Returns:
np.array: L*a*b*の配列 1x3, Nx3
"""
conv_2d = 0
if xyz_arr.ndim == 1:
xyz_arr = xyz_arr.reshape([1, len(xyz_arr)])
conv_2d = 1
lab_arr = np.array(
[
[
116.0 * lab_func(xyz[1] / white_xyz[1]) - 16.0,
500.0
* (lab_func(xyz[0] / white_xyz[0]) - lab_func(xyz[1] / white_xyz[1])),
200.0
* (lab_func(xyz[1] / white_xyz[1]) - lab_func(xyz[2] / white_xyz[2])),
]
for xyz in xyz_arr
]
)
if conv_2d == 1:
return lab_arr[0]
else:
return lab_arr
if __name__ == "__main__":
# --------------------------------------------------------------
# Excelからデータ取得
light_sp = read_excel("Spectrum.xlsx", "B2:AF2")[0]
cmf_sp_arr = read_excel("Spectrum.xlsx", "B3:AF5")
filter_sp_arr = read_excel("Spectrum.xlsx", "B6:AF8")
tongue_sp_arr = read_excel("Spectrum.xlsx", "B9:AF16")
# tcc_name = "handy_tcc"
# tcc_name = "takano_tcc"
tcc_name = "sugawara_tcc_mat8485"
tcc_sp_arr = read_excel(tcc_name + ".xlsx", "F2:AJ25") / 100.0
white_patch = (
0 if tcc_name == "handy_tcc" else 12
) # ホワイトバランスに用いる色票番号
# --------------------------------------------------------------
# 白色点の計算
white_xyz = np.array([light_sp.dot(cmf_sp) for cmf_sp in cmf_sp_arr])
k = 100.0 / white_xyz[1]
white_xyz *= k
print("白色点 XYZ=", white_xyz)
# --------------------------------------------------------------
# TCCのXYZ, L*a*b* 算出(正解値)
tcc_xyz_arr = np.array(
[
[
(light_sp * tcc_sp_arr[p]).dot(cmf_sp_arr[i]) * k
for i in range(len(cmf_sp_arr))
]
for p in range(len(tcc_sp_arr))
]
)
# print("tcc_xyz=")
# print(tcc_xyz)
tcc_lab_arr = xyz_to_lab(tcc_xyz_arr, white_xyz)
print("TCC L*a*b* 正解値=")
print(tcc_lab_arr)
# --------------------------------------------------------------
# カメラシミュレータのホワイトバランス
white_balance = np.array([0.9, 0.9, 0.9]) # ホワイトバランスのRGB値
raw_rgb = np.array(
[
(light_sp * tcc_sp_arr[white_patch]).dot(filter_sp_arr[i])
for i in range(len(filter_sp_arr))
]
)
# print("raw_rgb=", raw_rgb)
rgb_gain = white_balance / raw_rgb
# print("rgb_gain=", rgb_gain)
# --------------------------------------------------------------
# カメラシミュレータによる色票のRGB値算出
tcc_rgb_arr = np.array(
[
[
(light_sp * tcc_sp_arr[p]).dot(filter_sp_arr[i]) * rgb_gain[i]
for i in range(len(filter_sp_arr))
]
for p in range(len(tcc_sp_arr))
]
)
# print("TCC RGB=")
# print(tcc_rgb_arr)
# print(extend_rgb(tcc_rgb))
# --------------------------------------------------------------
# カメラシミュレータRGBからXYZ推定行列の算出
# tcc_rgb_input_arr = tcc_rgb_arr
tcc_rgb_input_arr = extend_rgb(tcc_rgb_arr)
conv_rgb_to_xyz = LinearRegression()
conv_rgb_to_xyz.fit(tcc_rgb_input_arr, tcc_xyz_arr)
# --------------------------------------------------------------
# カメラシミュレータRGB値からXYZ, L*a*b*の推定
tcc_xyz_pred_arr = conv_rgb_to_xyz.predict(tcc_rgb_input_arr)
# print("TCC XYZ 推定値=")
# print(tcc_xyz_pred_arr)
tcc_lab_pred_arr = xyz_to_lab(tcc_xyz_pred_arr, white_xyz)
print("TCC L*a*b* 推定値=")
print(tcc_lab_pred_arr)
# --------------------------------------------------------------
# TCC 正解値と推定値の色差計算
delta_e_arr = np.array(
[
np.linalg.norm(tcc_lab_pred_arr[i] - tcc_lab_arr[i], ord=2)
for i in range(len(tcc_lab_arr))
]
)
delta_e = np.mean(delta_e_arr)
print("TCC 推定色差 DeltaE=", delta_e)
# --------------------------------------------------------------
# 舌のXYZ, L*a*b* 算出(正解値)
tongue_xyz_arr = np.array(
[
[
(light_sp * tongue_sp_arr[p]).dot(cmf_sp_arr[i]) * k
for i in range(len(cmf_sp_arr))
]
for p in range(len(tongue_sp_arr))
]
)
# print("tongue_xyz=")
# print(tongue_xyz_arr)
tongue_lab_arr = xyz_to_lab(tongue_xyz_arr, white_xyz)
print("舌 L*a*b* 正解値=")
print(tongue_lab_arr)
# --------------------------------------------------------------
# カメラシミュレータによる舌のRGB値算出
tongue_rgb_arr = np.array(
[
[
(light_sp * tongue_sp_arr[p]).dot(filter_sp_arr[i]) * rgb_gain[i]
for i in range(len(filter_sp_arr))
]
for p in range(len(tongue_sp_arr))
]
)
print("舌 RGB=")
print(tongue_rgb_arr)
# --------------------------------------------------------------
# 舌のカメラシミュレータRGB値からXYZ, L*a*b*の推定
tongue_rgb_input_arr = extend_rgb(tongue_rgb_arr)
tongue_xyz_pred_arr = conv_rgb_to_xyz.predict(tongue_rgb_input_arr)
print("舌 XYZ 推定値=")
print(tongue_xyz_pred_arr)
tongue_lab_pred_arr = xyz_to_lab(tongue_xyz_pred_arr, white_xyz)
print("舌 L*a*b* 推定値=")
print(tongue_lab_pred_arr)
# --------------------------------------------------------------
# 舌 正解値と推定値の色差計算
delta_e_arr = np.array(
[
np.linalg.norm(tongue_lab_pred_arr[i] - tongue_lab_arr[i], ord=2)
for i in range(len(tongue_lab_arr))
]
)
delta_e = np.mean(delta_e_arr)
print("舌 推定色差 DeltaE=", delta_e)
fig = plt.figure(figsize=(15, 8))
ax_ab = fig.add_subplot(1, 2, 1, xlabel="a*", ylabel="b*")
plt.grid(True)
ax_ab.set_title(tcc_name + " a* - b*")
ax_ab.set_xlim(-40, 60)
ax_ab.set_ylim(-40, 60)
ax_ab.text(-35, 55, "red:GT, black:Estim, green:tcc")
for i in range(len(tongue_lab_arr)):
ax_ab.plot(
[tongue_lab_arr[i, 1], tongue_lab_pred_arr[i, 1]],
[tongue_lab_arr[i, 2], tongue_lab_pred_arr[i, 2]],
"-",
c="gray",
)
ax_ab.plot(tongue_lab_arr[:, 1], tongue_lab_arr[:, 2], "o", c="r")
ax_ab.plot(tongue_lab_pred_arr[:, 1], tongue_lab_pred_arr[:, 2], "o", c="k")
ax_ab.plot(tcc_lab_arr[:, 1], tcc_lab_arr[:, 2], "o", c="g")
ax_aL = fig.add_subplot(1, 2, 2, xlabel="a*", ylabel="L*")
plt.grid(True)
ax_aL.set_title(tcc_name + " a* - L*")
ax_aL.set_xlim(-40, 60)
ax_aL.set_ylim(0, 100)
ax_aL.text(-35, 95, "red:GT, black:Estim, green:tcc")
for i in range(len(tongue_lab_arr)):
ax_aL.plot(
[tongue_lab_arr[i, 1], tongue_lab_pred_arr[i, 1]],
[tongue_lab_arr[i, 0], tongue_lab_pred_arr[i, 0]],
"-",
c="gray",
)
ax_aL.plot(tongue_lab_arr[:, 1], tongue_lab_arr[:, 0], "o", c="r")
ax_aL.plot(tongue_lab_pred_arr[:, 1], tongue_lab_pred_arr[:, 0], "o", c="k")
ax_aL.plot(tcc_lab_arr[:, 1], tcc_lab_arr[:, 0], "o", c="g")
plt.show()