import numpy as np
import pandas as pd
from scipy import integrate
def calculate_xyz(wavelengths, intensities, cie_data_path='cie1931_2deg.csv'):
"""
计算光谱数据的 CIE XYZ 三刺激值
参数:
wavelengths -- 波长数组 (nm)
intensities -- 对应光强数组
cie_data_path -- CIE 1931标准观察者数据文件路径
返回:
XYZ -- 包含 X, Y, Z 值的 numpy 数组
"""
# 数据预处理:归一化光强
intensities = intensities / np.max(intensities)
# 插值到标准波长范围 (380-780nm, 1nm间隔)
new_wavelengths = np.linspace(380, 780, 401)
new_intensities = np.interp(new_wavelengths, wavelengths, intensities)
# 加载CIE 1931标准观察者数据
cmf = pd.read_csv(cie_data_path)
cmf_wavelengths = cmf['wavelength'].values
# 插值匹配光谱数据
x_bar = np.interp(new_wavelengths, cmf_wavelengths, cmf['x_bar'])
y_bar = np.interp(new_wavelengths, cmf_wavelengths, cmf['y_bar'])
z_bar = np.interp(new_wavelengths, cmf_wavelengths, cmf['z_bar'])
# 计算XYZ三刺激值
X = integrate.simps(new_intensities * x_bar, new_wavelengths)
Y = integrate.simps(new_intensities * y_bar, new_wavelengths)
Z = integrate.simps(new_intensities * z_bar, new_wavelengths)
# 归一化常数 k = 100/∫ȳ(λ)dλ
k = 100 / integrate.simps(y_bar, new_wavelengths)
return np.array([X, Y, Z]) * k
# 使用示例
if __name__ == "__main__":
# 示例数据(实际应从文件读取)
wavelengths = np.linspace(400, 700, 31) # 400-700nm, 31个点
intensities = np.exp(-0.5*((wavelengths-550)/50)**2) # 高斯分布
# 计算XYZ值
XYZ = calculate_xyz(wavelengths, intensities, 'cie1931_2deg.csv')
print(f"CIE XYZ 值: X={XYZ[0]:.4f}, Y={XYZ[1]:.4f}, Z={XYZ[2]:.4f}")
对于这个代码没有1931的文件如何直接使用