光学学报, 2020, 40 (5): 0520002, 网络出版: 2020-03-10   

基于非对称三角形仪器函数的最优加权表方法 下载: 945次

Optimum Weighting Table Method Based on Asymmetric Triangular Instrumental Function
作者单位
辽宁科技大学计算机与软件工程学院, 辽宁 鞍山114051
摘要
研究了仪器函数形状稍微偏离对称三角形情形下的三刺激值计算问题,探索了其对应的最优加权表。结果发现,最优权重可以通过求解三个系数矩阵的线性方程组来得到,且系数矩阵是对称正定三对角矩阵。将本课题组前期研究获得的最优加权表方法推广到仪器函数形状为非对称三角形的情形。仿真结果表明,对于10 nm和20 nm测量间隔,最优加权表方法的精度优于三点和五点校正方法;对于5 nm测量间隔,最好的计算方法是先对测量数据进行三点校正,再使用直接选取方法计算三刺激值。
Abstract
The problem of the calculation of the tristimulus value when the shape of the instrument function is slightly deviated from the symmetrical triangle is studied, and the corresponding optimal weighting table is explored. The results show that the optimal weight can be obtained by solving the linear equations of three coefficient matrices, and the coefficient matrices are symmetric positive definite tridiagonal matrices. The optimal weighting table method obtained in the previous study is extended to the case where the instrument function shape is asymmetric triangle. The simulation results show that the accuracy of the optimal weighting table method is better than that of the three-point and five-point calibration methods for the 10-nm and 20-nm measurement intervals. For the measurement interval of 5 nm, the best calculation method is to conduct three-point correction of the measurement data first, and then use the direct selection method to calculate the tri-stimulus value.

1 引言

三刺激值XYZ是色度学中最基本的空间,所有其他颜色空间都是经由XYZ空间转换而得到的,如国际照明委员会(CIE)的均匀颜色空间CIELAB,CIE色貌模型空间和基于色貌的均匀颜色空间CAM16-UCS[1]等。从理论上说,计算XYZ所涉及的被积函数都是波长(在可见光范围)的连续函数,但在实际应用中并不知道这些函数的确切解析表达式,因此如何准确计算三刺激值是色度学领域关注的焦点之一。一直以来,很多学者都致力于研究如何提高三刺激值计算的精度。常用的方法有ASTM E308-1985加权表法[2]、ASTM E308-1995加权表法[3]、ASTM E308-2018加权表法[4]、Oleari法[5-6]、最小二乘法[7]、最优权重法[8]以及通过校正测量反射比求取三刺激值的方法等。这些算法都是建立在仪器函数形状是面积为1的对称三角形的假设之上,但是工业应用中测色仪的仪器函数形状可能与理想假设稍微有偏差,从而基于理想的对称三角形形状的仪器函数导出的算法将会失效。

国际照明委员会(CIE)成立了一个技术委员会(TC2-60)[9-13],研究校正测量仪器函数的方法,同时给出了基于测量仪器函数的三点和五点校正法,但是相应的加权表方法鲜有报道。最优加权表方法是把反射比测量误差的校正都并入到加权表权重的构建中来,因此在计算三刺激值时直接使用测量数据会使计算非常简便。该方法尤其适用于在同一观察条件下计算多个物体色三刺激值,但依赖于仪器函数的形状。本文针对真实的仪器函数形状为稍微偏离对称三角形或非对称形状的情形,探究与其适配的最优权重法以计算三刺激值。

2 基于非对称三角形仪器函数的最优权重

国际照明委员会定义的反射物体颜色的三刺激值[14]

V=kabR(λ)E(λ)v-(λ)dλk=100/abE(λ)y-(λ)dλ,(1)

式中:λ为波长; V∈{X,Y,Z},X,Y,Z为三刺激值;E(λ)为标准照明体的相对光谱功率;R(λ)为光谱反射比; v-(λ)∈{ x-(λ), y-(λ), z-(λ)}, x-(λ), y-(λ), z-(λ)是CIE 1931或1964配色函数;(a,b)为可见光波长的范围;k为归一化因子。

由于被积函数并没有确切的解析表达式,在实际计算三刺激值时,CIE推荐采用1 nm波长间隔的函数值求和代替积分[14-16],即

V=kλ=a,step=1bR(λ)E(λ)v-(λ)k=100/λ=abE(λ)y-(λ),(2)

式中:step为步长。

尽管如此,(2)式的应用仍存在问题,原因是实际测量的反射比多数是10 nm或20 nm间隔。另外,记 R~i(Δλ)为反射比R(λ)在波长λi处的测量数据,且测量间隔为Δλ,其中i为测量点编号。则测量数据 R~i(Δλ)与真实反射比满足如下关系:

R~i(Δλ)=abPi(λ)R(λ)dλ,(3)

式中:Pi(λ)为测量仪器函数,满足

abPi(λ)dλ=1(4)

本文假定测量仪器函数Pi(λ)的表达式为

Pi(λ)=λ-λi-1sskew(Δλ)2,λi-1λ<λaλi+1-λ(2-sskew)·(Δλ)2,λaλλi+10,λ<λi-1orλ>λi+1,(5)

式中:λa为仪器函数在区间[λi-1,λi+1]内的一个点,仪器函数在该点达到最大值;sskew为非对称因子。

λa=λi-1+sskew·Δλ,(6)λi=a+i·Δλ,i=0,1,n,(7)n=(b-a)/Δλ(8)

通过非对称因子sskew并结合(6)式,控制仪器函数最高点所对应的横轴位置,进而控制仪器函数形状的倾斜程度。当sskew=1时,仪器函数Pi(λ)变为对称的仪器函数,如图1“+”线所示;当sskew偏离1时,Pi(λ)形状就变成了非对称的三角形形状,图1中的“·”线对应sskew=0.9,“*”线对应sskew=1.1的情况。在理想的情况下[4,16],sskew为1,但由于各种因素,如元器件质量、杂散光、装配过程、搬运过程等,真实的仪器函数形状会稍微偏离对称性。对于经常校正和维护的测色仪,一般sskew不会偏离1太远,故将sskew的范围定义为0.9~1.1,以0.01为步长,研究非对称仪器函数的细微变化。

图 1. 带宽为10 nm时Pi(λ)的模拟结果

Fig. 1. Simulated result of Pi(λ) when bandwidth is 10 nm

下载图片 查看所有图片

R~i(Δλ)为测量的反射比,最优权重表法[15,17-18]就是确定权重WV,i以计算三刺激值。

V'=i=0nWV,iR~(Δλ),(9)

式中:V'为准确三刺激值V的计算值。由(1)式可知,确定权重的原则是使V'V尽可能接近。由(1)式和(5)式可知,权重WV,i依赖于标准照明体、配色函数和仪器函数。

V'-V=i=0nWV,iR~i(Δλ)-abWV(λ)R(λ)dλ=abi=0nWV,iPi(λ)-WV(λ)R(λ)dλ=abDV(λ)R(λ)dλ,(10)

式中:

DV(λ)=i=0nWV,iPi(λ)-WV(λ),WV(λ)=k·E(λ)·v-(λ)(11)

理想情况下可通过选择WV,i使DV(λ)≡0,但在实际应用中,由于

(V'-V)2=abDV(λ)R(λ)dλ2abDV2(λ)dλabR2(λ)dλ,(12)

R2(λ)≤1,故

(V'-V)2(b-a)abDV2(λ)dλ(13)

F=F(WV,0(Δλ),WV,1(Δλ),,WV,i(Δλ),,WV,n(Δλ))=abD2(λ)dλ,(14)

因此,通过确定权重使V'V尽可能接近的问题就变成了求(14)式的极小值问题。由微积分原理可知,最优权重应满足 FWV,i=0或

FWV,i=2j=0nWV,jλi-1λi+1Pi(λ)Pj(λ)dλ-λi-1λi+1WV(λ)Pi(λ)dλ=0(15)

若令

ai,j=6Δλλi-1λi+1Pi(λ)Pj(λ)dλbi(V)=6Δλλi-1λi+1WV(λ)Pi(λ)dλ,(16)

则最优权重应满足

j=0nai,jWV,j=bi(V),i=0,1,,n(17)

由(16)式可知,当|i-j|≥2时,ai,j=0。对于aj,j,结合(5)式,有

aj,j=6Δλλj-1λj+1Pj(λ)Pj(λ)dλ=6Δλλj-1λa(λ-λj-1)2sskew2·(Δλ)4dλ+λaλj+1(λj+1-λ)2(2-sskew)2·(Δλ)4dλ=6Δλsskew2(Δλ)4·sskew3(Δλ)33+6Δλ(2-sskew)2(Δλ)4·(2-sskew)3(Δλ)33=2[sskew+(2-sskew)]=4(18)

图2图3分别给出了sskew<1和sskew>1时仪器函数的分布情况。可以看出,Pj-1(λ)与Pj(λ)相交的部分与横轴构成一个三角形, Pj(λ)和Pj+1(λ)相交的部分与横轴构成一个三角形,且这两个三角形的面积是一样的,因此有aj,j-1=aj-1,j

图 2. sskew<1时带宽函数的分布

Fig. 2. Distribution of bandpass function when sskew<1

下载图片 查看所有图片

图 3. sskew>1时带宽函数的分布

Fig. 3. Distribution of bandpass function when sskew>1

下载图片 查看所有图片

sskew∈(0.9,1.0)时,由图2可知:

aj,j-1=aj-1,j=6Δλλj-1λj+1Pj-1(λ)Pj(λ)dλ=6Δλλj-1λaλ-λj-1sskew·(Δλ)2·λj-λ(2-sskew)·(Δλ)2dλ+λaλjλj+1-λ(2-sskew)·(Δλ)2·λj-λ(2-sskew)·(Δλ)2dλ=6sskew2-sskew-sskew3+12+6(1-sskew)2(2-sskew)21-sskew3+12=2sskew2-6sskew+54-4sskew+sskew2(19)

sskew∈(1.0,1.1)时,由图3可知:

aj,j-1=aj-1,j=6Δλλj-1λj+1Pj-1(λ)Pj(λ)dλ=6Δλλj-1λaλ-λj-1sskew·(Δλ)2·λ-λj-2sskew·(Δλ)2dλ+λaλjλ-λj-1sskew·(Δλ)2·λj-λ(2-sskew)·(Δλ)2dλ=6(sskew-1)2sskew2sskew-13+12+6(2-sskew)sskew-2-sskew3+12=2sskew2-2sskew+1sskew2(20)

由(10)式可知,当R(λ)=c时,有

V'-V=i=0nWV,iR~i-abWV(λ)R(λ)dλ=abi=0nWV,iPi(λ)-WV(λ)R(λ)dλ=cabDV(λ)dλ=ci=0nWV,iPi(λ)-abWV(λ)dλ,(21)

式中:c为0和1之间的常数。

因此,若新的权重满足

j=0nWV,j=abWV(λ)dλ,(22)

V'=V。注意积分 abWV(λ)dλ是1 nm加权表的求和,并对应标准照明体的三刺激值。

由(17)式可知

i=0nbi(V)=i=0nj=0nai,jWV,j=j=0ni=0nai,jWV,j=(a0,0+a1,0)WV,0+j=1n-1(aj-1,j+aj,j+aj+1,j)WV,j+(an-1,n+an,n)WV,n(23)

由(19)、(20)式可知,aj,j-1, aj-1,j都与下标无关,因此aj-1,j+aj,j+aj+1,j是一个常数。令

d=aj-1,j+aj,j+aj+1,j,(24)

并调整a0,0an,n

a0,0=d-a0,1,an,n=d-an-1,n,(25)

则由(23)~(25)式有

i=0nbi(V)=dj=0nWV,j(26)

由(17)、(25)式求出的权重WV,j仍不满足(22)式,因此需要进一步调整,为此引入常数tV,即

tV=abWV(λ)dλ/i=0nbi(V),(27)

(17)式调整为

AuV=aV,(28)

式中:

A=gff4ff4ffg,(29)uV=WV,1WV,2WV,n,aV=d·tV·b0(V)b1(V)bn(V),(30)f=aj-1,j=aj+1,j,g=d-f(31)

注意矩阵An+1阶方阵,由(28)式得到的权重满足(22)式。由于系数矩阵A是对称的三对角矩阵,且是对称正定的,故方程组会有唯一解。另外,参数d·tVsskew发生微小变化,当sskew=1时乘积d·tV为1,而当sskew偏离1时,d·tV会稍微大于1。此外,还需要说明的是,当sskew=1时,(28)式正是文献[ 17]中对应的方程组(10)。

3 实验仿真

3.1 实验仿真准备

实验数据:1096个pantone色卡样本的1 nm间隔的反射比,测量波长范围为360~830 nm。通过这些标准的1 nm数据,可以采用(3)式进行模拟测量。

标准照明体:分别采用D65、A、D50、F02、F07、F11这六种标准照明体。根据CIE标准获得1 nm间隔、360~830 nm范围的六种标准照明体的相对光谱功率分布E(λ)。

标准观察者:采用10°视角的CIE1964标准色度观察者和2°视角的CIE1931标准色度观察者,根据CIE标准得到1 nm间隔、360~830 nm范围的相应标准观察者的颜色匹配函数 x-(λ), y-(λ), z-(λ)。

色差公式:选择最新的CIEDE2000色差公式[16]

统计度量:由于涉及大量数据,本文采用平均(Average)、最大(Max)和中间(Median)色差评价算法,每个指标越小越好。

3.2 模拟仿真

仿真过程如图4所示。

图 4. 仿真流程图

Fig. 4. Flow chart of simulation

下载图片 查看所有图片

首先,根据(2)式,使用标准反射比数据R(λ)、六种标准照明体的E(λ)可计算得到12组标准三刺激值XstdYstdZstd,其中每组有1096个色差。再根据(3)式,利用1 nm的真实反射比数据R(λ)计算指定间隔下380~780 nm范围内的测量数据 R~(Δλ)。测量间隔分别采用10 nm和20 nm。对每种间隔的仿真测量数据,分别选定标准照明体和配色函数,采用如下步骤。

1) 采用三点校正法和五点校正法校正 R~(Δλ),得到反射比 R~DO3(Δλ)R~DO5(Δλ)

2) 对校正后的反射比,采用3阶多项式插值方法进行插值获得1 nm间隔的数据,再结合1 nm间隔的标准照明体信息和配色函数计算三刺激值XDO3YDO3ZDO3XDO5YDO5ZDO5

3) 采用改进的最优权重表(LLR)法计算 R~(Δλ),即采用(16)~(31)式,得到三刺激值XLLRYLLRZLLR

4) 将XstdYstdZstd作为标准三刺激值,计算步骤2)、3)的XDO3YDO3ZDO3XDO5YDO5ZDO5XLLRYLLRZLLR与标准三刺激值之间的色差。

为了便于叙述分析,把六种标准照明体分为两类:一类是连续标准照明体,包括D65、A、D50;另一类是非连续标准照明体,包括F02、F07、F11,同时也将对应的两种配色函数的色差组合在一起。当固定了测量间隔和仪器函数形状(即选定sskew)时,在每种类型的标准照明体下,有1096×3×2个色差,分别取平均、最大和中间色差,以此衡量任何仪器函数形状下不同方法的效果优劣。

在测量数据间隔较大时,如10 nm和20 nm间隔,文献[ 2-4,8,15,17-18]中常考虑加权表方法。 本文首先研究10 nm和20 nm间隔的仿真结果。结果发现,无论用哪个指标进行评价,结论都是一样的。图5为10 nm间隔时三种方法的仿真结果。可以看到,对于连续标准照明体,选取21种不同的仪器函数,最优加权表法的最大色差均小于三点校正法和五点校正法;对于非连续标准照明体,仅有一种情况下最优加权表法的效果稍次于五点校正法,其余情况下最优加权表法均优于三点校正法和五点校正法。同时也可以看出,五点校正法的表现好于三点校正法。总之,在众多样本的色差中,即使是最大色差,改进的最优权重法的校正效果比其他两种方法更接近真实值,平均色差与中值色差更是如此。

图 5. 间隔为10 nm时三种校正法的最大色差比较。(a)连续标准照明体条件下;(b)非连续标准照明体条件下

Fig. 5. Maximum color differences of three correction methods with interval of 10 nm. (a) Under continuous standard light source condition; (b) under non-continuous standard light source condition

下载图片 查看所有图片

图6给出了测量间隔是20 nm时的仿真分析结果,可以看出,不管是连续标准照明体还是非连续标准照明体,最优加权表的表现是最好的,其次是五点校正,最差的是三点校正。另外,测量间隔是20 nm时的最大色差整体大于10 nm间隔下的。

图 6. 间隔为20 nm时三种校正法的最大色差比较。(a)连续标准照明体条件下;(b)非连续标准照明体条件下

Fig. 6. Maximum color differences of three correction methods with interval of 20 nm. (a) Under continuous standard light source condition; (b) under non-continuous standard light source condition

下载图片 查看所有图片

本文10 nm和20 nm间隔下的仿真结果与Li等[8]在仪器函数是对称三角形下的仿真结果一致。同时还发现,不管测量间隔是10 nm还是20 nm,当sskew=1时最优加权表方法表现最好。

最后考虑5 nm间隔时的仿真结果。在此情形下,CIE[16]推荐的是只要不涉及荧光测量,可以采用直接选择法计算三刺激值,即在5 nm间隔下利用所涉及的所有被积函数值的乘积求和形式计算三刺激值,该方法简称为DS法。在这里也把DS法纳入对比中。根据Li等[8]的结果,在应用DS方法之前,对5 nm间隔的测量反射比数据分别应用三点校正法(DO3+DS)和五点校正法 (DO5+DS)进行校正。图7给出了最大色差的比较结果。可以看出,误差均很小,都在10-3色差量级。另外,非对称的LLR方法表现最差(连续标准照明体)和次最差(非连续标准照明体)。表现最好的是CIE推荐的DS方法和三点或五点校正一起使用。这个结果也和Li等[8]在仪器函数是对称三角形下的仿真结果一致。根据同等表现,越简单越好的原则,对于5 nm间隔,建议先对测量数据进行三点校正,然后使用CIE推荐的DS方法计算三刺激值。

图 7. 间隔为5 nm时五种校正法的最大色差比较。(a)连续标准光源条件下;(b)非连续标准光源条件下

Fig. 7. Maximum color differences of five correction methods with interval of 5 nm. (a) Under continuous standard light source condition; (b) under non-continuous standard light source condition

下载图片 查看所有图片

需要说明的是,与10 nm和20 nm测量间隔下的结论不同,若采用平均或中值作为衡量指标,最好的方法是先采用五点校正的办法对测量数据进行校正,然后采用3阶多项式进行插值获得1 nm间隔的数据,最后采用1 nm乘积求和的方式计算三刺激值。

4 结论

针对仪器函数形状稍微偏离对称三角形的情形,研究了与指定仪器函数适配的最优权重表。通过引入非对称参数sskew控制仪器函数的倾斜程度,sskew的范围为[0.9,1.1]。该取值范围保证了仪器函数仅是稍微偏移,而当sskew=1时,仪器函数则是理想的对称三角形。

关于校正非对称仪器函数的测量过程,在传统的最优权重法的理论基础上,将非对称仪器函数的校正也放到了相应权重表的构建当中。研究结果表明,当仪器函数形状稍微偏离对称三角形时,对于10 nm和20 nm测量间隔的情形,推荐使用改进后的最优权重法计算三刺激值;对于5 nm测量间隔,建议先对测量数据进行三点校正,然后使用CIE推荐的DS方法计算三刺激值。

参考文献

[1] 徐海松. 颜色信息工程[M]. 2版. 杭州: 浙江大学出版社, 2015.

    Xu HS. Color information engineering[M]. 2nd ed. Hangzhou: Zhejiang University Press, 2015.

[2] American Society forTesting and Materials. Standard method for computing the colors of objects by using the CIE system: ASTM E308-13[S]. USA: ASTM, 2013.

[3] American Society forTesting and Materials. Standard practice for computing the colors of objects by using the CIE system: ASTM E308-17[S]. USA: ASTM, 2017.

[4] American Society forTesting and Materials. Standard practice for computing the colors of objects by using the CIE system: ASTM E308-18[S]. USA: ASTM, 2018.

[5] Oleari C. Spectral-reflectance-factor deconvolution and colorimetric calculations by local-power expansion[J]. Color Research & Application, 2000, 25(3): 176-185.

[6] Li C J, Oleari C, Melgosa M, et al. Methods for computing weighting tables based on local power expansion for tristimulus values computations[J]. Journal of the Optical Society of America A, 2011, 28(11): 2243-2252.

[7] Li C J, Luo M R, Wang G. Recent progress in computing weighting tables for calculating CIE tristimulus values[J]. Proceedings of SPIE, 2006, 6033: 603302.

[8] Li C J, Ronnier L M, Melgosa M, et al. Testing the accuracy of methods for the computation of CIE tristimulus values using weighting tables[J]. Color Research & Application, 2016, 41(2): 125-142.

[9] OhnoY. A flexible bandpass correction method for spectrometers[C]//AIC Color Conference. May 9-13, 2005, Granada, Spain. [S.l.: s.n.], 2005: 697- 700.

[10] Gardner J L. Bandwidth correction for LED chromaticity[J]. Color Research & Application, 2006, 31(5): 374-380.

[11] Woolliams E R, Baribeau R, Bialek A, et al. Spectrometer bandwidth correction for generalized bandpass functions[J]. Metrologia, 2011, 48(3): 164-172.

[12] Eichstädt S, Schmähling F, Wübbeler G, et al. Comparison of the Richardson-Lucy method and a classical approach for spectrometer bandpass correction[J]. Metrologia, 2013, 50(2): 107-118.

[13] Commission Internationale del'Eclairage. Effect of instrumental bandpass function and measurement interval on spectral quantities: CIE 214: 2014[S]. Vienna: CIE, 2014.

[14] Commission Internationale del'Eclairage. Colorimetry: CIE 15.2:1986[S]. 2nd ed. Vienna: CIE, 1986.

[15] Billmeyer F W, Fairman H S. CIE method for calculating tristimulus values[J]. Color Research & Application, 1987, 12(1): 27-36.

[16] Commission Internationale del'Eclairage. Colorimetry: CIE 015: 2018[S]. 4th ed. Vienna: CIE, 2018.

[17] 李长军, 崔桂华, 赵达尊, 等. 一种计算色度学三刺激值加权表的新方法[J]. 光学学报, 2003, 23(11): 1346-1353.

    Li C J, Cui G H, Zhao D Z, et al. A simple method for computing optimum weighting table for colorimetric tristimulus integration[J]. Acta Optica Sinica, 2003, 23(11): 1346-1353.

[18] Li C J, Luo M R, Rigg B. A new method for computing optimum weights for calculating CIE tristimulus values[J]. Color Research & Application, 2004, 29(2): 91-103.

李林璐, 王智峰, 李长军. 基于非对称三角形仪器函数的最优加权表方法[J]. 光学学报, 2020, 40(5): 0520002. Linlu Li, Zhifeng Wang, Changjun Li. Optimum Weighting Table Method Based on Asymmetric Triangular Instrumental Function[J]. Acta Optica Sinica, 2020, 40(5): 0520002.

本文已被 1 篇论文引用
被引统计数据来源于中国光学期刊网
引用该论文: TXT   |   EndNote

相关论文

加载中...

关于本站 Cookie 的使用提示

中国光学期刊网使用基于 cookie 的技术来更好地为您提供各项服务,点击此处了解我们的隐私策略。 如您需继续使用本网站,请您授权我们使用本地 cookie 来保存部分信息。
全站搜索
您最值得信赖的光电行业旗舰网络服务平台!