0   引言

在地面试验过程中, 超燃冲压发动机的主动冷却燃烧室燃烧气流的加热导致壁面结构温度急剧升高^[1], 其最高工作温度超过1000 ℃ ^[2-3], 壁面结构在热力载荷下会发生变形而产生热应力, 热应力过大将导致结构损坏. 因此, 燃烧室结构高温变形测量对燃烧室热变形匹配设计和结构强度评估十分关键.

燃烧室结构变形测量有接触式和非接触式2种方法. 应变片等接触式测量方法虽能保证单点测量精度, 但在高温条件下易脱落^[4-7], 无法满足燃烧室工作要求. 数字图像相关法(Digital Image Correlation, DIC) 是一种非接触变形测量方法, 其以散斑图案作为特征信息, 通过相似度函数运算并按照特定搜索路径匹配出变形前后特征点坐标位置，具有光路简单、抗干扰强、视场广、精度高等特点^[8-9]，在工程结构变形测量中得到应用^[10-12]，用于测量燃烧室结构高温变形有一定优势。但温度变化引起的像素灰度值变化可能会对匹配运算产生干扰，因此需要进一步分析。

对于燃烧室结构总变形测量，由于特征点像素位移常常过大，甚至超出数字图像相关搜索区域，引起误匹配，故需要改进数字图像相关法。加速鲁棒特征算法(Speeded Up Robust Features, SURF)具有很好的鲁棒性和尺度不变性，能够在光照改变和图像发生一定变化的情况下准确提取并匹配特征点^[13-16]，适合高温环境下结构变形初步匹配。

本文采用非接触方法开展了点火工作条件下燃烧室结构的图像采集，并利用SURF和改进的圆柱曲面数字图像相关法进行了结构总变形和局部应变计算及分析。

1   改进的数字图像相关法

1.1   DIC基本原理

数字图像相关法以物体表面的散斑纹理图案为特征，采集变形前后2幅图像进行对比，匹配图像子区相似度，以此进行特征点位移计算。常用的变形子区形状描述方式有刚性变形函数、一阶变形函数、二阶变形函数和高阶变形函数。针对燃烧室局部应变测量区域变形类型，采用一阶变形函数较为合理^[17]。在参考图像中选取一定大小的矩形参考子区，其中心为P(x₀, y₀)，依据一阶变形函数，变形图像必然存在一与参考子区相似度最高的子区图像，记为变形子区，其中心点坐标为P′(x′₀, y′₀)。参考子区和变形子区的相似度通常采用相关函数来度量，其基本原理如图 1所示。

图  1  变形前后的图像子区

Fig.  1  Subset region of image before and after deformation

下载: 全尺寸图片 幻灯片

采用零均值归一化互相关函数评价变形前后两子区的相似度^[18-19]，该函数在光照微弱变化环境中可降低干扰因素影响，其表达式为：

式中, r为图像子区半径，f(x, y)为参考图像中点(x, y)处的灰度值；g(x′, y′)为变形图像中点(x′, y′)处的灰度值；f_m和g_m为对应子区的灰度平均值。变形子区和参考子区的匹配度与相关系数C值正相关，C绝对值在[0, 1]之间，C=1代表二者完全匹配。

应变测量中可通过双三次样条插值获取亚像素位置像素值，并结合曲面拟合算法计算出亚像素位移结果，以提高应变测量精度。

1.2   大刚性位移微应变散斑匹配方法

采用高分辨率CCD相机采集燃烧室表面局部区域散斑图像，可提高应变测量精度，使燃烧室表面微小应变位移对应图像中的大像素位移; 但因此使得特征点位移区域过大，传统数字图像相关法搜索区域过大时的误匹配率大幅升高。限制数字图像相关法搜索范围需进行大刚性位移计算，以此结果代入数字图像相关法作为位移初值，再在该位置进行小范围数字图像相关搜索以得到准确的亚像素位移值。散斑包含大量黑白对比分明的角点特征，本文使用高效的SURF方法进行散斑图案角点特征匹配以求取大刚性位移。特征点匹配过程主要分3个阶段。

1) 特征点检测。SURF采用Hessian矩阵检测特征点。对图像I中的某一点P(x, y)，经过高斯滤波后具有尺度无关性的Hessian矩阵定义为：

式中，x为图像中特征点坐标，δ为特征点所在的空间尺度，L_xx(x, δ)、L_xy(x, δ)和L_yy(x, δ)为图像I中点(x, y)与高斯二阶偏导数的卷积。

实际运算中为降低计算时间，SURF算法使用盒式滤波器近似代替高斯二阶导数^[20]，用积分图像加速卷积，构造了一种Fast-Hessian矩阵，则Hessian矩阵行列式近似表示为:

式中，D_xx(x, δ)、D_xy(x, δ)和D_yy(x, δ)为盒式滤波器与图像中点(x, y)的卷积；w为权重系数，实际应用中常采用经验值0.9^[20]。

为初步定位特征点，遍历整幅图像各个尺度空间，对每个点的Hessian行列式的值与其周围3×3×3空间内其他26个点的值进行比较，采用非极大值抑制搜索极值点，并对局部极值点在尺度空间和图像空间中作插值计算^[21]，获得初步特征点位置及尺度值。

2) 特征点描述符。使用SURF确定每个特征点的主方向，使特征点描述符具有旋转不变性。以特征点为中心，构造半径为6δ的圆形区域，计算区域内各点沿x和y方向的Harr小波响应dx、dy，赋予响应值以不同的高斯权重系数，越靠近特征点的响应权重越大。将60°内的x和y方向Harr小波响应相加形成一个局部方向向量，遍历整个圆形区域，选择最长向量方向为该特征点的主方向。

以特征点为中心，沿主方向构造边长20δ的方形区域，并将其分为4×4个子区域，每个子区域中取25个采样点，在每个采样点分别求x和y方向的Harr小波响应dx、dy。在每一子区域中分别对Harr小波在x和y方向的响应dx、dy、|dx|、|dy|求和，得到1个四维向量ν=(∑dx, ∑dy, ∑| dx |, ∑|dy|)。而4×4个子区域里面都有1个四维向量，因此一共得到1个4×4×4的64维向量，这就是SURF算法的特征点描述符。Harr小波响应本身具有亮度不变性，通过特征矢量的归一化处理可实现对比度不变性。

3) 特征点匹配。特征点匹配是根据其描述符相似性与几何约束来进行的。获取SURF特征点后，优先采用k-d树在欧式空间搜索查找每个特征点的2个最近邻特征点。在这2个特征点中，如果最近距离与次近距离的比值小于某个比例阈值则接受这对匹配点。采用的欧式距离公式为:

式中，(x₁, x₂, …, x₆₄)、(x′₁, x′₂, …, x′₆₄)为待匹配的2个特征点的SURF特征向量。

在实际大刚性位移计算中，将待测区域划分若干大小均匀的网格区域，对待测点所在的网格子区进行SURF特征点匹配，提取匹配度最高的特征点计算变形前后的位移量，以此作为变形初值代入数字图像相关法进行最终的亚像素精确定位。

2   单目相机局部曲面测量校正

采用单目高分辨率CCD相机测量燃烧室局部圆柱曲面应变时，图像结果为平面图像，与实际曲面存在差异，需对其进行校正。

已知被测圆柱外表面半径R=100 mm，试验前固定相机安装位置，调整光轴正对被测圆柱中心，视距为480 mm。对相机进行像素当量标定，得到图像1个像素与实际长度的比例关系k₀。

如图 2所示，相机测量视场覆盖曲面MNKH，圆弧在所选曲面区域内。根据相机成像原理，圆弧在图像平面内为直线AB，曲面信息在成像平面内发生畸变成为平面M′N′K′H′。对测量区域内圆柱曲面进行校正：直线AB对应的像素长度即相机竖直方向像素高度L_pv=3072 pixel，其实际长度为L_v。取直线AB上任一点P，点P和相机光轴原点G连线相交圆弧于点P₀，其为点P在圆弧段对应的点，过点P作直线AB的垂线交圆弧于P′，此处以P′近似代替点P₀。

图  2  单目测量曲面校正

Fig.  2  Adjust cylindrical surface by single camera

下载: 全尺寸图片 幻灯片

在图像平面内，AB两点的横坐标与点P相同，其纵坐标表示采集图像的最上端和最下端，记A(x, 0) 和B(x, J)，J为局部应变测量相机竖直方向像素个数。则P(x, y)到u轴的像素距离d和实际距离d′分别为：

P′至u轴的实际距离为d′，结合圆柱曲面半径R可求得OP′相对u轴的夹角为：

P(x, y)对应曲面的实际坐标为P′(x, R, θ_P)。通过测量变形前后图像上任意两点的像素坐标值，即可计算其对应的平均应变。

综上，曲面应变计算公式为：

式中，ε为平均应变，(x₁, y₁)和(x₂, y₂)为变形前图像坐标系中两点的坐标值，(x₁′, y₁′)和(x₂′, y₂′)为变形后图像坐标系中两点的坐标值，θ_y1、θ_y2、θ_y1′、θ_y2′分别为两点变形前后的角度值。

3   燃烧室高温结构变形测量及分析

3.1   测量系统简介

为验证本文非接触测量方法在燃烧室高温工作环境下结构变形测量的可靠性，以某超燃冲压发动机主动冷却燃烧室为测量对象，开展结构高温变形非接触测量试验与计算研究。

燃烧室采用高温合金材料，总长约1.6 m，燃料为RP-3航空煤油。燃烧室点火工作时间20~100 s。采用点焊机在外壁面焊接K型热偶丝测量外壁温，沿轴向共布置15个测点。采用三维传热分析工具计算的燃烧室沿程外壁温为300~800 K。

燃烧室结构变形非接触测量试验系统主要由燃烧室模型、模型支架、相机支架、工业相机和工控机等组成，系统示意图如图 3所示，系统现场安装如图 4所示。燃烧室安装平台和光学实验平台均采用厚重稳固支架，且二者相互隔离。燃烧室安装过程中，在竖直和左右两侧方向进行固定，采用轴向滑移槽使燃烧室受热膨胀时能沿轴向移动。光学实验平台固定在地基上，并采用隔震垫隔离地基的振动。

图  3  测量系统安装示意图

Fig.  3  Schematic diagram of the measuring system

下载: 全尺寸图片 幻灯片

图  4  测量系统现场安装照片

Fig.  4  Photograph of the measuring system

下载: 全尺寸图片 幻灯片

采用一台2500万像素的Baumer的LXG250M型工业相机(#1)测量燃烧室左右两端轴向总变形，帧频为5帧/s。在燃烧室的入口和出口法兰分别制作十字特征标识。

安装#1相机时，调整其光轴角度，视场覆盖燃烧室左右两端且预留一定变形冗余空间。

采用一台1200万像素的Baumer的LXG120M型工业相机(#2)测量燃烧室中部圆柱表面局部应变，帧频为10帧/s。针对燃烧室中部圆弧曲面结构，设计了多自由度的相机支架，具有3个安装槽，本次试验中采用图 3中的安装槽②，相机与燃烧室壁面垂直距离为0.5 m。在局部待测量区域喷涂高温底漆，使用氧化铝粉末制作高温散斑图案，并制作若干特征标识。底漆和散斑厚度分别小于0.2和0.1 mm，在1000 K高温下不会产生脱落。

安装#2相机时，使其光轴正对待测圆柱中心线，根据燃烧室外部空间结构，选择视场大小为100 mm ×100 mm。架设工业聚光灯作为辅助光源，调整相机的曝光值至300 μs，从而降低燃烧室低频振动(约10 Hz)对图像采集的影响。

3.2   结构变形图像分析方法

燃烧室结构总变形采集图像分析方法：采用SURF方法跟踪全局序列图像中燃烧室两端特征标识的位移量，出口端轴向位移差值减去入口端轴向位移差值即为燃烧室结构轴向总变形量。

燃烧室局部应变采集图像分析方法：采用改进的圆柱曲面数字图像相关法对序列图像中的散斑图像和特征标识进行位移跟踪，两测点欧氏距离变化值和初始欧氏距离之比即为平均应变量。

3.3   结构高温变形测量与结果分析

燃烧室点火工作试验流程：1) 系统准备，相机和光源电源启动、采用网线将相机与工控机连接、启动测量软件并设置参数；2) 风洞系统开始试验前5 s，触发测量系统开始采集图像，以时间为序保存每帧图像；3) 试验开始，实时采集并存储燃烧室结构图像；4) 试验结束5 s后，触发相机停止采集图像；5) 存储数据，关闭测量系统，采用上述分析方法对采集图像进行结果分析。

试验开始前采集的燃烧室入口和出口法兰图像如图 5所示, 试验过程中采集的燃烧室入口和出口法兰图像如图 6所示。

图  5  试验前燃烧室入口和出口法兰采集图像

Fig.  5  Photograph of inlet and outlet before test

下载: 全尺寸图片 幻灯片

图  6  试验过程中燃烧室入口和出口法兰采集图像

Fig.  6  Photograph of inlet and outlet during test

下载: 全尺寸图片 幻灯片

根据采集的燃烧室入口和出口法兰图像，采用SURF方法分析得到的结构轴向总变形曲线如图 7所示。

图  7  结构总变形曲线

Fig.  7  Curve of axial total deformation

下载: 全尺寸图片 幻灯片

由图 7可知，燃烧室工作开始后，结构轴向总变形量随时间逐渐增大，在21 s达到峰值约为6.6 mm，在26 s后基本维持稳定值约为4.8 mm。

主动冷却燃烧室点火工作过程中，外壁面温度逐渐升高，最终达到热平衡状态。在稳定工作状态下，燃烧室外壁温沿轴向分布为316~820 K。利用燃烧室测量的沿程平均外壁温升，采用工程算法估算的燃烧室结构总变形量约为5.8 mm。可见，稳定工作状态下，燃烧室结构总变形的测量结果与工程估算结果比较吻合。

试验开始前采集的燃烧室局部图像如图 8所示，试验过程中采集的燃烧室局部图像如图 9所示。

图  8  试验前燃烧室局部采集图像

Fig.  8  Photograph of middle position before test

下载: 全尺寸图片 幻灯片

图  9  试验过程中燃烧室局部采集图像

Fig.  9  Photograph of middle position during test

下载: 全尺寸图片 幻灯片

根据燃烧室局部应变测量相机采集的图像，采用改进的圆柱曲面数字图像相关法计算得到的局部位移场和应变曲线分别如图 10和11所示。

图  10  位移场分布

Fig.  10  Displacement field distribution

下载: 全尺寸图片 幻灯片

图  11  局部应变曲线

Fig.  11  Curve of partial strain

下载: 全尺寸图片 幻灯片

由图 11可知，燃烧室结构局部应变曲线和轴向总变形曲线变化趋势一致。燃烧室开始工作后，局部应变随时间逐渐增大，在21 s达到峰值约为0.0075，在26 s后基本维持稳定值约为0.0049。

本次试验中，稳定工作状态下燃烧室局部应变测量部位的中心点温度约为600 K(该区域需要拍摄图像而无法布置测点，故假定该处与圆周对称部位的外壁温测点温度相同)。假定局部应变测量区域温度均匀且和该测点温度一致，采用工程算法估算的燃烧室结构局部平均应变约为0.0043。可见，稳定工作状态下，燃烧室结构局部平均应变测量结果与工程估算结果吻合较好。

4   结论

基于加速鲁棒特征算法和改进的圆柱曲面数字图像相关法，建立了燃烧室结构高温变形非接触测量方法，完成了主动冷却燃烧室在高温工作环境下结构总变形和局部应变测量及分析。

燃烧室结构高温变形测量结果表明，稳定工作条件下，结构轴向总变形约4.8 mm、局部平均应变约0.0049。测量结果与工程估算结果吻合较好，说明本文的非接触测量方法有效，能够支撑燃烧室热结构设计，测量数据能够用于三维结构强度数值计算的验证。

后续将以本文算法为基础，采用多目相机研究燃烧室外壁面结构三维变形，并使用更高像素的相机和增加测点数量的方法进一步提高测量和计算精度。