阶梯变迎角试验方式不仅风洞试验时间长,而且测量的模型迎角数量有限,提出连续变迎角试验方法[4-8],已经成为国际上常规测力试验的重要手段。前期连续变迎角试验显示该方法可节省常规高超声速测力试验时间约20%~40%,但试验过程中模型存在一定的谐振,谐振产生的惯性力会叠加到天平中,引入天平的测量误差[9]。
去除该测量误差常用的方法是采用加速度补偿[10-12],即在模型中特定位置安装加速度计测量振动加速度,基于加速度扣除此偏差。该方法对1~4 Hz的低频噪声处理效果不理想;此外,高超声速风洞试验模型一般较小,受模型安装空间、引线等限制,安装加速度计等额外传感器比较困难,故基于数据本身进行滤波补偿,提高测量精度具有重要意义。
文献[13-14]采用小波阈值函数及小波工具箱分析风洞连续信号,文献[15]采用基于db42小波的小波降噪方法处理连续变迎角信号。上述方法可以有效滤除噪声,但没有针对该低频振动信号的处理分析,直接针对连续变迎角数据中低频信号处理的文献较少。本文基于连续变迎角数据的近似对称、局部线性等特点,提出自适应分段拟合滤波方法。
1 自适应分段拟合滤波原理 1.1 连续变迎角数据特点 连续变迎角试验中模型存在振动,因而天平的测量结果中也叠加了该振动引起的天平输出。图 1为某次连续变迎角(3(°)/s的速度改变模型的迎角α)测力试验数据,基于阶梯变迎角试验数据处理方法(主要进行电压增益修正、基于天平校准系数矩阵等计算气动力系数、并作弹性角修正等)的处理结果y(真实信号,包含时域和频域)。
图 1 连续变迎角试验数据及其功率谱 Fig. 1 Data of continuous sweeping angle of attack test and its power spectrum |
图选项 |
信号在1~4 Hz区域中存在较大的波动,表明阶梯变迎角的数据处理方法直接对连续变迎角试验的数据处理存在不足。
振动数据具有一定的对称性,在一个或者数个周期的长度中,数据的偏差之和基本为0。图 2所示的阶跃响应为典型的振动信号,从其误差(反向)积分可以看出,一个振动周期内,其偏差积分基本为0。可借鉴振动信号近似对称整周期偏差积分为0的特点,使用测力试验数据自身(不采用加速度等数据),取单个或者对整数个周期长度的数据(实际中可能有多种振动频率、多次谐波),做线性或者2次拟合,以消除振动干扰。
图 2 振动数据误差特性 Fig. 2 Error characteristics of vibration data |
图选项 |
此外,在不考虑迎角-气动载荷剧烈变化(或存在奇异点)的情况下,测力试验的结果具有连续性以及局部的线性或近似线性。气动载荷剧烈变化的局部迎角范围,可以通过试验数据观测判定,并通过阶梯变迎角试验等其他方法获得相应的试验数据。
1.2 滤波原理 针对阶梯变迎角数据处理方法处理连续变迎角的数据结果中存在大约1~4 Hz的振动信号,基于振动信号的周期性、近似对称性、连续性和局部线性等特点,提出将数据分段拟合后再重组的方法。
本文方法的具体实现步骤如图 3所示,主要包括5个步骤:①基于方差确定分段长度;②二项式拟合;③确定下一次拟合起点;④加高斯窗对拟合数据重组;⑤平滑处理。
图 3 方法实现步骤 Fig. 3 Implementation steps of proposed method |
图选项 |
1) 分段长度确定。试验获得的待处理信号,在不同的迎角下因受力等状态不同,因而具有不同的振动周期,故用基于频谱分析的方法确定周期进而确定分段长度的方法不太合适,从图 1可看出上述特性。
由于信号中有某频率的周期振荡信号,如果信号的长度刚好为整数周期倍,一次或者二次拟合的数据正好两边偏差一致,是理想的结果,如果不是整数周期倍,则将往一边偏移,带来较大的偏差。
图 4(a)中的“var”为某试验数据“y”(与图 1中的y相同)取第1~n个数据点进行二项式拟合后的方差,按照MATLAB计算得到,如式(1)所示。具体地,取长度为n的数据段y,基于其对应横坐标α做二次拟合,拟合结果与y的差为残差,各α对应点残差的方差,即为var(n)。可以看出,当分段长度逐渐增加时,方差逐渐增大,在某时刻会有一定的下降然后再上升,即会在某个长度时刻形成极小值点(图 4(b)中横坐标约-3.298处),该长度对应的数据就是大约取到整数周期,即可用该极值点作为该次分段的终点。
(1) |
图 4 某试验数据及拟合方差 Fig. 4 Data of a test and its fitting variance |
图选项 |
为了减少计算量,可以设定最小的数据长度。为了避免特殊情况引起的数据不“收敛”,可以设定数据长度上限。
2) 二项式拟合。前期数据处理结果表明,进行多项式拟合时,取二次和三次、四次的差别不大,直线拟合偏差较大。为减少计算量,可选择二次多项式拟合。
3) 确定下一分段起点。如果将数据不重叠的部分分割成多段,数据重组时将出现台阶,因此将数据滑动的部分分解为多段。为了控制重叠量,下一分段的起点,选择为上一数据段数据的20%处,即如果上一数据段的起点为K,长度为L,则下一数据段的起点为K+0.2L。该系数0.2可根据重叠数量酌情调整,系数越小,重叠次数越多。
分段结果如图 5所示,分段过程中的方差变化情况如图 6所示,每一段的方差都是逐渐增大然后下降到某个(极小)值。图 7为数据分段拟合的结果。
图 5 数据分段 Fig. 5 Data segmentation |
图选项 |
图 6 长度搜索过程中的方差变化 Fig. 6 Variance change during length search |
图选项 |
图 7 分段拟合结果 Fig. 7 Results of piecewise fitting |
图选项 |
4) 数据重组。各数据段拟合完成后,通过高斯窗函数加权取平均,对整个数据段进行融合。由于多项式拟合后两端的数据偏差较大,高斯窗将两端的权重降低,使融合后的数据尽量接近真值。
5) 平滑滤波。数据重组后,分段处可能存在细小的阶梯,通过平滑滤波予以削减。通过选择合适的高斯窗,可以出现较小的台阶。
图 8中“原始值”为使用阶梯变迎角数据处理方法得到的结果,“自适应多项式平滑”为本文方法的处理结果,“残差”为上述两信号对应的差值。
图 8 滤波结果及对应残差 Fig. 8 Filtered results and corresponding residuals |
图选项 |
从滤波整体效果看,该方法滤波后的曲线明显比原始曲线平滑,且基本在原始曲线中心,具体表现为滤波残差在0附近上下波动。
数据段末尾存在较大的偏差,是因为数据段取到最后时,该数据段的长度并没有满足方差极小值点的条件,因而拟合出的曲线存在相对较大的“偏移”。但是该特点基本不影响工程应用,解决方案1是多取一段迎角数据,基于该方法处理后丢掉头尾,方案2是利用数据反转的方法,分别从前往后和从后往前取数据按照本文方法处理后再融合,确保每一段都取到方差极小。
2 测试信号误差分析 生成一个抛物线信号作为测试(原始)信号,以3个谐波和白噪声叠加模拟“噪声”信号,对本文方法进行测试。测试信号y0以及带噪声测试信号y1按照式(2)和式(3)计算得到。图 9给出了测试信号和图 1所示真实信号的时域和频域特征,测试信号的波形与功率谱与真实信号基本相当,表明测试信号具有一定代表性。
(2) |
(3) |
(4) |
图 9 测试信号 Fig. 9 Test signal |
图选项 |
测试结果如图 10所示,包含了测试信号、带噪声测试信号、滤波效果及滤波残差(残差等于滤波结果减去测试信号)。图 11为各信号的功率谱。
图 10 测试信号滤波结果及对应残差 Fig. 10 Filtered results and corresponding residuals of test signal |
图选项 |
图 11 测试信号及滤波后的功率谱 Fig. 11 Power spectrum of test signals and filtered results |
图选项 |
从功率谱看,1~4 Hz的波动信号大部分被滤除,残差为上下波动的原始信号。
表 1给出了滤波后的信号y2以及带噪声的测试信号y1,按照式(5)~式(10)计算的平均残差能量有效值powerk、平均残差能量有效值比例power_pk、最大残差max_dertak、最大(残差)比例max_derta_pk、残差均值dertak以及残差与最大值的比例derta_pk,其中yk分别为y1或y2,y0为测试信号,计算时丢掉了首尾各34个点。
表 1 试验数据分析结果 Table 1 Analyzed results of test data
y | powerk/% | power_pk/% | max_dertak/% | max_derta_pk/% | dertak/10-6 | derta_pk/10-6 |
y1 | 0.104 4 | 0.115 2 | 0.222 3 | 0.245 9 | -10.07 | -11.11 |
y2 | 0.035 7 | 0.039 4 | 0.083 3 | 0.092 0 | 7.91 | 8.73 |
34.24 | 34.24 | 37.40 | 37.40 | -78.52 | -78.52 |
表选项
(5) |
(6) |
(7) |
(8) |
(9) |
(10) |
从结果看,使用本文方法可以取得较好效果,具体表现如下:
1) 将测试信号的平均残差能量有效值比例从0.115 2%降低到0.039 4%,降低到了原来的约34%,振动信号去除明显。
2) 最大残差从原来的0.222 3%降至0.083 3%,能够有效提高常规测力试验精度。
3) 残差在均值附近波动,残差均值在滤波前后残差变化不大,比例在10×10-6左右,远小于试验精度需求,说明该滤波方法引入的“稳态偏差”很小,基本可以忽略。
3 标模验证与试验应用 图 12(a)为某高超通气标模在总压0.98 MPa,马赫数Ma=4.992条件下,测得的使用本文滤波方法滤波前(原始)和滤波后(滤波)的Cmz值和CA值及其残差(滤波前后数据之差),图 12(b)为图 12(a)的局部放大图(隐藏了残差数据),图 12(c)为图 12(a)的功率谱(其中CA功率谱放大了10倍)。从图 12(a)可看出残差在0上下波动;图 12(b)可看出滤波结果在原始曲线中央,图 12(b)显示滤波结果明显较平坦,说明稳态偏差较小;图 12(c)显示1~4 Hz的谱线基本被滤除,达到预期目的。
图 12 某高超标模测试结果 Fig. 12 Test results of a hypervelocity normal mode |
图选项 |
图 13为某型号试验的数据,其中包含连续变迎角实验数据基于阶梯变迎角方法处理后的“原始数据”U1,本文方法处理的“滤波结果”U2和“阶梯变迎角”试验数据U0。表 2给出了U0、U1、U2滤波前后的数据,以及相对于阶梯变迎角数据的偏差比例e1=(U1-U0)/U0和e2=(U2-U0)/U0。
图 13 某试验模型测试结果 Fig. 13 Test results based on some test model |
图选项 |
表 2 滤波效果分析 Table 2 Analyzed filtering results
α/(°) | 阶梯数据U0/V | 滤波前数据U1/V | 滤波后数据U2/V | 滤波前偏差e1/% | 滤波后偏差e2/% |
-6.45 | 0.959 8 | 0.955 8 | 0.955 2 | -0.41 | -0.47 |
-4.30 | 0.911 6 | 0.912 4 | 0.908 7 | 0.09 | -0.32 |
-3.24 | 0.878 7 | 0.880 4 | 0.881 1 | 0.19 | 0.28 |
-2.20 | 0.832 8 | 0.829 9 | 0.833 8 | -0.35 | 0.12 |
-1.79 | 0.803 1 | 0.807 2 | 0.805 8 | 0.52 | 0.33 |
-1.16 | 0.762 7 | 0.767 1 | 0.769 2 | 0.58 | 0.85 |
-0.13 | 0.734 0 | 0.735 9 | 0.735 3 | 0.27 | 0.19 |
0.28 | 0.758 8 | 0.756 9 | 0.761 3 | -0.25 | 0.33 |
0.91 | 0.805 4 | 0.809 5 | 0.804 7 | 0.52 | -0.08 |
1.97 | 0.805 2 | 0.812 6 | 0.812 1 | 0.93 | 0.86 |
3.02 | 0.809 4 | 0.809 3 | 0.812 8 | -0.01 | 0.41 |
4.08 | 0.851 6 | 0.851 9 | 0.851 6 | 0.04 | 0 |
6.20 | 0.904 2 | 0.902 9 | 0.902 8 | -0.13 | -0.14 |
表选项
结果显示,本文方法可以滤除一定的波动,总体上能够减小振动带来的偏差,且偏差满足常规测力试验的精度要求。
局部地,在变化较剧烈的局部范围内(如本实例中的1°~3°范围,表 2中迎角为0.91°、1.97°和3.02°的状态),偏差相对较大,说明该范围内可能存在比较“特殊”的气动力特性,需针对该范围开展进一步试验。采用阶梯变迎角,或者更慢更稳的连续变迎角试验,是较好的解决方案。
4 结束语 本文提出基于方差极小自适应确定分段长度的分段拟合滤波方法,可以较好地滤除连续变迎角试验数据中1~4 Hz的低频噪声,有效提高了连续变迎角试验数据处理精度,可推进连续变迎角试验的应用,提高风洞试验效率。下一步将开展小波等更多分析方法的对比,探索更好的试验数据分析处理方法。
参考文献
[1] | 舒海峰, 许晓斌, 孙鹏. 高超声速风洞多天平测力试验技术研究[J]. 实验流体力学, 2014, 28(4): 49-53. SHU H F, XU X B, SUN P. Technology investigation on force test with multi-balance in hypersonic wind tunnel[J]. Journal of Experiments in Fluid Mechanics, 2014, 28(4): 49-53. (in Chinese) |
[2] | 谢艳, 李平, 王瑞波, 等. 2.4米跨声速风洞连续变攻角测力试验技术研究[J]. 气动研究与发展, 2009, 19(3): 10-14. XIE Y, LI P, WANG R B, et al. An experimental investigation on force measurement of continuous angle of attack traverses in 2.4 m transonic wind tunnel[J]. Aerodynamics Research and Development, 2009, 19(3): 10-14. (in Chinese) |
[3] | 魏志, 谢艳, 吴军强, 等. 连续变攻角测力试验技术在大型暂冲式跨声速风洞中的应用[J]. 实验流体力学, 2011, 25(4): 99-102. WEI Z, XIE Y, WU J Q, et al. JApplication of continuous sweeping force measuring technology in large intermittent transonic wind tunnel[J]. Journal of Experiments in Fluid Mechanics, 2011, 25(4): 99-102. DOI:10.3969/j.issn.1672-9897.2011.04.019 (in Chinese) |
[4] | 黄辉, 黄昊宇, 凌忠伟, 等. Φ0.5米高超声速风洞连续变攻角测力试验数据处理方法研究[J]. 计算机测量与控制, 2019, 27(8): 281-285. HUANG H, HUANG H Y, LING Z W, et al. Research on data processing method of continuous variable angle of attack force in Φ0.5 m hypersonic wind tunne[J]. Computer Measurement & Control, 2019, 27(8): 281-285. (in Chinese) |
[5] | 唐志共, 许晓斌, 杨彦广, 等. 高超声速风洞气动力试验技术进展[J]. 航空学报, 2015, 36(1): 86-97. TANG Z G, XU X B, YANG Y G, et al. Research progress on hypersonic wind tunnel aerodynamic testing techniques[J]. Acta Aeronautica et Astronautica Sinica, 2015, 36(1): 86-97. (in Chinese) |
[6] | 张俊. 工程实用的飞行器低速风洞连续扫描试验技术研究[D]. 长沙: 国防科学技术大学, 2007. ZHANG J. Investigation of the continuous scan technology on engineer to aircraft in low speed wind tunnel[D]. Changsha: National University of Defense Technology, 2007(in Chinese). |
[7] | 唐乔乔, 张卫国, 刘忠华, 等. 8 m×6 m风洞特大迎角机构连续扫描试验技术研究与应用[J]. 实验流体力学, 2012, 26(2): 81-85. TANG Q Q, ZHANG W G, LIU Z H, et al. Research and application of the continuous scan technique to the high angle of attack equipment in 8 m×6 m wind tunnel[J]. Journal of Experiments in Fluid Mechanics, 2012, 26(2): 81-85. DOI:10.3969/j.issn.1672-9897.2012.02.018 (in Chinese) |
[8] | 张双喜, 盖文, 褚卫华, 等. 风洞试验连续变攻角控制策略[J]. 计算机测量与控制, 2014, 22(2): 390-396. ZHANG S X, GAI W, CHU W H, et al. Control strategy of uniform changing attack angle in wind tunnel[J]. Computer Measurement & Control, 2014, 22(2): 390-396. DOI:10.3969/j.issn.1671-4598.2014.02.025 (in Chinese) |
[9] | 孟宝清, 韩桂来, 姜宗林. 结构振动对大型激波风洞气动力测量的干扰[J]. 力学学报, 2016, 48(1): 102-110. MENG B Q, HAN G L, JIANG Z L. Theoretical investigation on aerodynamic force measurement interfered by structural vibrations in large shock tunnel[J]. Chinese Journal of Theoretical and Applied Mechanics, 2016, 48(1): 102-110. (in Chinese) |
[10] | 程忠宇, 陈宏, 张琦. 多加速度计振动分离惯性补偿测力技术[J]. 流体力学实验与测量, 1993, 13(4): 57-61. CHENG Z Y, CHEN H, ZHANG Q. Inertia compensation technology based on multi-accelerometer vibration separating[J]. Experiments and Measurements in Fluid Mechanics, 1993, 13(4): 57-61. (in Chinese) |
[11] | 吕金州, 张小庆, 陈光雄, 等. 基于惯性补偿的脉冲风洞测力天平瞬态研究[J]. 振动与冲击, 2018, 37(2): 216-222. LV J Z, ZHANG X Q, CHEN G X, et al. Transient simulation for dynamic output of force measuring balance in animpulse combustion wind tunnel based on intertia compensation[J]. Journal of Vibration and Shock, 2018, 37(2): 216-222. (in Chinese) |
[12] | 艾迪, 许晓斌, 王雄. 风洞天平动态特性多阶惯性补偿技术研究[J]. 实验流体力学, 2018, 32(4): 87-92. AI D, XU X B, WANG X. Investigation of wind tunnel balance dynamic characteristics' multi-order inertial compensation[J]. Journal of Experiments in Fluid Mechanics, 2018, 32(4): 87-92. (in Chinese) |
[13] | 张鹏, 谢艳, 孙宁, 等. 基于改进小波阈值函数的风洞连续信号降噪方法[J]. 计算机测量与控制, 2014, 22(4): 1300-1302. ZHANG P, XIE Y, SUN N, et al. Wind tunnel continuous signal de-noising method based on improved thresholding function[J]. Computer Measurement & Control, 2014, 22(4): 1300-1302. DOI:10.3969/j.issn.1671-4598.2014.04.103 (in Chinese) |
[14] | 张鹏, 刘晨雨, 曹宇晴. 基于Matlab GUI的风洞信号小波分析处理软件[J]. 兵工自动化, 2018, 37(1): 61-67. ZHANG P, LIU C Y, CAO Y Q. Analysis and processing software of wind tunnel signals on the basis of the Matlab GUI programming[J]. Ordnance Industry Automation, 2018, 37(1): 61-67. (in Chinese) |
[15] | 张鹏, 魏志, 王春, 等. 基于小波变换的风洞连续信号降噪分析[J]. 兵工自动化, 2013, 32(5): 63-67. ZHANG P, WEI Z, WANG C, et al. Analysis of wind tunnel continuous signal de-noising based on wavelet transformation[J]. Ordnance Industry Automation, 2013, 32(5): 63-67. (in Chinese) |