光谱色差与主分量特征融合的石窟表面风化智能评估

光谱色差与主分量特征融合的石窟表面风化智能评估 Intelligent evaluation of grotto surface weathering based on spectral chromatic aberration and principal component feature fusion 1 first-author 曹赤鹏（1995-），男，河北定州人，硕士研究生，主要从事多光谱图像处理、机器学习等方面的研究。E-mail：1309209924@qq.com http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870175&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870179&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870177&type= 22.01332855 32.00402832 曹赤鹏 （1995-），男，河北定州人，硕士研究生，主要从事多光谱图像处理、机器学习等方面的研究。E-mail： 1309209924@qq.com 1309209924@qq.com 1 corresp E-mail： hqwang@xauat.edu.cn E-mail： hqwang@xauat.edu.cn 王慧琴（1970-），女，山西长治人，博士，教授，博士生导师，1992年于西安冶金建筑学院获得学士学位，1997年于西安建筑科技大学获得硕士学位，2002年于西安交通大学获得博士学位，现为西安建筑科技大学信息与控制工程学院副院长，主要从事智能信息处理、信息理论与应用、信息技术与管理、数字建筑等研究。Email：hqwang@xauat.edu.cn http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870182&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870186&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870184&type= 22.01333618 32.00395966 王慧琴 （1970-），女，山西长治人，博士，教授，博士生导师，1992年于西安冶金建筑学院获得学士学位，1997年于西安建筑科技大学获得硕士学位，2002年于西安交通大学获得博士学位，现为西安建筑科技大学信息与控制工程学院副院长，主要从事智能信息处理、信息理论与应用、信息技术与管理、数字建筑等研究。Email： hqwang@xauat.edu.cn hqwang@xauat.edu.cn 1 2 2 2 西安建筑科技大学 信息与控制工程学院，陕西 西安 710055 西安建筑科技大学 信息与控制工程学院，陕西 西安 710055 School of Information and Control Engineering， Xi'an University of Architecture and Technology， Xi'an 710055， China School of Information and Control Engineering， Xi'an University of Architecture and Technology， Xi'an 710055， China 陕西省文物保护研究院，陕西 西安 710075 陕西省文物保护研究院，陕西 西安 710075 Shanxi Provincial Institute of Cultural Relics Protection， Xi’an 710075， China Shanxi Provincial Institute of Cultural Relics Protection， Xi’an 710075， China 针对单一光谱特征无法精细化表征石窟表面复杂风化区域风化类型及程度的问题，基于光谱分析与色度学理论，提出了一种光谱色差与主分量特征融合的石窟表面风化智能量化评估方法。重建石窟表面多光谱图像的反射光谱，计算每个像素点与基准点的色差，并使用主成分分析方法提取多光谱图像数据的主分量特征；然后融合石窟表面风化的光谱色差与主分量特征，对不同风化类型及程度进行表征；最后使用随机森林分类器，对石窟表面多光谱图像中每个像素点的风化程度进行智能评估。实验表明，融合光谱色差与主分量特征的方法对复杂风化区域不同风化类型及程度的表征能力优于单一光谱特征，石窟表面整体风化程度评估准确率为99.86%，Kappa系数为0.99，可实现复杂风化区域的精细化表征。 To overcome the problem of a single spectral feature not finely characterizing the type and degree of weathering in the complex weathering area on a grotto surface, this paper proposes an intelligent quantitative evaluation method for grotto surface weathering based on spectral analysis and colorimetric theory. First, we reconstruct the reflection spectrum of the multispectral image of the grotto surface, calculate the color difference between each pixel and reference point, and use the principal component analysis method to extract the principal component features from the multispectral image data. Then, we fuse the spectral color difference and principal component features of the weathering on the grotto surface to characterize different types and degrees of weathering. Finally, a random forest classifier is used to intelligently evaluate the weathering degree of each pixel in the grotto surface multispectral image. Experiments show that the method of fusing spectral chromatic aberration and principal component features performs better compared with a single spectral feature in characterizing different types and degrees of weathering in complex weathering areas. The evaluation accuracy of the overall grotto surface weathering degree is 99.86%, with a Kappa coefficient of 0.99. The proposed method can effectively realize a refined characterization of complex weathered areas. 石窟风化 光谱重建 主成分分析 特征融合 智能评估 cave weathering spectral reconstruction principal component analysis feature fusion intelligent evaluation 1 引　言 1 在大型石质文物表面风化等级评估方法研究中，不同风化类型及风化程度之间的表征问题一直是研究的热点。传统的风化程度检测主要依靠人工采集风化表面岩体的物理化学参数，通过计算各参数与风化类型及程度之间的相关性对不同风化类型及风化程度进行表征 ［ 1 1 ］ 。采集过程工作量巨大，并且采样点的规划极易受人为主观因素影响，难以满足大型石质文物表面复杂风化区域的精细化表征。 光谱成像技术 ［ 2 2 - 4 4 2-4 ］ 迅速发展，广泛应用于壁画颜料识别 ［ 5 5 ］ 、油画颜料信息提取 ［ 6 6 ］ 、石质文物表面风化程度检测 ［ 7 7 ］ 等工作。在石质文物表面风化程度检测方法研究中，国内外众多学者将多光谱成像技术和光谱分析技术相结合，通过采集石质文物表面风化离散光谱通道的光谱信息和空间信息对不同风化类型及风化程度区域进行稀疏表征。反射光谱不仅可以反映物质的状态特征，还反映了表面颜色、表面粒度等物理特征 ［ 8 8 ］ 。柴勃隆 ［ 9 9 ］ 等将多光谱成像技术与石窟考古研究相结合，利用壁画表面不同物质反射光谱特征的不同，成功对莫高窟壁画在可见光范围内肉眼无法辨识的彩色区域进行了修复。周霄 ［ 10 10 ］ 等利用光谱成像技术采集了云冈石窟表面风化区域的光谱反射率数据，利用光谱反射率幅值的差异对风化程度进行表征，有效区分出强风化和弱风化区域。但反射光谱曲线特征变化较为平缓，光谱反射率数据存在较大冗余，不能实现数据的快速处理，并且对于石窟表面复杂风化区域，光谱反射率幅值特征较为单一，难以达到精细化表征，需要探索一种有效的表征方法。 颜色是对可见光特征的表征，可见光的波长决定了人类眼睛所感知的颜色。对于石窟而言，不同风化类型及风化程度区域的颜色由其表面的反射光谱决定。不同风化区域的颜色差异真实地反映了石窟表面不同风化类型及风化程度之间的差别。主成分分析方法 ［ 11 11 ］ 是对多波段数据的线性变换，将原始光谱数据转换到一个数据维度较少的新空间，通过计算使光谱数据的差异达到最大，对于增强石窟表面风化光谱信息含量、隔离噪声及减少数据维数非常有效 ［ 12 12 ］ 。基于此，本文提出了一种基于光谱色差与主分量特征融合的石窟表面风化智能评估方法。通过重建石窟表面风化光谱反射率，计算多光谱成像数据中每个像素点与基准点的色差；并对石窟表面原始多光谱成像数据进行主成分分析，提取主分量特征数据，将光谱色差数据与主分量特征数据进行融合，建立不同风化类型及风化程度表征数据库；利用随机森林算法（Random Forest，RF）对石窟表面复杂风化区域每个像素点进行智能评估，为大型石质文物表面复杂风化区域精细化评估提供可行性方案。 2 石窟表面不同风化区域光谱色差特性分析 1 石窟表面在长期风化作用下其表面物质、表面颜色发生改变。表面岩体中的钙、镁等氢氧化物在空气中二氧化碳的长期作用下，使得岩体表面颜色变为白色，风化程度越深颜色表征越为明显 ［ 13 13 ］ 。 以黄色砂砾岩岩体石窟为例，利用多光谱成像系统对石窟表面风化区域进行多光谱成像数据采集，如 图1 图1 所示为石窟表面不同风化类型及风化程度纯净风化区域的多光谱成像平均反射光谱数据。随盐析风化程度的加深，其所对应的反射光谱的幅值越高；随积尘风化的加深，其所对应的反射光谱幅值越低。反射光谱不仅反映石窟表面岩体物质成分的改变，同时还反映不同风化类型及风化程度区域间颜色的差异。为了准确分析石窟表面风化类型及风化程度的改变对颜色变化的影响，分别对不同的风化类型及风化程度区域的反射光谱进行CIE1976L*a*b*、CIE1931XYZ以及RGB颜色空间的线性变换，通过各参数值的变化来反映石窟表面不同风化类型及程度区域间表面颜色的变化情况。 表1 表1 为石窟表面风化区域的CIE1976L*a*b*颜色空间、CIE1931XYZ颜色空间以及R，G，B各参数随不同风化类型及风化程度的变化情况，其中变量 <math id="M1"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">s</mi><mi mathvariant="normal">t</mi><mi mathvariant="normal">r</mi><mi mathvariant="normal">o</mi><mi mathvariant="normal">n</mi><mi mathvariant="normal">g</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870199&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870196&type= 8.97466660 3.72533321 ， <math id="M2"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">w</mi><mi mathvariant="normal">e</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">k</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870202&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870201&type= 8.21266651 3.72533321 ， <math id="M3"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">s</mi><mi mathvariant="normal">l</mi><mi mathvariant="normal">i</mi><mi mathvariant="normal">g</mi><mi mathvariant="normal">h</mi><mi mathvariant="normal">t</mi><mi mathvariant="normal">l</mi><mi mathvariant="normal">y</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870209&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870206&type= 9.82133293 3.72533321 ， <math id="M4"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">d</mi><mi mathvariant="normal">u</mi><mi mathvariant="normal">s</mi><mi mathvariant="normal">t</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870214&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870211&type= 7.36600018 3.72533321 ， <math id="M5"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">b</mi><mi mathvariant="normal">e</mi><mi mathvariant="normal">n</mi><mi mathvariant="normal">c</mi><mi mathvariant="normal">h</mi><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">r</mi><mi mathvariant="normal">k</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870220&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870218&type= 12.69999981 3.72533321 分别表示强盐析风化，弱盐析风化，微盐析风化，积尘风化，基准点。盐析风化类型中风化程度1强于风化程度2，积尘风化类型中风化程度2强于风化程度1。通过对比 表1 表1 中各参数值可以发现，随盐析风化程度的加深，石窟表面岩体中钙、镁等的氢氧化物溶液与空气中的二氧化碳反应生成白色的碳酸钙或碳酸镁等的结晶，导致表面岩体颜色变白，使得岩体表面亮度逐渐增加，三种颜色空间的各参数值呈现逐步上升趋势。对于积尘风化区域，由于石窟表面酥化，产生的碎屑掺杂着尘土覆盖于石窟表面，使得积尘风化区域的明亮度变低，相对于盐析风化区域其RGB各通道值最低。不同的风化类型及风化程度区域间其表面颜色的差异可以用光谱色差来表征。以基准点L*a*b*值为标准，计算不同风化类型及风化程度间的光谱色差值，用 <math id="M6"><mi mathvariant="normal">Δ</mi><mi mathvariant="normal">E</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870226&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870223&type= 4.99533367 2.79399991 表示。在国际颜色标准CIE1976L*a*b*的评价体系中，色差的判断标准如 表2 表2 所示。不同风化类型及风化程度区域与基准间的光谱色差值如 表3 表3 所示。 Weathering type <math id="M7"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">s</mi><mi mathvariant="normal">t</mi><mi mathvariant="normal">r</mi><mi mathvariant="normal">o</mi><mi mathvariant="normal">n</mi><mi mathvariant="normal">g</mi><mn>1</mn></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870233&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870231&type= 8.72066593 3.21733332 <math id="M8"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">s</mi><mi mathvariant="normal">t</mi><mi mathvariant="normal">r</mi><mi mathvariant="normal">o</mi><mi mathvariant="normal">n</mi><mi mathvariant="normal">g</mi><mn>2</mn></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870241&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870238&type= 8.80533314 3.21733332 <math id="M9"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">w</mi><mi mathvariant="normal">e</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">k</mi><mn>1</mn></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870247&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870244&type= 8.04333401 3.21733332 <math id="M10"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">w</mi><mi mathvariant="normal">e</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">k</mi><mn>2</mn></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870254&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870251&type= 8.04333401 3.21733332 <math id="M11"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">s</mi><mi mathvariant="normal">l</mi><mi mathvariant="normal">i</mi><mi mathvariant="normal">g</mi><mi mathvariant="normal">h</mi><mi mathvariant="normal">t</mi><mi mathvariant="normal">l</mi><mi mathvariant="normal">y</mi><mn>1</mn></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870259&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870257&type= 9.48266602 3.21733332 <math id="M12"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">s</mi><mi mathvariant="normal">l</mi><mi mathvariant="normal">i</mi><mi mathvariant="normal">g</mi><mi mathvariant="normal">h</mi><mi mathvariant="normal">t</mi><mi mathvariant="normal">l</mi><mi mathvariant="normal">y</mi><mn>2</mn></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870265&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870262&type= 9.48266602 3.21733332 <math id="M13"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">d</mi><mi mathvariant="normal">u</mi><mi mathvariant="normal">s</mi><mi mathvariant="normal">t</mi><mn>1</mn></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870271&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870268&type= 7.28133297 3.21733332 <math id="M14"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">d</mi><mi mathvariant="normal">u</mi><mi mathvariant="normal">s</mi><mi mathvariant="normal">t</mi><mn>2</mn></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870276&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870275&type= 7.28133297 3.21733332 <math id="M15"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">b</mi><mi mathvariant="normal">e</mi><mi mathvariant="normal">n</mi><mi mathvariant="normal">c</mi><mi mathvariant="normal">h</mi><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">r</mi><mi mathvariant="normal">k</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870706&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870704&type= 11.17599964 3.21733332 L* 46.4 43.7 42.6 37.7 33.0 28.0 8.0 7.9 0 a* 3.5 4.7 3.5 3.1 3.6 2.5 0.2 0.2 0 b* 8.9 10.5 8.3 8.7 9.0 7.6 7.9 7.8 0 X 15.7 13.9 13.0 10.0 7.7 5.5 0.8 0.8 0 Y 15.6 13.6 12.9 9.9 7.5 5.5 0.9 0.9 0 Z 9.8 8.1 8.1 6.0 4.4 3.2 0.3 0.3 0 R 121 117 111 98 88 74 27 27 0 G 108 100 98 87 75 64 23 23 0 B 95 85 86 74 63 55 10 10 0 <math id="M16"><mi mathvariant="normal">Δ</mi><mi>E</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870291&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870288&type= 4.31799984 2.37066650 Chromatic aberration degree identification 0~0.5 Trace 0.5~1.5 Lightweight 1.5~3.0 Can feel 3.0~6.0 obvious 6.0~12.0 tremendous <math id="M17"><mo>≥</mo></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870297&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870294&type= 3.21733332 2.79399991 12.0 huge Color difference type <math id="M18"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">s</mi><mi mathvariant="normal">t</mi><mi mathvariant="normal">r</mi><mi mathvariant="normal">o</mi><mi mathvariant="normal">n</mi><mi mathvariant="normal">g</mi></mrow></msub><mo>-</mo><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">b</mi><mi mathvariant="normal">e</mi><mi mathvariant="normal">n</mi><mi mathvariant="normal">c</mi><mi mathvariant="normal">h</mi><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">r</mi><mi mathvariant="normal">k</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870307&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870304&type= 23.19866753 3.21733332 <math id="M19"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">w</mi><mi mathvariant="normal">e</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">k</mi></mrow></msub><mo>-</mo><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">b</mi><mi mathvariant="normal">e</mi><mi mathvariant="normal">n</mi><mi mathvariant="normal">c</mi><mi mathvariant="normal">h</mi><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">r</mi><mi mathvariant="normal">k</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870315&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870311&type= 22.52133369 3.21733332 <math id="M20"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">s</mi><mi mathvariant="normal">l</mi><mi mathvariant="normal">i</mi><mi mathvariant="normal">g</mi><mi mathvariant="normal">h</mi><mi mathvariant="normal">t</mi><mi mathvariant="normal">l</mi><mi mathvariant="normal">y</mi></mrow></msub><mo>-</mo><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">b</mi><mi mathvariant="normal">e</mi><mi mathvariant="normal">n</mi><mi mathvariant="normal">c</mi><mi mathvariant="normal">h</mi><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">r</mi><mi mathvariant="normal">k</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870321&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870318&type= 23.96066666 3.21733332 <math id="M21"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">d</mi><mi mathvariant="normal">u</mi><mi mathvariant="normal">s</mi><mi mathvariant="normal">t</mi></mrow></msub><mo>-</mo><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">b</mi><mi mathvariant="normal">e</mi><mi mathvariant="normal">n</mi><mi mathvariant="normal">c</mi><mi mathvariant="normal">h</mi><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">r</mi><mi mathvariant="normal">k</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870328&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870324&type= 21.75933456 3.21733332 <math id="M22"><mi mathvariant="normal">Δ</mi><mi mathvariant="normal">E</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870344&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870342&type= 4.23333359 2.37066650 76 45.19 38.81 29.12 11.24 <math id="M23"><mi mathvariant="normal">Δ</mi><mi mathvariant="normal">E</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870344&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870342&type= 4.23333359 2.37066650 94 45.19 38.81 29.12 11.24 <math id="M24"><mi mathvariant="normal">Δ</mi><mi mathvariant="normal">E</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870344&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870342&type= 4.23333359 2.37066650 2 000 32.37 27.01 19.57 8.2 通过 表3 表3 可以看出，由于石窟表面不同风化程度区域表面颜色的不同，导致不同风化程度区域与基准间的色差存在较大差异。盐析风化作用越强与基准间的色差值越高，强盐析风化与弱盐析风化之间色差为6.38、弱盐析风化与微盐析风化之间色差为9.69、微盐析风化与积尘风化之间色差为17.88，在色差判断标准中其差异均处于很大或者非常大水平。根据以上分析，光谱色差可以对石窟表面不同风化类型及风化程度进行有效表征，为光谱色差在石窟表面风化智能评估提供理论支撑与数据基础。 3 相关原理 1 3.1　基于多光谱成像数据的石窟表面光谱反射率重建 2 石窟表面风化多光谱成像数据由多波段组成，具有较强的离散特性，使得光谱色差计算时存在较大误差。光谱反射率重建可以消除多光谱成像数据的离散特性，重建后的光谱信息可以充分表征石窟不同风化区域表面物质及颜色差异，精确的重现不同风化类型及风化程度区域的颜色信息。 在前期光谱重建理论研究与实验基础之上，本文采用多核支持向量回归方法，重建石窟表面多光谱成像的光谱反射率数据。与伪逆法、多项式法、单核支持向量回归等传统重建方法相比，多核支持向量回归将局部核函数和全局核函数相结合，通过局部核函数对测试点较近函数值的强学习能力与全局核函数对测试点较远函数值的强泛化能力，能够有效解决传统重建模型泛化能力差、重建精度低的问题，并在光谱重建颜色复原实验中获得较高的精度 ［ 14 14 ］ 。本文采用局部核函数Poly核函数与全局核函数Cauchy核函数相结合的方法，引入到支持向量回归中，利用多核支持向量回归对石窟表面不同风化类型及风化程度的光谱反射率进行重建。回归函数为： <math display="block" id="M25"><mi>R</mi><mo stretchy="false">(</mo><msub><mrow><mi>x</mi></mrow><mrow><mi>i</mi></mrow></msub><mo stretchy="false">)</mo><mo>=</mo><munder><mo movablelimits="false">∑</mo><mrow><msub><mrow><mi>x</mi></mrow><mrow><mi>i</mi></mrow></msub></mrow></munder><mrow><mo stretchy="false">(</mo><msub><mrow><mi>a</mi></mrow><mrow><mi>i</mi></mrow></msub><mo>-</mo><msubsup><mrow><mi>a</mi></mrow><mrow><mi>i</mi></mrow><mrow><mi mathvariant="normal">*</mi></mrow></msubsup><mo stretchy="false">)</mo></mrow><mo>·</mo><msub><mrow><mi>K</mi></mrow><mrow><mi>c</mi><mi>p</mi></mrow></msub><mo stretchy="false">(</mo><msub><mrow><mi>x</mi></mrow><mrow><mi>i</mi></mrow></msub><mo>,</mo><mi>y</mi><mo stretchy="false">)</mo><mo>+</mo><mi>b</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870358&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870354&type= 55.37200165 7.45066643 ， 其中： <math id="M26"><mi>R</mi><mo stretchy="false">(</mo><msub><mrow><mi>x</mi></mrow><mrow><mi>i</mi></mrow></msub><mo stretchy="false">)</mo></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870364&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870361&type= 8.89000034 4.48733330 为波长为 <math id="M27"><msub><mrow><mi>λ</mi></mrow><mrow><mi>i</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870370&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870367&type= 2.53999996 3.72533321 处的不同风化类型及风化程度的光谱反射率重建结果； <math id="M28"><msub><mrow><mi>K</mi></mrow><mrow><mi>c</mi><mi>p</mi></mrow></msub><mo stretchy="false">(</mo><msub><mrow><mi>x</mi></mrow><mrow><mi>i</mi></mrow></msub><mo>,</mo><mi>y</mi><mo stretchy="false">)</mo></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870375&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870373&type= 15.40933418 4.74133301 表示重建光谱反射率的核函数； <math id="M29"><mi>i</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870400&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870397&type= 1.10066664 2.79399991 表示需要重建的光谱响应值个数，取值范围为 <math id="M30"><mfenced separators="|"><mrow><mi>i</mi><mo>=</mo><mn mathvariant="normal">1,2</mn><mo>,</mo><mo>⋅</mo><mo>⋅</mo><mo>⋅</mo><mo>,</mo><mi>n</mi></mrow></mfenced></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870387&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870383&type= 23.11400032 3.97933316 ； <math id="M31"><msub><mrow><mi>x</mi></mrow><mrow><mi>i</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870571&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870568&type= 2.79399991 3.72533321 表示第 <math id="M32"><mi>i</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870400&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870397&type= 1.10066664 2.79399991 个测试样本光谱响应数据； <math id="M33"><mi>y</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870406&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870404&type= 1.86266661 3.72533321 表示训练样本光谱反射率数据； <math id="M34"><mi>b</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870412&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870409&type= 1.86266661 2.79399991 表示测量误差。通过多核核函数对重建算法进行优化，多核核函数构造如 式（2） 式（2） 所示： <math display="block" id="M35"><msub><mrow><mi>K</mi></mrow><mrow><mi>c</mi><mi>p</mi></mrow></msub><mo stretchy="false">(</mo><msub><mrow><mi>x</mi></mrow><mrow><mi>i</mi></mrow></msub><mo>,</mo><mi>y</mi><mo stretchy="false">)</mo><mo>=</mo><mfrac><mrow><mn mathvariant="normal">1</mn></mrow><mrow><mrow><mrow><msup><mrow><mfenced open="‖" close="‖" separators="|"><mrow><mi>x</mi><mo>-</mo><msub><mrow><mi>x</mi></mrow><mrow><mi>i</mi></mrow></msub></mrow></mfenced></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msup></mrow><mo>/</mo><mrow><mi>σ</mi><mo>+</mo><mn mathvariant="normal">1</mn></mrow></mrow></mrow></mfrac><mo>·</mo><mfenced open="[" close="]" separators="|"><mrow><msup><mrow><mfenced separators="|"><mrow><msup><mrow><mi>x</mi></mrow><mrow><mi>T</mi></mrow></msup><msub><mrow><mi>x</mi></mrow><mrow><mi>i</mi></mrow></msub><mo>+</mo><mn mathvariant="normal">1</mn></mrow></mfenced></mrow><mrow><mi>d</mi></mrow></msup></mrow></mfenced></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870419&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870415&type= 69.93466949 10.75266743 ， 其中： <math id="M36"><mi>x</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870425&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870423&type= 2.03200006 2.79399991 表示石窟表面待重建区域不同风化类型及风化程度的多光谱相机响应输入值； <math id="M37"><mi>σ</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870435&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870431&type= 1.69333339 2.79399991 表示柯西核函数的宽度，作为各支持向量之间的关联系数； <math id="M38"><mi>d</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870442&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870439&type= 2.20133328 2.79399991 表示多项式核函数阶数。则石窟表面风化多光谱成像数据中每个像素点的光谱反射率重建结果可以表示为 <math id="M39"><msup><mrow><mfenced separators="|"><mrow><msub><mrow><mi>R</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub><mo stretchy="false">(</mo><mi>x</mi><mo stretchy="false">)</mo><mo>,</mo><msub><mrow><mi>R</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub><mo stretchy="false">(</mo><mi>x</mi><mo stretchy="false">)</mo><mo>,</mo><mo>⋅</mo><mo>⋅</mo><mo>⋅</mo><mo>,</mo><msub><mrow><mi>R</mi></mrow><mrow><mi>i</mi></mrow></msub><mo stretchy="false">(</mo><mi>x</mi><mo stretchy="false">)</mo></mrow></mfenced></mrow><mrow><mi mathvariant="normal">T</mi></mrow></msup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870450&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870446&type= 43.85733414 5.58799982 。 3.2　石窟表面风化像素级光谱色差计算 2 利用支持向量回归对石窟表面不同风化类型及风化程度区域每个像素点的光谱反射率进行重建；将重建后的光谱反射率，从光谱空间转换到CIE1976L*a*b*颜色空间，计算每个像素点的L*，a*，b*参数值。在CIE1976L*a*b*颜色空间中，颜色通过三种对比的刺激值进行表示，分别为黑-白、绿-红、蓝-黄。L*表示黑/白刺激值，当L*为0时，表示黑色，当L*为100时表示白色；a*表示绿/红刺激值，当a*为正值时表示红的程度，当a*为负值时表示绿的程度；b*表示蓝/黄刺激值，当b*为正值时表示黄的程度，当b*为负值时表示蓝的程度。CIE1976L*a*b*颜色空间的各参数值是由CIE1931XYZ颜色空间通过线性变换而来。石窟表面风化多光谱成像数据计算 X、Y、Z 三刺激值的计算过程如式（3）~ 式（7） 式（7） 所示： <math display="block" id="M40"><mi>C</mi><mi>I</mi><mi>E</mi><mn mathvariant="normal">1931</mn><mi>X</mi><mi>Y</mi><mi>Z</mi><mo>=</mo><mfenced open="{" close="" separators="|"><mrow><mtable columnalign="left"><mtr><mtd><mi>X</mi><mo>=</mo><mi>k</mi><mstyle displaystyle="true"><munderover><mo largeop="true">∑</mo><mrow><mn mathvariant="normal">380</mn></mrow><mrow><mn mathvariant="normal">780</mn></mrow></munderover></mstyle><mrow><mi>φ</mi><mo stretchy="false">(</mo><mi>λ</mi><mo stretchy="false">)</mo></mrow><mo>·</mo><mover><mrow><mi>x</mi></mrow><mo>¯</mo></mover><mo stretchy="false">(</mo><mi>λ</mi><mo stretchy="false">)</mo><mi mathvariant="normal">Δ</mi><mi>λ</mi></mtd></mtr><mtr><mtd><mi>Y</mi><mo>=</mo><mi>k</mi><mstyle displaystyle="true"><munderover><mo largeop="true">∑</mo><mrow><mn mathvariant="normal">380</mn></mrow><mrow><mn mathvariant="normal">780</mn></mrow></munderover></mstyle><mrow><mi>φ</mi><mo stretchy="false">(</mo><mi>λ</mi><mo stretchy="false">)</mo></mrow><mo>·</mo><mover><mrow><mi>y</mi></mrow><mo>¯</mo></mover><mo stretchy="false">(</mo><mi>λ</mi><mo stretchy="false">)</mo><mi mathvariant="normal">Δ</mi><mi>λ</mi></mtd></mtr><mtr><mtd><mi>Z</mi><mo>=</mo><mi>k</mi><mstyle displaystyle="true"><munderover><mo largeop="true">∑</mo><mrow><mn mathvariant="normal">380</mn></mrow><mrow><mn mathvariant="normal">780</mn></mrow></munderover></mstyle><mrow><mi>φ</mi><mo stretchy="false">(</mo><mi>λ</mi><mo stretchy="false">)</mo></mrow><mo>·</mo><mover><mrow><mi>z</mi></mrow><mo>¯</mo></mover><mo stretchy="false">(</mo><mi>λ</mi><mo stretchy="false">)</mo><mi mathvariant="normal">Δ</mi><mi>λ</mi></mtd></mtr></mtable></mrow></mfenced></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870457&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870454&type= 66.20932770 29.63333130 ， <math display="block" id="M41"><mtable><mtr><mtd><mi>C</mi><mi>I</mi><mi>E</mi><mn mathvariant="normal">1976</mn><mi>L</mi><mi mathvariant="normal">*</mi><mi>a</mi><mi mathvariant="normal">*</mi><mi>b</mi><mi mathvariant="normal">*</mi><mo>=</mo></mtd></mtr><mtr><mtd><mfenced open="{" close="" separators="|"><mrow><mtable columnalign="left"><mtr><mtd><mi>L</mi><mi mathvariant="normal">*</mi><mo>=</mo><mn mathvariant="normal">116</mn><mi>f</mi><mfenced separators="|"><mrow><mfrac><mrow><mi>Y</mi></mrow><mrow><msub><mrow><mi>Y</mi></mrow><mrow><mi>N</mi></mrow></msub></mrow></mfrac></mrow></mfenced><mo>-</mo><mn mathvariant="normal">16</mn></mtd></mtr><mtr><mtd><mi>a</mi><mi mathvariant="normal">*</mi><mo>=</mo><mn mathvariant="normal">500</mn><mfenced separators="|"><mrow><mi>f</mi><mfenced separators="|"><mrow><mfrac><mrow><mi>X</mi></mrow><mrow><msub><mrow><mi>X</mi></mrow><mrow><mi>N</mi></mrow></msub></mrow></mfrac></mrow></mfenced><mo>-</mo><mi>f</mi><mfenced separators="|"><mrow><mfrac><mrow><mi>Y</mi></mrow><mrow><msub><mrow><mi>Y</mi></mrow><mrow><mi>N</mi></mrow></msub></mrow></mfrac></mrow></mfenced></mrow></mfenced></mtd></mtr><mtr><mtd><mi>b</mi><mi mathvariant="normal">*</mi><mo>=</mo><mn mathvariant="normal">200</mn><mfenced separators="|"><mrow><mi>f</mi><mfenced separators="|"><mrow><mfrac><mrow><mi>Y</mi></mrow><mrow><msub><mrow><mi>Y</mi></mrow><mrow><mi>N</mi></mrow></msub></mrow></mfrac></mrow></mfenced><mo>-</mo><mi>f</mi><mfenced separators="|"><mrow><mfrac><mrow><mi>Z</mi></mrow><mrow><msub><mrow><mi>Z</mi></mrow><mrow><mi>N</mi></mrow></msub></mrow></mfrac></mrow></mfenced></mrow></mfenced></mtd></mtr></mtable></mrow></mfenced></mtd></mtr></mtable></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870464&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870460&type= 44.45000076 39.03133011 ， <math display="block" id="M42"><mtable><mtr><mtd><mi>f</mi><mfenced separators="|"><mrow><mfrac><mrow><mi>X</mi></mrow><mrow><msub><mrow><mi>X</mi></mrow><mrow><mi>N</mi></mrow></msub></mrow></mfrac></mrow></mfenced><mo>=</mo></mtd></mtr><mtr><mtd><mfenced open="{" close="" separators="|"><mrow><mtable columnalign="left"><mtr><mtd><msup><mrow><mfenced separators="|"><mrow><mfrac><mrow><mi>X</mi></mrow><mrow><msub><mrow><mi>X</mi></mrow><mrow><mi>N</mi></mrow></msub></mrow></mfrac></mrow></mfenced></mrow><mrow><mfrac><mrow><mn mathvariant="normal">1</mn></mrow><mrow><mn mathvariant="normal">3</mn></mrow></mfrac></mrow></msup><mtext>                      </mtext><mfrac><mrow><mi>X</mi></mrow><mrow><msub><mrow><mi>X</mi></mrow><mrow><mi>N</mi></mrow></msub></mrow></mfrac><mo>></mo><mn mathvariant="normal">0.008</mn><mtext> </mtext><mn mathvariant="normal">856</mn></mtd></mtr><mtr><mtd><mn mathvariant="normal">7.787</mn><mfenced separators="|"><mrow><mfrac><mrow><mi>X</mi></mrow><mrow><msub><mrow><mi>X</mi></mrow><mrow><mi>N</mi></mrow></msub></mrow></mfrac></mrow></mfenced><mo>+</mo><mfrac><mrow><mn mathvariant="normal">16</mn></mrow><mrow><mn mathvariant="normal">116</mn></mrow></mfrac><mtext>      </mtext><mfrac><mrow><mi>X</mi></mrow><mrow><msub><mrow><mi>X</mi></mrow><mrow><mi>N</mi></mrow></msub></mrow></mfrac><mo>≤</mo><mn mathvariant="normal">0.008</mn><mtext> </mtext><mn mathvariant="normal">856</mn></mtd></mtr></mtable></mrow></mfenced></mtd></mtr></mtable></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870471&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870467&type= 61.63733673 34.88266373 ， <math display="block" id="M43"><mi>f</mi><mfenced separators="|"><mrow><mfrac><mrow><mi>Y</mi></mrow><mrow><msub><mrow><mi>Y</mi></mrow><mrow><mi>N</mi></mrow></msub></mrow></mfrac></mrow></mfenced><mo>=</mo><mfenced open="{" close="" separators="|"><mrow><mtable columnalign="left"><mtr><mtd><msup><mrow><mfenced separators="|"><mrow><mfrac><mrow><mi>Y</mi></mrow><mrow><msub><mrow><mi>Y</mi></mrow><mrow><mi>N</mi></mrow></msub></mrow></mfrac></mrow></mfenced></mrow><mrow><mfrac><mrow><mn mathvariant="normal">1</mn></mrow><mrow><mn mathvariant="normal">3</mn></mrow></mfrac></mrow></msup><mtext>                      </mtext><mfrac><mrow><mi>Y</mi></mrow><mrow><msub><mrow><mi>Y</mi></mrow><mrow><mi>N</mi></mrow></msub></mrow></mfrac><mo>></mo><mn mathvariant="normal">0.008856</mn></mtd></mtr><mtr><mtd><mn mathvariant="normal">7.787</mn><mfenced separators="|"><mrow><mfrac><mrow><mi>Y</mi></mrow><mrow><msub><mrow><mi>Y</mi></mrow><mrow><mi>N</mi></mrow></msub></mrow></mfrac></mrow></mfenced><mo>+</mo><mfrac><mrow><mn mathvariant="normal">16</mn></mrow><mrow><mn mathvariant="normal">116</mn></mrow></mfrac><mtext>      </mtext><mfrac><mrow><mi>Y</mi></mrow><mrow><msub><mrow><mi>Y</mi></mrow><mrow><mi>N</mi></mrow></msub></mrow></mfrac><mo>≤</mo><mn mathvariant="normal">0.008856</mn></mtd></mtr></mtable></mrow></mfenced></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870476&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870475&type= 77.97799683 23.53733444 ， <math display="block" id="M44"><mtable><mtr><mtd><mi>f</mi><mo stretchy="false">(</mo><mfrac><mrow><mi>Z</mi></mrow><mrow><msub><mrow><mi>Z</mi></mrow><mrow><mi>N</mi></mrow></msub></mrow></mfrac><mo stretchy="false">)</mo><mo>=</mo></mtd></mtr><mtr><mtd><mfenced open="{" close="" separators="|"><mrow><mtable columnalign="left"><mtr><mtd><msup><mrow><mo stretchy="false">(</mo><mfrac><mrow><mi>Z</mi></mrow><mrow><msub><mrow><mi>Z</mi></mrow><mrow><mi>N</mi></mrow></msub></mrow></mfrac><mo stretchy="false">)</mo></mrow><mrow><mfrac><mrow><mn mathvariant="normal">1</mn></mrow><mrow><mn mathvariant="normal">3</mn></mrow></mfrac></mrow></msup><mtext>                      </mtext><mfrac><mrow><mi>Z</mi></mrow><mrow><msub><mrow><mi>Z</mi></mrow><mrow><mi>N</mi></mrow></msub></mrow></mfrac><mo>></mo><mn mathvariant="normal">0.008</mn><mtext> </mtext><mn mathvariant="normal">856</mn></mtd></mtr><mtr><mtd><mn mathvariant="normal">7.787</mn><mo stretchy="false">(</mo><mfrac><mrow><mi>Z</mi></mrow><mrow><msub><mrow><mi>Z</mi></mrow><mrow><mi>N</mi></mrow></msub></mrow></mfrac><mo stretchy="false">)</mo><mo>+</mo><mfrac><mrow><mn mathvariant="normal">16</mn></mrow><mrow><mn mathvariant="normal">116</mn></mrow></mfrac><mtext>      </mtext><mfrac><mrow><mi>Z</mi></mrow><mrow><msub><mrow><mi>Z</mi></mrow><mrow><mi>N</mi></mrow></msub></mrow></mfrac><mo>≤</mo><mn mathvariant="normal">0.008</mn><mtext> </mtext><mn mathvariant="normal">856</mn></mtd></mtr></mtable></mrow></mfenced></mtd></mtr></mtable></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870482&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870479&type= 61.55266953 29.80266571 ， 其中： X，Y，Z 分别是CIE1931XYZ颜色空间的三刺激值； <math id="M45"><msub><mrow><mi>X</mi></mrow><mrow><mi>n</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870488&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870485&type= 3.80999994 3.72533321 ， <math id="M46"><msub><mrow><mi>Y</mi></mrow><mrow><mi>n</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870493&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870490&type= 4.06400013 3.72533321 ， <math id="M47"><msub><mrow><mi>Z</mi></mrow><mrow><mi>n</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870500&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870497&type= 3.47133350 3.72533321 表示标准白板的三刺激值，并通过 <math id="M48"><msub><mrow><mi>Y</mi></mrow><mrow><mi>n</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870505&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870502&type= 4.06400013 3.72533321 值为100进行标准化处理。通过颜色空间的线性变换得到了石窟表面不同风化类型及风化程度区域内每个像素点的L*，a*，b*参数值。并以基准点为计算基准对石窟表面不同风化类型及风化程度区域的色差进行计算，用 <math id="M49"><mi mathvariant="normal">Δ</mi><msub><mrow><mi>E</mi></mrow><mrow><mi>i</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870510&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870507&type= 5.84200001 3.72533321 表示，则计算公式如下： <math display="block" id="M50"><mi mathvariant="normal">Δ</mi><msub><mrow><mi>E</mi></mrow><mrow><mi>i</mi></mrow></msub><mo>=</mo><mo stretchy="false">(</mo><msubsup><mrow><mi>L</mi></mrow><mrow><mn mathvariant="normal">0</mn></mrow><mrow><mi mathvariant="normal">*</mi></mrow></msubsup><mo>-</mo><msubsup><mrow><mi>L</mi></mrow><mrow><mi>i</mi></mrow><mrow><mi mathvariant="normal">*</mi></mrow></msubsup><msup><mrow><mo stretchy="false">)</mo></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msup><mo>+</mo><mo stretchy="false">(</mo><msubsup><mrow><mi>a</mi></mrow><mrow><mn mathvariant="normal">0</mn></mrow><mrow><mi mathvariant="normal">*</mi></mrow></msubsup><mo>-</mo><msubsup><mrow><mi>a</mi></mrow><mrow><mi>i</mi></mrow><mrow><mi mathvariant="normal">*</mi></mrow></msubsup><msup><mrow><mo stretchy="false">)</mo></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msup><mo>+</mo><mo stretchy="false">(</mo><msubsup><mrow><mi>b</mi></mrow><mrow><mn mathvariant="normal">0</mn></mrow><mrow><mi mathvariant="normal">*</mi></mrow></msubsup><mo>-</mo><msubsup><mrow><mi>b</mi></mrow><mrow><mi>i</mi></mrow><mrow><mi mathvariant="normal">*</mi></mrow></msubsup><msup><mrow><mo stretchy="false">)</mo></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msup><mtext>   </mtext><mo stretchy="false">(</mo><mi mathvariant="normal">i</mi><mo>=</mo><mn mathvariant="normal">1,2</mn><mo>,</mo><mo>⋅</mo><mo>⋅</mo><mo>⋅</mo><mo>,</mo><mi>n</mi><mo stretchy="false">)</mo></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870516&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870513&type= 90.93199921 4.65666676 ， 其中： <math id="M51"><msubsup><mrow><mi>L</mi></mrow><mrow><mn mathvariant="normal">0</mn></mrow><mrow><mi mathvariant="normal">*</mi></mrow></msubsup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870522&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870519&type= 3.55599999 3.97933316 ， <math id="M52"><msubsup><mrow><mi>a</mi></mrow><mrow><mn mathvariant="normal">0</mn></mrow><mrow><mi mathvariant="normal">*</mi></mrow></msubsup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870529&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870526&type= 3.04800010 3.97933316 ， <math id="M53"><msubsup><mrow><mi>b</mi></mrow><mrow><mn mathvariant="normal">0</mn></mrow><mrow><mi mathvariant="normal">*</mi></mrow></msubsup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870537&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870533&type= 3.04800010 3.97933316 表示基准点的CIE1976L*a*b*值； <math id="M54"><msubsup><mrow><mi>L</mi></mrow><mrow><mi>i</mi></mrow><mrow><mi mathvariant="normal">*</mi></mrow></msubsup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870545&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870541&type= 3.55599999 3.97933316 ， <math id="M55"><msubsup><mrow><mi>a</mi></mrow><mrow><mi>i</mi></mrow><mrow><mi mathvariant="normal">*</mi></mrow></msubsup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870551&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870548&type= 3.04800010 3.97933316 ， <math id="M56"><msubsup><mrow><mi>b</mi></mrow><mrow><mi>i</mi></mrow><mrow><mi mathvariant="normal">*</mi></mrow></msubsup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870557&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870554&type= 3.04800010 3.97933316 表示石窟表面不同风化类型及风化程度区域每个像素点的CIE1976L*a*b*值。 3.3　石窟表面风化多光谱成像数据主成分分析 2 对石窟表面多光谱成像数据进行降维，可以有效降低多光谱成像数据的冗余，提高智能评估方法的计算效率。光谱成像数据降维的常用方法主要分为非线性和线性降维，其中非线性降维以核主成分分析（Kernel Principal Component Analysis， KPCA） ［ 15 15 ］ 、核独立成分分析（Kernel Independent Component Analysis， KICA） ［ 16 16 ］ 为主，通过核函数完成非线性映射，以实现对非线性数据降维，同时保留了原始光谱数据的主要信息 ［ 17 17 ］ ；线性降维以主成分分析（Principal Component Analysis， PCA） ［ 18 18 ］ 、独立成分分析（Independent Component Analysis， ICA） ［ 19 19 ］ 为主，通过线性空间变换将高维数据空间投影到低维空间，即保留原始数据的大部分方差信息又消除了原始数据成分间的相关分量 ［ 20 20 ］ 。由于独立成分分析方法适用于数据的解混，且石窟表面风化光谱数据变化较为平缓，不同风化类型及风化程度间的反射光谱特性曲线线性可分，所以本文采用主成分分析方法对石窟表面多光谱成像数据进行线性降维处理。首先以统计计算的方式获得石窟表面风化多光谱成像数据中各波段间的相关性，从中提取特征向量和特征值，再以二者为基础进行正交变换。主成分分析方法对石窟表面风化多光谱成像数据处理的关键在于计算转换矩阵。通过转换矩阵将石窟表面风化多光谱成像数据从高维空间映射到低维空间，用较少的独立的主分量图像全面综合反映有用信息 ［ 21 21 ］ 。特征值的大小也反映着数据的变化程度，特征值越大说明该像素点的反射光谱响应数据变化程度越大，对于分类的评估结果影响程度越高，在多光谱成像数据中越重要。将石窟表面多光谱成像数据按相同像素点位置的相机响应值作为向量，以此建立石窟表面风化多光谱成像数据矩阵。原始的石窟表面风化多光谱成像数据由16个通道组成，每个通道的图片大小为1 276×1 028，将16个通道的光谱成像数据表示为16个大小为1×1 311 728的向量集合 <math id="M57"><mi>X</mi><mo>=</mo><mo stretchy="false">[</mo><msub><mrow><mi>x</mi></mrow><mrow><mn mathvariant="normal">1</mn></mrow></msub><mo>,</mo><msub><mrow><mi>x</mi></mrow><mrow><mn mathvariant="normal">2</mn></mrow></msub><mo>⋅</mo><mo>⋅</mo><mo>⋅</mo><msub><mrow><mi>x</mi></mrow><mrow><mi>i</mi></mrow></msub><msup><mrow><mo stretchy="false">]</mo></mrow><mrow><mi mathvariant="normal">T</mi></mrow></msup></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870564&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870559&type= 28.02466965 4.65666676 ，每个向量 <math id="M58"><msub><mrow><mi>x</mi></mrow><mrow><mi>i</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870571&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870568&type= 2.79399991 3.72533321 有1 311 728个反射光谱数据构成， <math id="M59"><msub><mrow><mi>x</mi></mrow><mrow><mi>i</mi></mrow></msub><mo>=</mo><mo stretchy="false">[</mo><msub><mrow><mi>x</mi></mrow><mrow><mi>i</mi><mn mathvariant="normal">1</mn></mrow></msub><mo>,</mo><msub><mrow><mi>x</mi></mrow><mrow><mi>i</mi><mn mathvariant="normal">2</mn></mrow></msub><mo>,</mo><mo>⋅</mo><mo>⋅</mo><mo>⋅</mo><mo>,</mo><msub><mrow><mi>x</mi></mrow><mrow><mi>i</mi><mi>n</mi></mrow></msub><mo stretchy="false">]</mo></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870576&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870573&type= 31.32666588 4.48733330 ，并对每一列光谱数据进行零均值化处理，即每一列的反射光谱数据减去该行数据的均值 <math id="M60"><mover><mrow><msub><mrow><mi>x</mi></mrow><mrow><mi>n</mi></mrow></msub></mrow><mo>¯</mo></mover></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870583&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870580&type= 3.21733332 4.23333359 ，然后计算 <math id="M61"><msub><mrow><mi>x</mi></mrow><mrow><mi>i</mi><mi>j</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870590&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870586&type= 3.30200005 3.72533321 的协方差矩阵 <math id="M62"><mi>C</mi><mi>o</mi><mi>v</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870595&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870592&type= 6.26533365 2.79399991 ，由协方差矩阵求其特征值 <math id="M63"><msub><mrow><mi>λ</mi></mrow><mrow><mi>n</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870601&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870599&type= 2.87866688 3.72533321 及其特征向量，计算公式为： <math display="block" id="M64"><mi>C</mi><mi>o</mi><mi>v</mi><mo>=</mo><mfrac><mrow><mn mathvariant="normal">1</mn></mrow><mrow><mi>i</mi><mo>-</mo><mn mathvariant="normal">1</mn></mrow></mfrac><mstyle displaystyle="true"><munderover><mo largeop="true">∑</mo><mrow><mi>i</mi><mo>=</mo><mn mathvariant="normal">1</mn></mrow><mrow><mn mathvariant="normal">16</mn></mrow></munderover></mstyle><mrow><mo stretchy="false">(</mo><msub><mrow><mi>x</mi></mrow><mrow><mi>i</mi><mi>j</mi></mrow></msub><mo>-</mo><mover><mrow><msub><mrow><mi>x</mi></mrow><mrow><mi>n</mi></mrow></msub></mrow><mo>¯</mo></mover><mo stretchy="false">)</mo><mo stretchy="false">(</mo><msub><mrow><mi>x</mi></mrow><mrow><mi>i</mi><mi>j</mi></mrow></msub><msup><mrow><mo>-</mo><mover><mrow><msub><mrow><mi>x</mi></mrow><mrow><mi>n</mi></mrow></msub></mrow><mo>¯</mo></mover><mo stretchy="false">)</mo></mrow><mrow><mi mathvariant="normal">T</mi></mrow></msup></mrow></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870606&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870603&type= 58.58933258 8.38199997 ， <math display="block" id="M65"><mi mathvariant="normal">d</mi><mi mathvariant="normal">e</mi><mi mathvariant="normal">t</mi><mo stretchy="false">(</mo><msub><mrow><mi>λ</mi></mrow><mrow><mi>n</mi></mrow></msub><mi>I</mi><mo>-</mo><mi>C</mi><mi>o</mi><mi>v</mi><mo stretchy="false">)</mo><mo>=</mo><mn mathvariant="normal">0</mn></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870613&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870609&type= 31.15733528 4.48733330 ， <math display="block" id="M66"><mo stretchy="false">(</mo><msub><mrow><mi>λ</mi></mrow><mrow><mi>n</mi></mrow></msub><mi>I</mi><mo>-</mo><mi>C</mi><mi>o</mi><mi>v</mi><mo stretchy="false">)</mo><msub><mrow><mi mathvariant="normal">A</mi></mrow><mrow><mi>n</mi></mrow></msub><mo>=</mo><mn mathvariant="normal">0</mn></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870619&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870615&type= 30.64933395 4.48733330 ， 其中：将特征值矩阵记为 <math id="M67"><mi>Λ</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870623&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870622&type= 2.53999996 2.79399991 ，由特征值所对应的特征向量组成的正交矩阵记为 <math id="M68"><msub><mrow><mi mathvariant="normal">A</mi></mrow><mrow><mi>n</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870629&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870626&type= 4.06400013 3.64066648 ，将特征值从大到小排列，将其中对应的前 n 个特征向量组成转换矩阵 <math id="M69"><mi>T</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870634&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870632&type= 2.53999996 2.79399991 ，将主成分集表示为 <math id="M70"><mi>Y</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870640&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870637&type= 2.87866688 2.79399991 ，则主成分提取公式为： <math display="block" id="M71"><mi>Y</mi><mo>=</mo><mi>T</mi><mi>X</mi></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870646&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870644&type= 13.12333298 2.79399991 . 4 实验及结果分析 1 4.1　实验部分 2 利用海洋光学公司的SpectroCamVIS型号多光谱成像系统对陕西省延安市清凉山万佛寺万佛窟表面风化岩体进行多光谱图像数据采集。对于复杂风化区域，波长采集范围400~940 nm，波段间隔20~60 nm，共采集16个波段，分辨率大小为1 276×1 028 pixel，如 图2 图2 所示。首先对清凉山万佛寺万佛窟表面风化多光谱图像数据进行处理，将大小为1 276×1 028×3的主分量特征数据与大小为1 276×1 028×1的光谱色差特征数据进行重组，重组后的特征数据大小为1 276×1 028×4；然后按像素点位置提取该像素点的光谱色差数据和主分量特征数据，建立大小为1 311 728×4的石窟表面风化待测样本数据库；最后对石窟表面不同风化类型及风化程度纯净风化区域进行标记，按标记位置提取光谱色差与主分量特征融合数据，并将该像素点所对应的不同风化类型及风化程度标签标记于特征数据矩阵最后一列，建立数据大小为34 860×5的训练及测试样本数据库。 为验证光谱色差与主分量融合数据对不同风化类型及风化程度表征的适用性及有效性，利用随机森林算法、K-近邻算法（K Near Neighbor， KNN）、BP神经网络（Back Propagation Neural Networks， BP）、径向基神经网络（Radial Basis Function， RBF）分别对光谱色差数据、主分量特征数据、光谱色差数据与主分量融合数据以及16通道原始反射光谱数据进行对比实验。为降低不同风化类型及风化程度数据间的不平衡对模型训练的影响，将融合表征数据库的10%用于预测，剩余部分采用5折交叉验证方法，轮流选取每种风化类型及风化程度数据的20%用于测试，80%用于训练，循环5次。 通过四种实验的预测准确率、Kappa系数及整体风化程度预测评估结果与该像素区域实际风化情况对比，对评估结果进行精度评价。石窟表面风化智能评估方法技术框图如 图3 图3 所示。 4 . 1 . 1　石窟表面风化光谱色差计算 3 为准确还原石窟表面不同风化区域的颜色信息，减少离散多光谱成像数据对不同风化类型及风化程度间光谱色差计算的影响，利用支持向量回归方法对石窟表面多光谱成像数据进行光谱反射率重建。不同风化类型及风化程度区域的光谱反射率重建后的光谱特性曲线如 图4 图4 所示。 对石窟表面不同风化区域的重建光谱特征数据进行颜色空间的线性变换，计算每个像素点的L*，a*，b*参数值；以基准点L*a*b*值为标准计算石窟表面不同风化类型及风化程度区域间的光谱色差。 在国际色差评价标准下，强盐析风化与弱盐析风化区域光谱色差值为5.83，颜色差异较为明显；强盐析风化与微盐析风化、积尘风化之间的光谱色差值均大于12，颜色差异非常大；弱盐析风化与微盐析风化之间的光谱色差值为6.85，颜色差异很大；弱盐析风化与积尘风化区域光谱色差值均大于12，颜色差异非常大；微盐析风化与积尘风化之间光谱色差值为8.33，颜色差异很大。不同风化类型及风化程度间的光谱色差计算结果如 表4 表4 所示。利用光谱色差特性可以对石窟表面不同风化类型及风化程度区域进行有效表征。 Color difference type <math id="M72"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">s</mi><mi mathvariant="normal">t</mi><mi mathvariant="normal">r</mi><mi mathvariant="normal">o</mi><mi mathvariant="normal">n</mi><mi mathvariant="normal">g</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870691&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870689&type= 7.87400007 3.21733332 <math id="M73"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">w</mi><mi mathvariant="normal">e</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">k</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870694&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870692&type= 7.19666624 3.21733332 <math id="M74"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">s</mi><mi mathvariant="normal">l</mi><mi mathvariant="normal">i</mi><mi mathvariant="normal">g</mi><mi mathvariant="normal">h</mi><mi mathvariant="normal">t</mi><mi mathvariant="normal">l</mi><mi mathvariant="normal">y</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870698&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870695&type= 8.63599968 3.21733332 <math id="M75"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">d</mi><mi mathvariant="normal">u</mi><mi mathvariant="normal">s</mi><mi mathvariant="normal">t</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870703&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870701&type= 6.43466663 3.21733332 <math id="M76"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">b</mi><mi mathvariant="normal">e</mi><mi mathvariant="normal">n</mi><mi mathvariant="normal">c</mi><mi mathvariant="normal">h</mi><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">r</mi><mi mathvariant="normal">k</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870706&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870704&type= 11.17599964 3.21733332 <math id="M77"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">s</mi><mi mathvariant="normal">t</mi><mi mathvariant="normal">r</mi><mi mathvariant="normal">o</mi><mi mathvariant="normal">n</mi><mi mathvariant="normal">g</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870691&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870689&type= 7.87400007 3.21733332 0 5.83 12.68 21.01 74.05 <math id="M78"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">w</mi><mi mathvariant="normal">e</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">k</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870694&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870692&type= 7.19666624 3.21733332 5.83 0 6.85 5.83 68.22 <math id="M79"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">s</mi><mi mathvariant="normal">l</mi><mi mathvariant="normal">i</mi><mi mathvariant="normal">g</mi><mi mathvariant="normal">h</mi><mi mathvariant="normal">t</mi><mi mathvariant="normal">l</mi><mi mathvariant="normal">y</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870698&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870695&type= 8.63599968 3.21733332 12.68 6.85 0 8.33 61.37 <math id="M80"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">d</mi><mi mathvariant="normal">u</mi><mi mathvariant="normal">s</mi><mi mathvariant="normal">t</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870703&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870701&type= 6.43466663 3.21733332 21.01 15.18 8.33 0 53.04 <math id="M81"><msub><mrow><mi>W</mi></mrow><mrow><mi mathvariant="normal">b</mi><mi mathvariant="normal">e</mi><mi mathvariant="normal">n</mi><mi mathvariant="normal">c</mi><mi mathvariant="normal">h</mi><mi mathvariant="normal">m</mi><mi mathvariant="normal">a</mi><mi mathvariant="normal">r</mi><mi mathvariant="normal">k</mi></mrow></msub></math> http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870706&type= http://html.publish.founderss.cn/rc-pub/api/common/picture?pictureId=21870704&type= 11.17599964 3.21733332 74.05 68.22 61.37 53.04 0 4 . 1 . 2　主分量数据特征提取 3 利用主成分分析方法对原始多光谱成像数据的16个波段图像进行降维处理。 图5 图5 表示主成分分析方法对石窟表面风化多光谱成像数据进行处理的各波段特征值分布情况，横轴代表波段值，纵轴代表特征值。 通过特征值数量随波段的变化结果可以看出，石窟表面风化多光谱成像的数据信息主要集中在第一波段和第二、三波段，其余波段所含的信息较少。 图6 图6 （a）表示石窟表面风化多光谱成像数据经过主成分分析第一主成分分量（PC1）图像， 图6 图6 （b）表示第二主成分分量（PC2）图像， 图6 图6 （c）表示第三主成分分量（PC3）图像， 图6 图6 （d）表示前三主成分分量合成的伪彩色图像。通过对比可以看出第一主成分具有的信息最多，第三主成分具有的信息较少，噪声占比上升。 前3个主分量贡献率及累计贡献率如 表5 表5 所示。第一主成分分量具有的图像信息十分丰富，其信息贡献率达到总信息量的95.41%，第二主成分分量和第三主成分分量贡献率较低分别为3.30%和0.53%，但包含了较多的石窟表面的岩性信息，有利于对不同风化类型及风化程度的区分识别。 Principal component Eigenvalues Contribution rate Cumulative contribution rate PC1 16757.90 95.41% 95.41% PC2 579.80 3.30% 98.71% PC3 92.21 0.53% 99.24% 前三个主成分分量累计贡献率为99.24%。利用主成分分析方法对多光谱成像数据进行处理，即降低了数据维度又保留了石窟表面风化反射光谱的主要信息。 4.2　石窟表面风化整体评估结果分析 2 将石窟表面风化光谱色差数据、主分量特征数据、反射率幅值表征数据、光谱色差与主分量特征融合表征数据分为训练、测试以及预测训练样本。通过RF、KNN、BP、RBF四种算法分别对不同风化类型及风化程度的石窟表面风化特征数据进行学习，并保存训练模型，通过预测准确率及Kappa系数对四种风化程度表征方法的有效性进行对比分析。 表6 表6 记录了四种评估方法对石窟表面不同风化类型及风化程度特征数据的分类准确率实验结果。 Method Accuracy evaluation Color difference - principal component Reflectance spectrum principal component Color difference 3 RF Training accuracy 100.00% 99.91% 98.52% 93.03% Prediction accuracy 99.86% 98.49% 91.34% 70.26% Kappa coefficient 0.99 0.98 0.89 0.53 3 KNN Training accuracy 99.79% 94.99% 88.60% 75.43% Prediction accuracy 97.50% 91.69% 82.68% 66.84% Kappa coefficient 0.96 0.90 0.69 0.42 3 BP Training accuracy 63.93% 63.19% 46.13% 36.84% Prediction accuracy 50.81% 48.00% 42.19% 32.51% Kappa coefficient 0.28 0.25 0.21 0.18 3 RBF Training accuracy 96.12% 95.95% 87.52% 70.08% Prediction accuracy 94.73% 92.88% 85.96% 65.84% Kappa coefficient 0.92 0.91 0.75 0.39 实验表明，在四种分类算法中光谱色差与主分量融合表征数据的预测准确率均较高，相比于光谱反射率幅值表征方法、主分量表征方法、色差表征方法，光谱色差与主分量特征融合的表征方法对石窟表面不同风化类型及风化程度的表征能力更强。其中在随机森林算法中的预测准确率最高，高达99.86%，Kappa系数为0.99，与评估结果较好的16通道反射光谱幅值表征方法相比，光谱色差与主分量特征融合表征方法使评估准确率提高了1.37%。利用训练好的随机森林算法模型，对石窟表面整体风化类型及风化程度分布情况进行评估，石窟表面风化多光谱成像数据中，待评估区域的总像素点数为1 311 728 points。在不同风化类型及风化程度表征数据下，随机森林算法对石窟表面每个像素点的预测评估结果如 图7 图7 所示。将不同风化类型及风化程度使用不同颜色进行可视化标记，红色表示强盐析风化、蓝色表示弱盐析风化、青色表示微盐析风化、绿色表示积尘风化。 其中光谱色差与主分量特征融合表征方法对石窟表面整体风化程度的评估结果为：强盐析风化区域像素点数为67 625 points（5.16%），弱盐析风化区域像素点数为368 921 points （28.12%），微盐析风化区域像素点数为359 457 points（27.40%），积尘风化区域像素点数为515 725 points（39.32%）。 对于图像中部，弱盐析风化与强盐析风化混杂区域，由于表面风化类型同为盐析风化，其物质状态特征相似，颜色差异较小，使得主分量特征难以达到精细化表征。 将光谱色差数据与主分量特征数据进行融合表征可以有效减少同种风化不同风化程度间的误划分，通过实验结果对比可以得出，光谱色差与主分量特征融合方法对石窟表面不同风化类型及风化程度的表征能力较强，其评估结果与实际风化状况最为接近。但对于图像底部弱盐析风化与积尘风化混杂区域，中部以弱盐析风化为主；对于底部突起部分的边缘，其表面物质主要为积尘，风化类型以积尘风化为主。由于此区域的风化类型及风化程度是由多种风化因素共同作用的结果，此区域上部为大面积的盐析风化区域，盐析风化产物随石窟表面的酥化脱落，覆盖于下部突起部分的表面，使得此区域的表面的物质状态特征与盐析风化区域相近，表面的反射光谱特性与弱盐析风化区域的反射光谱形态及幅值差异较小。利用随机森林算法对下部突起部分的风化类型及风化程度进行评估时，光谱色差数据根据颜色的差异将此部分划分为积尘风化，而主分量特征数据及其与光谱色差融合表征数据将大部分划分为弱盐析风化，光谱反射率幅值表征数据将此区域的风化类型进行了折中处理。对于颜色差异较大，物质状态特征相似的复杂风化区域，光谱色差特征难以达到精细化表征。主分量特征具有石窟表面物质光谱特性及部分岩性信息，可以减小不同风化类型及风化程度区域因颜色差异而带来的误划分，使得突起部分其表面风化类型主要以弱盐析风化为主，边缘以积尘风化为主，与实际表面物质状态及风化状况基本一致。光谱色差与主分量特征融合数据对于复杂风化区域不同风化类型及风化程度区域，可以达到精细化表征。 图7 图7 （d）为只用主分量特征数据进行训练，对石窟表面风化类型及风化程度整体表征的实验结果。对于弱盐析风化与微盐析风化、弱盐析风化与强盐析风化混杂区域，由于石窟表面物质状态特征相似，主分量特征对同种风化类型不同风化程度的表征能力较弱，使得评估结果与实际风化状况差异较大。 图7 图7 （e）为光谱色差对不同风化类型及风化程度评估结果，由于色差数据较为单一，不能有效表征石窟表面风化区域的物质状态特征。在通过随机森林算法进行分类时，弱盐析风化与积尘风化混杂区域、微盐析风化与积尘风化混杂区域，出现大面积的误划分，石窟表面整体风化程度预测的准确率降低。16通道的光谱反射率幅值特征数据对不同风化类型及风化程度的表征能力虽然较好，但强盐析风化与弱盐析风化混杂区域、弱盐析风化与积尘风化混杂区域的风化状况与实际风化状况及光谱色差与主分量融合表征评估结果对比分析表明，光谱色差与主分量融合的特征数据对风化状况的表征能力更强，对复杂风化区域的评估结果更加精细、准确。 5 结　论 1 本文提出了一种基于光谱色差与主分量特征融合的石窟表面风化智能量化评估方法。实验表明，光谱色差数据可以对同风化类型不同风化程度区域进行有效表征；主分量特征可以有效反映石窟表面风化的物质状态特征及部分岩性信息。融合光谱色差与主分量特征对石窟表面复杂风化区域的风化类型及风化程度进行表征的方法，有效解决了对于复杂风化区域，单一光谱特征表征能力弱的问题，使得随机森林算法对石窟表面整体风化程度评估的准确率提高了1.37%。光谱色差与主分量融合数据方法对于石窟表面复杂风化区域的风化类型及风化程度可以达到精细化表征，其评估结果更加准确。 参考文献 QIN Y ， WANG Y H ， LI L L ， et al . Experimental weathering of weak sandstone without direct water participation by using sandstone from the Yungang grottoes in Datong， China ［J］. Rock Mechanics and Rock Engineering ， 2016 ， 49 （ 11 ）： 4473 - 4478 . doi: 10.1007/s00603-016-1003-3 http://dx.doi.org/10.1007/s00603-016-1003-3 Experimental weathering of weak sandstone without direct water participation by using sandstone from the Yungang grottoes in Datong， China HE X ， LIU Y J ， BECKETT P ， et al . Hybrid color filters for multispectral imaging ［J］. Advanced Theory and Simulations ， 2020 ， 3 （ 11 ）： 2000137 . doi: 10.1002/adts.202000137 http://dx.doi.org/10.1002/adts.202000137 Hybrid color filters for multispectral imaging 易定容 ， 孔令华 ， 赵艳丽 ， 等 . 窄带快照式多光谱成像色彩还原方法 ［J］. 光谱学与光谱分析 ， 2021 ， 41 （ 1 ）： 183 - 187 . 窄带快照式多光谱成像色彩还原方法 YI D R ， KONG L H ， ZHAO Y L ， et al . Color recovery method for snapshot narrow band spectral imaging technology ［J］. Spectroscopy and Spectral Analysis ， 2021 ， 41 （ 1 ）： 183 - 187 . （in Chinese） Color recovery method for snapshot narrow band spectral imaging technology 高泽东 ， 高洪兴 ， 朱院院 ， 等 . 快照式光谱成像技术综述 ［J］. 光学 精密工程 ， 2020 ， 28 （ 6 ）： 1323 - 1343 . doi: 10.3788/ope.20202806.1323 http://dx.doi.org/10.3788/ope.20202806.1323 快照式光谱成像技术综述 GAO Z D ， GAO H X ， ZHU Y Y ， et al . Review of snapshot spectral imaging technologies ［J］. Opt. Precision Eng. ， 2020 ， 28 （ 6 ）： 1323 - 1343 . （in Chinese） . doi: 10.3788/ope.20202806.1323 http://dx.doi.org/10.3788/ope.20202806.1323 Review of snapshot spectral imaging technologies BRUNETTI B ， MILIANI C ， ROSI F ， et al . Non-invasive investigations of paintings by portable instrumentation： the MOLAB experience ［J］. Topics in Current Chemistry ， 2016 ， 374 （ 1 ）： 1 - 35 . doi: 10.1007/s41061-015-0008-9 http://dx.doi.org/10.1007/s41061-015-0008-9 Non-invasive investigations of paintings by portable instrumentation： the MOLAB experience DE MANINCOR N ， MARCHIORO G ， FIORIN E ， et al . Integration of multispectral visible-infrared imaging and pointwise X-ray fluorescence data for the analysis of a large canvas painting by Carpaccio ［J］. Microchemical Journal ， 2020 ， 153 ： 104469 . doi: 10.1016/j.microc.2019.104469 http://dx.doi.org/10.1016/j.microc.2019.104469 Integration of multispectral visible-infrared imaging and pointwise X-ray fluorescence data for the analysis of a large canvas painting by Carpaccio 曹赤鹏 ， 王慧琴 ， 王可 ， 等 . 基于多光谱成像和随机森林算法的石窟表面风化智能评估方法 ［J］. 光学学报 ， 2020 ， 40 （ 22 ）： 201 - 211 . doi: 10.3788/AOS202040.2230001 http://dx.doi.org/10.3788/AOS202040.2230001 基于多光谱成像和随机森林算法的石窟表面风化智能评估方法 CAO CH P ， WANG H Q ， WANG K ， et al . Intelligent evaluation method of grottoes surface weathering based on multispectral imaging and random forest algorithm ［J］. Acta Optica Sinica ， 2020 ， 40 （ 22 ）： 201 - 211 . （in Chinese） . doi: 10.3788/AOS202040.2230001 http://dx.doi.org/10.3788/AOS202040.2230001 Intelligent evaluation method of grottoes surface weathering based on multispectral imaging and random forest algorithm 蒋晓东 ， 曹建劲 ， 李易安 ， 等 . 近红外光谱技术在红砂岩文物风化研究中的应用 ［J］. 光谱学与光谱分析 ， 2011 ， 31 （ 8 ）： 2102 - 2105 . doi: 10.3964/j.issn.1000-0593(2011)08-2102-04 http://dx.doi.org/10.3964/j.issn.1000-0593(2011)08-2102-04 近红外光谱技术在红砂岩文物风化研究中的应用 JIANG X D ， CAO J J ， LI Y A ， et al . Application of near-infrared spectrum technology to research of weathering of red sandstone relics ［J］. Spectroscopy and Spectral Analysis ， 2011 ， 31 （ 8 ）： 2102 - 2105 . （in Chinese） . doi: 10.3964/j.issn.1000-0593(2011)08-2102-04 http://dx.doi.org/10.3964/j.issn.1000-0593(2011)08-2102-04 Application of near-infrared spectrum technology to research of weathering of red sandstone relics 柴勃隆 ， 范宇权 ， 苏伯民 ， 等 . 莫高窟194窟漫漶壁画多光谱成像调查与研究 ［J］. 敦煌研究 ， 2010 （ 6 ）： 74 - 76 . doi: 10.3969/j.issn.1000-4106.2010.06.014 http://dx.doi.org/10.3969/j.issn.1000-4106.2010.06.014 莫高窟194窟漫漶壁画多光谱成像调查与研究 CHAI B L ， FAN Y Q ， SU B M ， et al . Investigating and studying the fading away wall-paintings of Cave 194 at Mogao Grottoes with the UV-induced Fluorescence Photography ［J］. Dunhuang Research ， 2010 （ 6 ）： 74 - 76 . （in Chinese） . doi: 10.3969/j.issn.1000-4106.2010.06.014 http://dx.doi.org/10.3969/j.issn.1000-4106.2010.06.014 Investigating and studying the fading away wall-paintings of Cave 194 at Mogao Grottoes with the UV-induced Fluorescence Photography 周霄 ， 高峰 ， 张爱武 ， 等 . VIS/NIR高光谱成像在中国云冈石窟砂岩风化状况分布研究中的进展 ［J］. 光谱学与光谱分析 ， 2012 ， 32 （ 3 ）： 790 - 794 . doi: 10.3964/j.issn.1000-0593(2012)03-0790-05 http://dx.doi.org/10.3964/j.issn.1000-0593(2012)03-0790-05 VIS/NIR高光谱成像在中国云冈石窟砂岩风化状况分布研究中的进展 ZHOU X ， GAO F ， ZHANG A W ， et al . Advance in the study of the powdered weathering profile of sandstone on China Yungang grottoes based on VIS/NIR hyperspectral imaging ［J］. Spectroscopy and Spectral Analysis ， 2012 ， 32 （ 3 ）： 790 - 794 . （in Chinese） . doi: 10.3964/j.issn.1000-0593(2012)03-0790-05 http://dx.doi.org/10.3964/j.issn.1000-0593(2012)03-0790-05 Advance in the study of the powdered weathering profile of sandstone on China Yungang grottoes based on VIS/NIR hyperspectral imaging 李杰 ， 李晓龙 ， 唐秋华 ， 等 . 基于PCA的时间分辨油荧光光谱分析及优化 ［J］. 光学 精密工程 ， 2017 ， 25 （ 4 ）： 352 - 358 . doi: 10.3788/OPE.20172504.0884 http://dx.doi.org/10.3788/OPE.20172504.0884 基于PCA的时间分辨油荧光光谱分析及优化 LI J ， LI X L ， TANG Q H ， et al . Spectral analysis and optimization of time-resolved oil fluorescence based on PCA ［J］. Opt. Precision Eng. ， 2017 ， 25 （ 4 ）： 352 - 358 . （in Chinese） . doi: 10.3788/OPE.20172504.0884 http://dx.doi.org/10.3788/OPE.20172504.0884 Spectral analysis and optimization of time-resolved oil fluorescence based on PCA 屈东胜 ， 刘连伟 ， 樊宏杰 ， 等 . 基于主成分分析方法的光电系统抗干扰性能评估研究 ［J］. 激光与红外 ， 2020 ， 50 （ 6 ）： 712 - 717 . doi: 10.3969/j.issn.1001-5078.2020.06.011 http://dx.doi.org/10.3969/j.issn.1001-5078.2020.06.011 基于主成分分析方法的光电系统抗干扰性能评估研究 QU D SH ， LIU L W ， FAN H J ， et al . Evaluation method of infrared system countermeasures performance based on principal components analysis ［J］. Laser & Infrared ， 2020 ， 50 （ 6 ）： 712 - 717 . （in Chinese） . doi: 10.3969/j.issn.1001-5078.2020.06.011 http://dx.doi.org/10.3969/j.issn.1001-5078.2020.06.011 Evaluation method of infrared system countermeasures performance based on principal components analysis 王来贵 ， 丁盛鹏 ， 何慧娟 ， 等 . 冻融循环作用下含结核砂岩风化特征实验研究 ［J］. 工程地质学报 ， 2018 ， 26 （ 3 ）： 611 - 619 . 冻融循环作用下含结核砂岩风化特征实验研究 WANG L G ， DING SH P ， HE H J ， et al . Experimental study on cracking characteristics of concretion sandstone with freezing-thawing cycles ［J］. Journal of Engineering Geology ， 2018 ， 26 （ 3 ）： 611 - 619 . （in Chinese） Experimental study on cracking characteristics of concretion sandstone with freezing-thawing cycles 赵丽娟 ， 王慧琴 ， 王可 ， 等 . 基于多核支持向量回归的光谱反射率重建方法 ［J］. 液晶与显示 ， 2018 ， 33 （ 12 ）： 1008 - 1018 . doi: 10.3788/YJYXS20183312.1008 http://dx.doi.org/10.3788/YJYXS20183312.1008 基于多核支持向量回归的光谱反射率重建方法 ZHAO L J ， WANG H Q ， WANG K ， et al . Spectral reflectance reconstruction based on multi-kernel support vector regression ［J］. Chinese Journal of Liquid Crystals and Displays ， 2018 ， 33 （ 12 ）： 1008 - 1018 . （in Chinese） . doi: 10.3788/YJYXS20183312.1008 http://dx.doi.org/10.3788/YJYXS20183312.1008 Spectral reflectance reconstruction based on multi-kernel support vector regression GAN L ， XIA J S ， DU P J ， et al . Multiple feature kernel sparse representation classifier for hyperspectral imagery ［J］. IEEE Transactions on Geoscience and Remote Sensing ， 2018 ， 56 （ 9 ）： 5343 - 5356 . doi: 10.1109/tgrs.2018.2814781 http://dx.doi.org/10.1109/tgrs.2018.2814781 Multiple feature kernel sparse representation classifier for hyperspectral imagery LI J Y ， LI S M ， CHEN X H ， et al . The hybrid KICA-GDA-LSSVM method research on rolling bearing fault feature extraction and classification ［J］. Shock and Vibration ， 2015 ， 2015 ： 1 - 9 . doi: 10.1155/2015/512163 http://dx.doi.org/10.1155/2015/512163 The hybrid KICA-GDA-LSSVM method research on rolling bearing fault feature extraction and classification 王卓薇 ， 罗鉴鹏 ， 李学识 ， 等 . 结合径向基函数和KPCA的食用油太赫兹光谱特征提取方法 ［J］. 光谱学与光谱分析 ， 2020 ， 40 （ 2 ）： 391 - 396 . 结合径向基函数和KPCA的食用油太赫兹光谱特征提取方法 WANG ZH W ， LUO J P ， LI X SH ， et al . Edible oil terahertz spectral feature extraction method combining radial basis function and KPCA ［J］. Spectroscopy and Spectral Analysis ， 2020 ， 40 （ 2 ）： 391 - 396 . （in Chinese） Edible oil terahertz spectral feature extraction method combining radial basis function and KPCA WOLD S ， ESBENSEN K ， GELADI P . Principal component analysis ［J］. Chemometrics and Intelligent Laboratory Systems ， 1987 ， 2 （ 1/2/3 ）： 37 - 52 . doi: 10.1016/0169-7439(87)80084-9 http://dx.doi.org/10.1016/0169-7439(87)80084-9 Principal component analysis XIAO F L ， LU C ， ZHAO D J ， et al . Independent component analysis and graph theoretical analysis in patients with narcolepsy ［J］. Neuroscience Bulletin ， 2019 ， 35 （ 4 ）： 743 - 755 . doi: 10.1007/s12264-018-0307-6 http://dx.doi.org/10.1007/s12264-018-0307-6 Independent component analysis and graph theoretical analysis in patients with narcolepsy 唐勇波 ， 桂卫华 ， 彭涛 ， 等 . PCA和KICA特征提取的变压器故障诊断模型 ［J］. 高电压技术 ， 2014 ， 40 （ 2 ）： 557 - 563 . doi: 10.13336/j.1003-6520.hve.2014.02.031 http://dx.doi.org/10.13336/j.1003-6520.hve.2014.02.031 PCA和KICA特征提取的变压器故障诊断模型 TANG Y B ， GUI W H ， PENG T ， et al . Transformer fault diagnosis model based on PCA and KICA feature extraction ［J］. High Voltage Engineering ， 2014 ， 40 （ 2 ）： 557 - 563 . （in Chinese） . doi: 10.13336/j.1003-6520.hve.2014.02.031 http://dx.doi.org/10.13336/j.1003-6520.hve.2014.02.031 Transformer fault diagnosis model based on PCA and KICA feature extraction 廖一鹏 ， 陈诗媛 ， 杨洁洁 ， 等 . NSST域改进ORB的泡沫流动特征提取及加药状态识别 ［J］. 光学 精密工程 ， 2020 ， 28 （ 12 ）： 2684 - 2699 . doi: 10.37188/OPE.20202812.2684 http://dx.doi.org/10.37188/OPE.20202812.2684 NSST域改进ORB的泡沫流动特征提取及加药状态识别 LIAO Y P ， CHEN SH Y ， YANG J J ， et al . Dosing status identification and froth flow feature extraction based on improved ORB in NSST domain ［J］. Opt. Precision Eng. ， 2020 ， 28 （ 12 ）： 2684 - 2699 . （in Chinese） . doi: 10.37188/OPE.20202812.2684 http://dx.doi.org/10.37188/OPE.20202812.2684 Dosing status identification and froth flow feature extraction based on improved ORB in NSST domain 10.37188/OPE.20212910.2444 O433.4 A 1004-924X（2021）10-2444-12 曹赤鹏,王慧琴,王可等.光谱色差与主分量特征融合的石窟表面风化智能评估[J].光学精密工程,2021,29(10):2444-2455. 曹赤鹏,王慧琴,王可等.光谱色差与主分量特征融合的石窟表面风化智能评估[J].光学精密工程,2021,29(10):2444-2455. CAO Chi-peng,WANG Hui-qin,WANG Ke,et al.Intelligent evaluation of grotto surface weathering based on spectral chromatic aberration and principal component feature fusion[J].Optics and Precision Engineering,2021,29(10):2444-2455. CAO Chi-peng,WANG Hui-qin,WANG Ke,et al.Intelligent evaluation of grotto surface weathering based on spectral chromatic aberration and principal component feature fusion[J].Optics and Precision Engineering,2021,29(10):2444-2455. 国家重点研发计划项目（No.2019YFC1520500）；陕西省自然科学基金项目（No.2021JM-377）；陕西省科技厅科技合作项目（No.2020KW-012）；陕西省教育厅智库项目（No.18JT006）；西安市科技局高校人才服务企业项目（No.GXYD10.1） 2021-02-09 2021-03-08 2021-10-15 2021-11-11 光学精密工程 10 29