分段最小二乘法(共8篇)
分段最小二乘法 篇1
在工程实践与科学实验中, 常常需要从一组带噪声的试验观测数据 ( xi, yi) ; i = 1, 2, …, n中找出自变量x与因变量y之间隐含的函数关系, 数据拟合[1]是一种常用的处理方法。其中多项式曲线拟合又是一种较常用的数据拟合方法。当数据点较多时, 多项式阶数太低, 拟合精度和效果不太理想。要提高拟合精度和效果就需要提高曲线阶数, 但阶数太高又带来计算上的复杂性及其他方面的不利。因此, 如果只采用一种多项式曲线函数拟合较多的数据点, 难以取得较好的拟合精度和效果。为有效的解决上述问题, 一般采用分段曲线拟合, 在每段区间上进行局部最小二乘拟合[2—4]。
传统的分段曲线拟合根据主观经验和绘制数据散点图来确定拟合的经验函数和分段点。文献[5] 提出分段区间重合的拟合方法, 由每4个数据点决定一个三次曲线, 但分段区间太密, 不适用于密集的数据拟合。文献[6]提出的多项式基函数的全局连续拟合方法, 只限于2个分段区间。文献[7]提出多分段区间全局连续的曲线拟合方法, 但基函数只限于一次多项式。文献[5—7]提出的方法都是根据主观经验来确定经验函数和分段区间, 然后进行拟合。文献[8]引入拟合误差限度 εmax, 在误差大于该限度的点处分段, 但该限度 εmax不容易界定。
提出一种自动分段多项式曲线拟合方法, 该方法在实际运用中具有如下有点: 1提供几种不同的经验函数, 根据不同经验函数拟合的数据和实测数据的误差的方差, 自动确定较优的经验函数; 2自动确定数据的分段区间, 在各个区间进行最小二乘拟合。通过数值实验, 证明该方法有较高的拟合精度。
1经验函数的确定及拟合步骤
1. 1经验函数的确定
数据拟合方法有很多, 例如对数曲线拟合, 反函数曲线拟合, 二次曲线拟合, 三次曲线拟合, 幂函数曲线拟合, 指数曲线拟合等。一般先观察散点图来确定曲线的类型, 不过散点图都是相关关系的粗略表示, 有时候散点图可能与几种曲线都很接近, 这时建立相应的经验函数可能都是合理的, 但由于选择不同的曲线, 得到同一个问题的多个不同经验函数, 怎样从这些经验函数中选择最优的一个。文献[1] 提出的, 用几种函数进行拟合, 计算历史数据点实测值和拟合值的误差平方和最小的为最优经验函数这一方法, 可能存在这样的问题: 误差平方和最小, 但误差波动较大, 即一些点误差很小, 一些点误差相对较大。针对这种情况, 本文提出一种新的确定经验函数的方法, 用几种不同的函数进行拟合, 从中选取最优的经验函数, 最优经验函数确定的条件如下:
1) 历史数据点误差为正和误差为负的个数之差小于适应性参数, 在本文试验中选用的适应性参数为3。
2) 计算误差的方差, 方差最小的为最优的拟合函数。
方差是各个数据与平均数之差的平方和的平均数。通俗点讲, 就是各点与平均数偏离的程度, 来衡量一批数据相较于平均数的波动的大小。方差的数学描述为
式 ( 1) 中把x-视为第一个历史点的误差, xi为后面各历史点的误差, 该公式即可理解为各点误差相较于第一个点的误差的波动大小, 波动越小, 说明各个点的误差相差越小, 分布较均匀, 那么该函数即为最优经验函数。经验函数确定的流程图如图1所示。
1. 2自动分段拟合的步骤
本文选取三种较为典型的回归方程类型
分段拟合步骤如下, 分段拟合流程图如图2所示。
1) 根据1. 1节提出的方法拟合所有的历史数据, 在上述的三种拟合函数中选择最优的拟合函数。
2) 由第1) 步选取的经验函数, 计算历史数据点拟合值和实测值的误差, 计算历史数据点误差绝对值的均值S。
3) 比较各历史点误差与误差均值S, 若连续3个点的误差绝对值大于均值S, 则从第一个误差大于均值的点处分段, 执行第4) 步, 若不存在拟合值和实测值误差大于均值的点, 或者这种点只有一个或两个, 则不进行分段, 执行第5) 步。
4) 从分段点处到历史数据最后一个点重新拟合, 选择最优的函数, 重复以上步骤。
5) 根据以上步骤得到的分段, 以及各段选出的最优经验函数对历史数据进行拟合。
2数值实验
2. 1自动分段拟合实验
为检验所论述的自动分段多项式曲线拟合方法, 采用两组数据, 厦门高崎机场1990 ~ 2012年的起降架次的数据[9]和1990 ~ 2004年中国人均GDP的数据来进行拟合。
提出的自动分段多项式曲线拟合方法对厦门高崎机场1990 ~ 2012年起降架次统计数据和1990 ~ 2004年中国人均GDP数据进行拟合, 用JAVA编程实现, JFree Chart绘制曲线。表1和图3是厦门高崎机场起降架次的数据与经过分段曲线拟合后的数据的比较, 表2和图4是中国人均GDP数据与拟合后的数据的比较。
对于厦门高崎机场起降架次的数据, 提出的该方法自动将数据分为3段, 即数据1990 ~ 1994, 1994 ~ 2003, 2003 ~ 2012, 这三段数据选取的最优经验函数都是指数函数形式。对于中国人均GDP的数据, 程序自动将数据分为三段, 即数据段1990 ~ 1995, 1995 ~ 1999, 1999 ~ 2004, 这三段数据选取的最优经验函数分别是多项式函数形式, 直线形式和多项式函数形式。由表1、表2和图3、图4可以看出, 提出的自动分段拟合方法拟合精度较高, 能够很好的拟合历史数据。
2. 2对比实验
本文选取1990 ~ 2004年中国人均GDP数据, 分别用本文提出的自动分段拟合方法和传统直线拟合、指数拟合、多项式拟合方法对该组数据进行实验, 对比实验各种方法的相对误差如表3所示。
由表3可以看出, 自动分段拟合方法相对误差比传统拟合方法相对误差小, 且相对误差波动比传统方法波动小。自动分段拟合方法不仅在拟合精度上比传统方法高, 而且分段方便, 不需用户根据主观经验和绘制散点图来分段, 方便用户使用。
3结论
提出的自动分段曲线拟合方法采取对历史数据进行拟合, 自动选取较优的经验函数, 自动进行分段, 使得在进行数据拟合时, 不需要人为的绘制出历史数据散点图来选取经验函数和分段。通过编程实现和数值实验验证, 该方法不仅便于在计算机上编程实现, 而且拟合效果较好。
摘要:针对传统的分段曲线拟合方法在选择拟合函数和确定分段区间时经验成分较多的不足, 提出一种自动分段多项式曲线拟合方法, 根据误差方差和误差均值, 自动确定经验函数和分段区间。通过实际数据的检验, 验证了该方法的拟合效果。
关键词:数据拟合,分段拟合,多项式曲线,最小二乘法
参考文献
[1] 吴宗敏.散乱数据拟合的模型、方法和理论.北京:科学出版社, 2007
[2] 杨维, 张晓明.利用最小二乘法进行回归分析及经验公式的确定.沈阳工业大学学报, 1991;13 (2) :1—6
[3] 侯超钧, 曾艳姗, 吴东庆, 等.全局连续的分段最小二乘曲线拟合方法.重庆师范大学学报 (自然科学版) , 2011;28 (6) :44—48
[4] Grama L, Rusu C.Phase approximation by divide-and-conquer piecewise linear fitting of gain.International Symposium on Signals, Circuits and Systems, 2009:1—4
[5] 蔡山, 张浩, 陈洪辉, 等.基于最小二乘法的分段三次曲线拟合方法研究.科学技术与工程, 2007;7 (3) :352—355
[6] 刘晓莉, 陈春梅.基于最小二乘原理的分段曲线拟合方法.伊犁教育学院学报, 2004;17 (3) :132—134
[7] Gluss B.An alternative method for continuous line segment curve-fitting.Information and Control, 1964;7 (2) :200—206
[8] 王成山, 黄碧斌, 李鹏, 等.燃料电池复杂非线性静态特性模型简化方法.电力系统自动化, 2011;35 (7) :64—69
[9] 刘强, 朱敏, 王小维, 等.灰色神经模型在空中交通流量预测中的应用.全国第19届计算机技术与应用学术会议, 乐山, 2008
分段最小二乘法 篇2
【关键词】BP神经网络 最小二乘法 盾构故障
一、故障预测常用方法种类及特点
常用的故障预测种类有以下几种:
(一)曲线拟合法:本方法简单实用,好理解,较多的应用在监控机械装置的系统,缺点是预测误差结果比较大。
(二)卡尔曼滤波法:计算量小、预测精度高,主要用于线性系统,在非线性系统中需要进行扩展,但其模型的不确定性很差。
(三)灰色预测:该方法需将两个部分进行置换,即“随机过程”转换为“灰色过程”,把“随机量”处理成“灰色量”,然后在理论系统的模型GM(1,1)中进行处理。
(四)神经网络:具有学习记忆功能,能很好解决非线性系统问题,不对预测模型做限制,同时还能把以前的历史数据完整映射到未来的数据库中,可广泛用在故障检测中。实际使用广泛的神经网络预测有:BP神经网络和自组织特征映射网络。
二、最小二乘法与BP神经网络相结合的预测方法
基于本文采用最小二乘法与BP神经网络相结合的方法来预测故障,因此需要对所选取的数据样本集进行直线最小二乘法的拟合,然后用BP神经网络算法进行计算,最后将这些值作为样本集,输入到BP神经网络模型来进行训练预测。
(一)BP神经网络模型的构建
在构建神经网络之前,需将数据进行归一化处理,将归一化的数据样本构建成BP神经网络模型。在本论文的神经网络模型中,训练函数,学习函数,传读函数,性能函数,训练次数,训练误差,学习率分别为:traingdx,learngdm,logsig,mse,1000次,0.0001,0.08.输入层,输出层及隐层分别为:6,1,5。
(二)最小二乘法的直线拟合过程
根据最小二乘法的理论[4]及BP神经网络的模型可知,要预测故障,既是对出现故障的时间序列的预测,而预测时间序列X={Xt︱t=0,1,2…m},实际上就是利用已知的{xk,xk-1…xk-n+1}来求xk+1。
将得到的数据样本再通过BP神经网络的算法对样本计算,得到其对应的各个参变量值,即给出Xk,得到期望值N(Xk),其中k=n-1,n…m,然后把这些值作为样本集,用在BP神经网络中来训练整个网络。
三、盾构机PU电流序列预测应用
(一)单一的BP神经网络的预测结果
单一的BP神经网络采用以1小时为间隔所提取得 100个PU电流数据作为样本集对进行训练,采用输入连续的6个PU电流值,紧接着输出相连续的一个PU电流值,即输入与输出层节点数比例为6:1,隐层节点数为5,将网络进行收敛,当达到预定精度以后,预测其结果如图1中曲线。
(二)最小二乘法及BP神经网络相结合的预测结果
用上面的样本集数据,对最小二乘法与BP神经网络的合成的样本进行训练,并进行网络收敛,当达到预定精度,预测其结果如图2中曲线所示。
(三)两者结果比较
图1和图2中的两条预测线分别代表PU实际电流变化曲线(实线)以及用神经网络预测的电流变化曲线(虚线)。可以看到,在图2中,这两种曲线的变化趋势比较接近,仅有极少的时间序列点有较大偏差,而图1中的两条曲线则在多个时间序列段出现较大的偏差,两条曲线走势明显偏离,所以,可以认为,采用合成的BP神经网络算法预测到的PU电流时间序列精度较高,比采用单一的BP神经网络更能客观反映盾构机的实际电流变化情况。虽然有个边点出现误差,也是未考虑机组负荷的影响。因此,整体上说合成的BP神经网络算法对盾构机的PU电流还是具有很高的预测能力。
图1 图2
四、总结
上述内容主要采用最小二乘法与BP神经网络结合的方法对盾构机的故障进行预测,将二者结合起来,可以充分利用二者的优点,优势互补。笔者将其应用到了实际的工作中,起到了良好的预测效果。
参考文献:
[1]蒋瑜,杨雪,阮启明.机械设备故障规律及运行趋势预测方法综述[J].机电一体化,2001(3):14-17.
[2]陈敏泽,周东华.动态系统的故障预报技术[J].控制理论与应用,2003,20(6):820-822.
[3]苏春华,罗雷,海军,梅检民,肖云魁.基于BP神经网络的汽车发动机寿命预测[J].军事交通学院学报,2009,11(4):50-51.
基于最小二乘法离场航迹构造方法 篇3
飞机离场航迹是飞机起飞过程的形象化体现。目前,飞行航迹的构造方法通常采用两类方式,第一类是在飞行计划确定及气象资料完整的情况下,结合飞行动力学和运动学模型正向推导的方法。飞机起飞航迹计算方法研究提出了对不同机型统一的离场航迹构造方法,该方法主要针对离场航迹剖面进行了构造,缺少对水平投影面构造的方法;基于ANP数据库的飞机起飞仿真研究是基于详细的飞行资料和性能参数的前提下,提出了飞机离场剖面航迹构造方法;离场航迹降噪优化设计的多目标智能方法是一种利用航段飞行特征约束求解离场航迹的方法。第二类是在拥有较为准确的雷达位置信息点的情况下对雷达数据去噪,拟合出最佳函数匹配从而得到平滑的航迹。经纬仪目标交汇测量及航迹曲线拟合文中提出根据不同时刻的坐标,用最小二乘法对目标航迹进行拟合,从而推测下一时刻的位置速度及加速度;三维航迹的B样条曲线拟合算法利用B样条曲线的几何性质,解决了飞行器三维航迹拟合中的边界条件等约束问题。第二类多用于飞机离场结束后航迹的拟合。上述提出的两类方法用于离场航迹的构造存在以下三种问题:一是由于飞机离场属于低空飞行,雷达捕捉飞行器在低空飞行的位置信息不准确,飞机离场的雷达点相比于真实点误差较大,且飞机离场的方式不同,导致无法单一的利用拟合离场雷达数据的方法确定离场航迹。二是现有方法多为对离场剖面航迹进行构造,忽略由于离场方式的不同导致水平面航迹存在误差。三是在离场数据、飞行计划及气象资料缺失的情况下,无法对飞机离场航迹进行构造的问题。由此本文为了解决上述问题,采用最小二乘法结合两种离场方式特征,提出离场航迹逆向构造方法。
相关工作
飞机离场过程是指飞机高于起飞表面450m(1500ft),并完成从起飞到航路形态的转变,达到规定的速度和爬升梯度。飞机起飞过程包括起飞场道滑跑阶段和起飞航道阶段。起飞航迹依据飞机的构型、发动机的推力状态、对爬升梯度的要求等分为第一爬升阶段、第二爬升阶段、第三平飞阶段、最后爬升阶段等四个阶段。本文忽略平飞过程,把此过程与第四阶段融合在一起下文统称第三阶段,分别对第一、二及三阶段分别进行构造。
最小二乘法是通过最小化误差的平方和求得待定系数从而寻找数据的最佳函数匹配。假设多项式:
由最小二乘法确定待定系数a0,a1,(43)am,设数据点权为1,令
得方程组:
该方程组称为多项式拟合的法方程,令:
基于最小二乘法飞机离场航迹逆向构造方法
飞行高度在350m以上的雷达监测位置信息点较为准确,本文把航迹投影到水平面和剖面分别构造,并利用准确数据逆推离场航迹。首先根据两种离场方式特征判别离场方式,其次利用350m以上的雷达监测位置信息点分别拟合两种离场方式的第三阶段航迹水平面投影,最后计算各阶段关键参数确定滑跑航迹及运动过程。不得不提在离场数据、飞行计划、气象资料缺失前提下。起飞过程模拟难度很大,为了确保模拟的航迹准确性,不可避免需要根据飞机的离场规则,假设合理数值。
P0(x0,y0,z0)点为已知飞机在跑道滑跑的加速始点,此点的各方向速度均为V0=0,P1(x1,y1,z1)点为飞机的离地点,其中z1已知为地平面高度,P2(x2,y2,z2)为直线离场时高度为时的坐标点,转弯离场时为转弯点,P3(x3,y3,z3)为已知雷达数据第一个点。
基于贝叶斯推理的飞机离场方式估计方法
直线离场方式特征:特征一,起始离场航线与跑道中线方向角度相差小于15°;特征二,离场航迹偏于跑道中线一侧而在DER(离场末端)的横向距离不大于300m。但只要实际可能,离场航线就应与跑道中线延长线一致。转弯离场方式特征:离场过程中出现离场角度要求大于15°的转弯,并且规定在飞机起飞离场到达DER标高之上才允许转弯,在此之前为直线飞行。
由于已知雷达数据的高度值不同分为以下两种情况,利用不同的离场特征分别作为判别离场方式的依据。
(1)现有雷达数据高度值在450m以上,飞机可能已经结束离场并立即发生转弯,由雷达数据得到的航线与跑道中线的夹角大于15°不能说明转弯发生在离场过程中,因此不能用特征一判断,利用特征二更加适合。利用特征二区分时,计算雷达数据第一个点与跑道直线的水平距离如果小于300m则为直线,否则为转弯。
(2)现有雷达数据高度值在450m以下,飞机还未结束离场,利用特征一就可以判断。已知离场结束后航迹的雷达数据,前几个点组成的航迹是顺沿离场航迹的方向产生,继承了离场方向变化趋势。利用特征一区分时,利用前四个数据点进行一次的最小二乘法拟合方程,跑道直线方程,由正切公式计算两直线夹角,则,当时为转弯离场,否则为直线离场。
爬升阶段构造方法
坐标由公式得:
公式为:
起飞场道阶段构造方法
在实际计算起飞滑跑距离时,可以将地面滑跑期间发动机的推力取平均值并视为常数,同样换算摩擦系数也取其平均值并看作常数,飞机地面滑跑距离的近似估算公式:
综上飞机离场航迹由此确定。
方法验证
本实验以首都机场为例,采用首都机场提供的2013年4、5月的飞机飞行数据及机场地理位置信息,随机抽取1000条离场航迹进行实验。首先,依据构造方法对数据进行处理拟合出离场航迹,利用拟合出的航迹计算出飞机离场滑跑距离及转弯高度,查看计算出的上述参数是否符合飞机离场规定,从而验证方法的有效性。利用基于贝叶斯推理的飞机离场方式选择评分函数计算分数,验证飞机离场方式选择方法的准确性。
为了大致了解离场航迹线的特点,利用雷达数据做出三维散点图。进而根据航迹线水平面的投影的特征对航迹线进行分类。航迹线水平面投影类型分为两类,第一类为近似顺延跑道延长线如图2(a)。第二类为偏离跑道延长线,在跑道的一侧与跑道延长线形成较大夹角如图2(b)。从对雷达数据初步分析可以看出本文对离场模型的分类假设是正确的。
图3(a)和(b)分别为直线离场实验前350m以上原始雷达数据三维曲线图和实验后效果图。图2(c)和(d)分别为转弯离场实验前350m以上原始雷达数据三维曲线图和实验后效果图。通过实验效果前后对比可以看出利用本文方法构造出的离场航迹与实际情况较为一致。
由表1可以看出转弯高度均允许转弯高度值120m以上。滑跑距离在均合理范围以内。
由统计学办法得出P(line)及P(arc)。利用飞机离场方式选择方法从直线离场航迹中筛选出符合直线起飞特征的航迹及符合转弯起飞特征的航迹,从转弯离场航迹中筛选符合直线起飞特征的航迹及符合转弯起飞特征的航迹,利用统计学办法得出P(w/line)、、P(w/arc)及。最后利用评分方法计算评分。评分结果均在85分以上表明离场方式选择方法的准确度较高。
由图1和表2可以看出本文对离场模型的分类假设是正确的,离场方式选择方法准确度较高。由图2和表1可以看出本文构造出的离场航迹与实际情况一致。由此验证了基于最小二乘法离场航迹构造方法的有效性及准确性。
结语
本文提出了一种基于最小二乘法离场航迹逆向构造方法。新方法的可行性已在多次实验中得到验证。并利用基于贝叶斯推理的飞机离场航迹选择评分函数计算评分,评分结果理想,验证了离场方式的选取方法的有效性。新方法不仅可以利用于解决在飞行计划和气象资料缺失的前提下,无法对飞机离场航迹进行构造的问题,还可以为场间雷达与空中雷达连接方案的制定提供有利参考。
分段最小二乘法 篇4
随着微机技术不断应用到无损检测领域,同时无损检测本身的新方法和新技术也不断出现,从而使无损检测技术在工业材料的检测及质量管理中成为一个不可缺少的重要环节。钢铁材料电磁检测技术是无损快速检测技术领域中的一种实用技术,其目的是利用钢铁件内部结构异常或缺陷存在所引起的磁导率变化来进行各种钢材、零部件的硬度、含碳量测定及钢铁件混料分选,以便改进制造工艺,提高产品质量。本文将C8051F020单片机用于钢材无损检测分选仪的设计中,一方面提高了工件检测精度和速度,另一方面缩小了仪器体积,便于检测人员携带。
1 电磁检测原理
电磁检测是从使用电流的磁性方法中分离出来的,电磁方法也称磁感应法,它是以电磁感应为基础,利用交变磁场直接作用于铁磁材料本身,并通过对感应电流或电压振幅、相位的科学分析,从而成功地完成钢铁材料的性能测试、质量检查和监控,进行混料分选,完成材料的硬度、表面淬硬层或覆盖层深度的测定。该方法采用低频电流,通常不超过400Hz,比较常用的是60Hz~100Hz。电磁检测系统通常由3部分组成,即交变电流的检测线圈(探头)、检测电流的仪器和被检的金属工件。
电磁检测钢铁质量的磁导率检测法按电源类型可分为直流法和交流法两种,但第一种方法速度慢,不易实现自动化,而交流磁导率法是目前国内外应用最广泛的一种方法。交流磁导率法又分为中强磁场下的磁导率法和弱磁场下的初始磁导率法,二者区别是前者在实用中工件磁化区域多处于巴克豪森跃迁区,故对供电电源要求很高,弱磁场下的初始磁导率法避免了上述缺陷。所谓初始磁导率法,即工件在磁畴畴壁的可逆位移区域被磁化的磁导率法[1]。钢铁件的成分和组织结构不同,必然导致其磁性能和机械性能不同,即钢铁件成分和组织结构与其磁性能和机械性能是直接相关的。本文采用初始磁导率法,用电磁分选仪测出检测线圈的感应电势或感应电流,便可求得其相对磁导率,从而根据钢铁件的成分和机械性能与磁导率之间的相关性,间接测量出其硬度、含碳量,并可分选钢铁混料。
2 钢材无损检测分选仪设计
2.1 硬件设计
该仪器的硬件整体设计见图1,以通用型C8051F020单片机为核心,采用单片机片内ADC模块,通过模拟通道读入外部简单处理过的差分反馈信号,通过数字端口读入键盘的数据,通过数字端口输出实时的数据信号以驱动LCD显示,通过数字端口输出声光报警控制信号及外部设备控制信号,通过数字端口设置模拟开关选通,用以控制励磁电流的大小及频率高低。其特点在于:应用了C8051F020型单片机,将ADC模块、信号放大模块集成,整体封装采用TQFP贴片形式,大大缩小了占用空间,提高了整机的移动性能;在仪器的前向通道设计中,采用高性能仪表运算放大器及相应的噪声抑制措施对探头反馈信号进行处理后直接送单片机,减少外围电路,同时减小了信号的失真,保证了信号的完整性;键盘采用按钮加IN4148的组合形式,给出类似BCD码的数据形式,一改传统的扫描形式,不仅提高了键盘的反应速度,而且减少了单片机资源的浪费,提高了整体速度与精度;分选仪在分析处理现实数据的同时还能根据信号及设定的程序控制外部设备。为了增强分选仪的实用性,励磁信号的强度及频率必须能够变化,单片机通过MAX335模拟开关控制励磁信号强度及频率,考虑到稳定及避免占用MCU资源,励磁信号的频率由晶振及分频器CD4060提供。
2.2 软件设计
采用C8051F020的软件系统示意图见图2[2]。C8051F020的U-EC2适配器通过USB端口和PC主机连接,它的另一端插入目标板的标准JTAG接口上,JTAG接口与PC主机之间利用单片机的TMS、TCK、TDI和TDO线,通过在线通信协议进行通信。Silabs公司的集成开发环境和ICP(在线编程)功能可对C8051F系列单片机的FLASH存储器进行编程,利用该软件工具,用户可以在Windows操作系统下实现书写程序、编译、下载、电路内仿真及调试程序等功能,可对硬件和软件进行实时测试。
3 分段最小二乘法
现有分选仪的设计多数建立在钢材的材质特性参数与其磁导率参数成线性关系这一假设的基础上,采用两点法标定仪器,从而实现对未知钢材材质特性的检测,但这只是钢铁件材质-磁导率关系曲线的一种简单近似,其精度不高。其实钢铁件的材质特性与磁导率特别是初始磁导率之间往往存在双值甚至多值关系,例如初始磁导率与回火硬度在某区域往往呈现出“N”型关系。为了提高检测精度,采用基于分段最小二乘原理的数据拟合法,即根据一组测量数据,确定变量之间的函数关系。
3.1 直线的数据拟合
设钢材工件磁导率为u,硬度为h,含碳量为c,当工件的磁导率和工件的硬度在某区域呈直线时,采用基于分段最小二乘原理的直线数据拟合方法。给定一组基点u1
H1(u)=a1+b1u u11≤u
H2(u)=a2+b2u u1N1≤u
…
Hk(u)=ak+bku uk-1Nk-1≤u≤ukNk。
在整个拟合区间上就得到了k条直线段,将k条直线段连接起来,该函数就成为整个区间上的数据拟合函数,这样就得到了钢材工件的磁导率与硬度的函数表达式,同理也可求出磁导率与含碳量之间的函数表达式。
3.2 曲线的数据拟合
该仪器在许多实际检测问题中,变量之间内在的关系并不像上面所述的那样简单,呈线性关系,而是呈复杂的曲线关系。所以仅考虑用一种回归函数拟合其数据点,难以得到最佳的拟合曲线,要根据实测数据点的分布图,确定分段数目以及相应拟合曲线类型,然后再应用最小二乘法原理求得各分段拟合函数。
假设由实验得到一组数据点(ui,ci),i=1,2,…,n,根据散点图发现其前m个点较接近直线,后n-m个点呈现二次曲线关系,故分2段来进行数据拟合。前m个点用一次多项式拟合,后n-m个点用二次多项式拟合,即有:
要保证在分段点(um,cm)连续光滑,存在两个约束条件:
aum+b=k0uundefined+k1um+k2 ,
a=2k0um+k1 。
由此可推导出:
b=k2-k0uundefined。 (3)
令误差的平方和δ为最小二乘估计量,结合式(1)、式(2)和式(3)得到:
undefined。
根据最小二乘原理,当使误差的平方和δ达到最小时,曲线拟合达到最佳,因为误差都是实数,它们的平方和是正数,要求这些正数的和尽可能地小,便保证了这些误差的绝对值尽可能地小。对k0、k1、k2求偏微商,再使偏微商等于零,可求得当δ最小时k0、k1、k2的数值,进一步得到a、b的数值,从而确定拟合曲线的函数。
当然,如果由实验得到的散点图更复杂时,可以考虑分为M段,并且选用m次的多项式来拟合,然后在段点处选择约束条件,最终经过数学计算得到所需的函数关系式。
4 实验与结论
该实验制备了20件相同材质、不同硬度且尺寸均为长10cm的Φ8mm圆柱形碳钢试样,从中选取10件有代表性的作为标准试件来对仪器进行标定,另10件作为被测对象进行检测用,以检验不同算法对工件材质特性参数测量结果的影响。在检测过程中,标准试件的正确选取关系到检测的成败,为此,在仪器设计中增加了标准试件选取功能,以有效地剔除无效试样,不同算法下的测量结果见表1。
HRC
检测结果表明,采用分段最小二乘法检测的误差要比采用单一最小二乘法检测的误差小,因而将分段最小二乘拟合曲线的方法引入到检测过程中,确实能提高钢材无损检测分选仪的分选精度和效率。
参考文献
[1]张迎新.C8051F系列SOC单片机原理及应用[M].北京:国防工业出版社,2005.
[2]潘琢金,施国君.C8051FXXX高速SOC单片机原理及应用[M].北京:北京航空航天大学出版社,2002.
[3]李健生.基于模式识别的电磁无损检测系统的研究[D].哈尔滨:哈尔滨理工大学,2005:5-13.
[4]张东林.分段最小二乘曲线拟合[J].沈阳大学学报,1994(2):80-82.
分段最小二乘法 篇5
目前,我国10~35 k V系统以中性点不接地或经消弧线圈接地(称为小接地电流系统)的方式为主。中性点不接地方式的显著优点是,系统在发生单相接地故障(约占60%以上)时,可以短时继续运行。但是,此时由于电容电流的存在,就有可能引发间歇性弧光接地、铁磁谐振、中性点位移过电压等,甚至发展为事故。随着城市建设及城市电网改造,城市变电站大量采用了电力电缆送电线路,而在相对落后的农村地区,虽然农网改造已经大大地改善了当地的供电条件,但仍存在配电系统供电半径过大等问题。这些情况都会造成变电站不接地系统电容电流过大,弧光过电压时有出现,严重影响了系统和人身安全以及对用户的可靠供电。
电网电容电流的大小是决定是否装设消弧线圈以及消弧线圈调谐的依据,为了掌握电容电流增长变化情况以及正确设定消弧线圈补偿的档位,准确测量系统电容电流是非常必要的。
电容电流的测量可采用直接法和间接法。直接法由于对试验人员和配网系统均存在很多不安全的因素,现在几乎不再采用了。目前广泛采用的是分相外加电容法和变压器中性点外加电容法等间接法。但是这些测量方法都要接触到高压一次侧,且存在操作繁琐、准备工作耗时长、测量工作效率低等缺点[1]。
基于上述原因,本文介绍一种采用参数辨识原理测量电容电流的新方法。新方法采用注入信号法与递推式最小二乘法辨识相结合的方案,从PT二次侧进行采样,通过输入和输出间的数学模型,辨识计算出电网电容电流。
1 辨识算法的测量原理
1.1 注入信号法
首先,需要明确电容电流辨识的关键——信号注入法与最小二乘法辨识之间的关系。最小二乘法是一种系统辨识方法,用在这里是为了辨识电力系统的接地电容。然而,只有在一个系统模型中包括进去电力系统的接地电容,才能运用最小二乘法进行包括接地电容在内的系统参数的辨识。这里面的电力系统模型就是在信号注入法的实现方案中得到体现的。信号注入法运用到电压互感器PT,就包括了对电压互感器PT回路的电路等效与建模。运用这个电路模型作为目标系统模型,用注入信号角频率ω作为系统输入,用测量到的回路阻抗Z和相角θ作为系统输出,那么最小二乘法就能够对此系统发挥作用,有效地辨识到电压互感器PT回路的电路参数,得到接地电容,进而求得电容电流[2]。具体方案设计如图1所示。
图中WA、WB、WC分别为电压互感器(PT)三相的高压绕组,二次绕组Wa、Wb、Wc组成开口三角形;CA、CB、CC为三相导线对地电容。若在PT开口三角端注入一恒定电流i0,就会有3个大小相等、相位相同的电流i1、i2、i3从PT的高压侧流出,这3个电流将分别在PT三相的绕组电阻R、漏抗XL和导线对地电容上产生压降。根据图2所示的PT等值电路,可得:
一般地,三相PT的参数(绕组电阻R和漏抗XL)是对称的,而且三相导线对地电容CA、CB、CC也是基本相等的,因此三相电流i1、i2、i3分别在三相PT与导线对地电容中产生的压降是基本相同的,即uA=uB=uC,这时在PT开口三角端就可以测到一零序电压u0。
将式(2)代入式(1),得:
式(3)中:i0、u0为输入和输出变量;R、L和C为待辨识参数。这样,问题就转化为参数辨识问题[3]。
1.2 递推式最小二乘法
递推式最小二乘法的计算方法是这样的。一般最小二乘法计算系统参数的公式为
那么,在进行第N+1次数据观测后,有
经过推导可以将递推算法的公式总结为:
其中:Pi=(XiTXi)-1,γi+1=(1+XTi+1Pi Xi+1)-1,KN+1=γN+1PNxN+1为增益矩阵[4]。
式(3)描述了系统输入i0,系统输出u0与系统参数R,L,C的关系。但是其中电气量的微分,积分算子,不易进行数值化和系统实现[5]。所以,这里对公式(3)进行了变形。首先,对式(3)由时域变换到频域,微分算子d/dt变换为jω,积分算子∫dt变换为1/jω,得到:
式中:是PT开口端的电压频域变量;是PT开口端的电流频域变量;ω是注入信号的角频率。接着,对上式移到左边,公式两边取虚部,进一步变形为
其中:为阻抗模,θ为阻抗相角,n1、n2分别是电压互感器PT的高、低压绕组匝数。
这里信号注入法的测量对象是ω,U0和I0。通过测量PT开口端特定频点处的U0、I0就可以求得Z和sinθ。获得PT等效回路的电路模型,测出了ω,Z和sinθ,接下来的工作就是运用最小二乘法对电路模型中的L,C系统参数进行辨识,得到电力系统接地电容。递推式最小二乘法辨识的目标是获得系统参数L,C,记为Φ=(L,C),系统输入为注入信号的角频率ω和其倒数1/ω,记为X=(ω,1/ω),系统输出为等效阻抗的模Z和相角θ,记为Y=Zsinθ。递推式最小二乘法的算法流程如图3所示。
递推过程是这样的:在PT开口三角侧注入角频率ω逐渐变化的电流信号。在10~200 Hz的范围里,经过大量仿真实验和计算比较,最终选择扫频范围10~181 Hz之间的20个频点,各频点间隔9 Hz,作为注入信号的频率。不断地读取系统输入X和相应的输出Y,i时刻读取的一对输入、输出记为(Xi,Yi)。每次读取一对(Xi,Yi),用最小二乘法计算出当前迭代的系统参数iΦ。如此循环,直到系统参数Φn与系统参数Φn+1之差的绝对值小于一个非常小的正数ε,那么迭代结束,取得系统辨识参数Φn=(Ln,Cn)。最后,计算出电网电容电流Ic=3ωCU相。
2 电容电流辨识算法实现
最小二乘法辨识电容电流的实现主要是完成辨识算法程序,并将其加载到DSP芯片中。辨识电容电流程序的目标器件选用TMS320LF2407芯片,采用C语言进行开发,不仅能够方便地实现算法的编程,而且也能够达到较高的辨识速度,不会影响电容电流辨识的实时性。在前面设计的辨识算法的基础上,可以比较容易地制定出辨识电容电流的程序流程,如图4所示。
其中,程序运行的终止条件可以根据系统要求和实际情况而定。如果对辨识电力系统电容电流的反应时间有规定,那么可以设定迭代次数不超过N。经过多次仿真实验和计算,这里N一般可取10,就可以使系统辨识结果具有较高精度。如果对辨识电力系统电容电流的精度有较高要求,那么可以设置辨识系统参数Φn和下一次的参数Φn-1之差的绝对值小于一个非常小的正数ε,那么迭代结束。一般地,ε可以取10-4,所用的辨识迭代次数通常在10次以内。另外,最小二乘法程序中涉及到对矩阵的较多操作,比如辨识算法中PiXi+1是n×n方阵和1×n向量相乘,(I-Ki+1XTi+1)P是两个n×n方阵的乘积,Pi=(XiTXi)-1需要进行矩阵求逆运算。程序中对各种类型的矩阵运算进行了单独编程,并且用函数进行实现,从而能将辨识过程进行明确划分,更好地执行辨识流程。
3 电容电流辨识算法仿真
这里,我们给出CCS2.20中TMS320LF2407芯片完成电容电流辨识仿真的过程和结果。CCS2.20是TI推出的用于开发其DSP芯片的集成开发环境。它采用Windows风格界面,集编辑、编译、链接、软件仿真、硬件调试及实时跟踪等功能于一体,极大地方便了DSP程序的设计与开发。对电压互感器PT回路的建模采用Matlab实现,仿真并获得PT开口端的等效回路的各次电流Ii和电压Ui,从而计算得到各个角频率相应的电路阻抗值Zsinθ。根据各次电流Ii和电压Ui,在CCS2.20中仿真运行递推式最小二乘法算法程序,然后我们就可以在CCS2.20中设置断点,监测到程序辨识到的接地电容,从而验证DSP的电容电流辨识程序的正确性。
在PT开口端输入10~181 Hz,频率间隔为9Hz,幅度为1 A,相位为0°的20个电流信号Ii。在每输入一个电流信号之后,根据PT开口端测量到的电压值,就可以计算等效阻抗值Zsinθ。各个角频率对应的电路阻抗值Zsinθ见表1所示。
然后,把电流注入法检测到的各次电路等效阻抗Z和相角θ作为输入在CCS2.20中运行递推式最小二乘法程序,检测接地电容C。表1中列出了辨识的结果。由表1可见系统对电感L的辨识误差从第三迭代开始都在1%以内,而对电容C的辨识误差则在第一次辨识就达到1%以内。
由此,我们可以看到,这里实现的递推式最小二乘法辨识电容电流具有较高的辨识效率,而且相比一般最小二乘法,可以使系统计算开销降低,辨识速度提高。
4 结论
本文实现的接地电容电流辨识的实现方案具有一定创新性。这里辨识算法首先对电力系统加入了电压互感器PT,并且建立了电力系统的等效模型和系统方程。辨识所需的电力系统状态的获取是通过信号注入法获得的,而对系统参数的辨识中使用的是递推式最小二乘法。本文在接地电容电流辨识中创新性地对PT等效回路的电路方程进行了变换,使得其中的系统输入、输出变量的采集、获取比较容易实现。此外,本文电路系统参数的求解运用了递推式最小二乘法,它的特点是能够实时、在线辨识系统参数,满足消弧线圈实时补偿电力系统电容电流的要求。
参考文献
[1]孙结中,蒋学军.配电系统电容电流测试仪的研制[J].广西电力,2006(2):1-4.SUN Jie-zhong,JIANG Xue-jun.Development of Capacitive Current Tester in Distribution System[J].Guangxi Electric Power,2006(2):1-4.
[2]赵正军,姜新宇.信号注入法在配电网电容电流测量中的研究[J].广东电力,2005,17(6):25-28.ZHAO Zheng-jun,JIANG Xin-yu.Research of Signal Injection Method used in the Measurement of Capacitive Current in Distribution Network[J].Guangdong Electric Power,2005,17(6):25-28.
[3]徐宁寿.系统辨识技术及其应用[M].北京:机械工业出版社,2003.XU Ning-shou.System Identification Technology and Its Application[M].Beijing:China Machine Press,2003.
[4]陈自宽,母国光.相关数据集的最小二乘处理方法[J].数据采集与处理,2007,11(1):66-68.CHEN Zi-kuan,MU Guo-guang.Least Square Method Sets to Deal with the Related Data[J].Data Acquisition and Processing,2007,11(1):66-68.
分段最小二乘法 篇6
MLS法是对传统加权最小二乘法的进一步推广, 克服了最小二乘法的不足, 突出了加权函数的紧支性, 据此建立的拟合函数, 具有局部近似的特点。MLS法的基本理论本文不再赘述详见参考文献[4]。本文针对具体二维和三维算例分析了MLS的局部近似性, 为其他研究领域中引入该方法解决复杂问题提供了一种新的思路。
2 算法流程
结合移动最小二乘法的特点, 给出了基于MATLAB软件编写的MLS法的曲线曲面拟合程序, 程序设计流程如图1所示。
3 权函数的影响分析
由图2所示:权函数为恒定常数, 样本点的权函数值均为1, 在求解过程中系数a (x) 不再是自变量x的函数, 而是在整个区域内a (x) 为常数, 可以看出系数a (x) 与权函数的选取有关, 因此, 曲线拟合就退化成最小二乘法的曲线拟合。这样权函数形成的近似为线性拟合。
由图3所示:权函数为紧支常数, 只有在中心点x (28) 5.0权函数存在非零区域, 非零区域尺寸即为节点权函数支持域尺寸, 这样在任意一点x的支持域内均有两个节点对其有贡献 (n (28) m (28) 2) , 权函数为1, 其他节点的权函数为零, 权函数构造的近似就变成了分段线性插值, 类似有限元法的函数逼近。由图4所示:权函数为紧支且光滑的函数 (三次样条函数) , 计算点x的支持域内包含的节点数大于基函数的个数 (n (29) m) , 为使局部近似误差最小, 加权平方和J (x) 对a (x) 求偏导, 最终得到MLS的近似位移, 这样的求解过程用于EFG方法中, 即使是线性基函数, 近似值也是连续光滑的, 因为MLS近似继承了权函数的连续光滑性。因此, 本部分的分析, 选取三次样条函数。
4 曲线拟合
采用常数基函数]1[, 拟合结果如图5所示。
采用一次基函数[1x], 拟合结果如图6所示。
采用二次基函数[1xx2]拟合结果如图7所示。
由图7可以看出, MLS拟合函数及其一阶导数拟合的效果比较好, 二阶导数拟合呈现稍微的震荡。为了更加精确的分析拟合结构, 引入相对误差范数公式:
精确解公。式中:N为样本点个数, yN和yE分别为对应的数值解和
在移动最小二乘法中基函数的选取也有一定的规律, 从以上分析, 结合表1可以看出:同一个函数采用不同基函数, 对计算结果有直接的影响, 线性基的基函数的次数越高, 逼近效果越好且数值计算的相对误差也越小, 二次基函数能够得到良好的拟合效果。
5 曲面拟合
采用MLS法拟合以下函数:
设定区域在区域W上均匀布置8 (8个节点, 图8显示了精确值的曲面分布图, 利用MATLAB中Griddata函数插值得到的。为进行曲面拟合, 采用具有C0 (W) 连续的高斯型权函数 (b=30.) 和线性基函数[1xy], 在区域上均匀划分25 (25个规则的节点网格, 由MLS法公式, 求出各节点的近似函数值, 绘制网格曲面, 如图10所示, 从图中可以看出, 用移动最小二乘法拟合的曲面比较光滑, 并与Griddata插值形成的曲面相比, 如图9所示, 曲面光滑准确, 从而验证该方法的有效性。
为了衡量拟合的质量, 定义拟合的相对误差:
其中:zEXACT (x, y) 是精确解, zMLS (x, y) 是使用移动最小二乘法计算值, 本文得到的拟合误差为:L=2.1526, 这说明MLS拟合具有比较高的精度, 虽然MLS拟合的计算量稍微大一些, 和其精确性相比, 是可以接受的。
6 结语
移动最小二乘法虽然不具备传统有限元或边界元形函数所具有的插值特征, 即不满足Kronecker函数的属性, 给边界条件处理带来了困难, 但是它采用了较低的基函数通过选取适当的权函数, 获得具有较高连续性和相容性的形函数, 数值精度较高, 数据简单等特点, 是其他数值近似方法无法比拟的, 对以后的数值近似的研究具有重要的参考价值。
参考文献
[1]王能超.数值分析简明教程[M].北京:高等教育出版社, 2003.
[2]常青, 余湘娟.软土次固结变形特性的试验研究[J].南水北调与水利科技, 2008 (1) :148-150.
[3]张荣.利用VB实现水泵性能曲线的自动绘制[J].水利科技与经济, 2006 (1) :62-64.
分段最小二乘法 篇7
由于岩土介质都是不均匀、不连续的,在其中传播的波动现象也非常复杂,波的反射与折射、介质体内的内摩擦导致能量不断被吸收而发生能量散耗现象,使得地震波向外传播的过程是一个指数不断衰减的过程。
由于爆源的复杂性,以及传播介质体物理力学性质和地质结构构造的多样性,使得爆破地震波具有随机波非重复性的特性,因此,爆破地震波是一种随机波。爆破地震波不仅在振幅上随时间做复杂的变化,而且波的频率构成和持续时间等也随环境条件、爆心距、爆破规模以及地层地质等的影响表现为极为复杂的现象。
1 回归分析原理
利用《爆破安全规程》预测爆破地震动强度的萨道夫斯基的经验公式,并利用最小二乘法对其进行回归分析,其公式的具体形式是:
对式(1)两边同取自然对数得:
得到:
式(3)中x是可控或能够精确测量的变量,若使x取一组不完全相同的值xi(i=1,2,3,…,n),进行计算分析,就得到与之对应的观测值yi(i=1,2,3,…,n)。(x1,y1),(x2,y2),…,(xn,yn)就是一组样本。显然可以通过这组样本估计y的数字特征。其中数学期望μ=μ(x)是重要的数字特征,函数μ(x)称为y对x的回归。如果y与x的关系是线性的,则称为一元线性回归。一元线性回归的数学模型:
其中,εi(i=1,2,3,…,n)为一组相互独立且服从于同一正态分布N(0,σ2)的随机变量。在上述条件下,y是服从正态分布N(βo+βx,σ2)的随机变量,在一元回归模型下有:
称它为y对x的线性回归方程,其图形称为回归直线。其中βo和βx称为回归方程的回归系数。求βo,βx的估计量bo和bx,也就是要确定一条经验回归直线y°=bo+bx来近似表达y与x的关系以下利用最小二乘法进行分析
对于每一个xi,都可由式y°=bo+bxi确定一个回归值y°i,回归值与观察值之差yi-y°i=yi-bo-bxi刻画了与回归直线的偏离程度。显然,对于所有的xi,若y°i与yi的偏离越小,则就认为回归直线和所有试验点拟合的越好。显然,全部观察值yi与回归值y°i的偏差平方和刻画了全部观察值与回归直线的偏离程度。
最小二乘法是使得偏差平方和Q(bo,b)达到最小的一种确定bo和b的方法。因此求bo和b就成为求Q(bo,b)最小值的问题。由于Q(bo,b)是bo和b的非负二次函数,所以它存在最小值。根据微积分学中的极值原理,要求的估计值bo和b是下列方程组的解:
式(6)方程组称为正规方程组,它可转化为:
解方程组得:
2 监测数据分析
以下利用所测得的数据对公式进行回归分析,建立爆破地震振动速度与爆心距之间的线性关系。选取某地下工程施工时对现场进行的部分监测数据。
根据表1中的数据,利用最小二乘法回归分析的原理进行线性回归,求得:K=186.49,α=1.434 7。故可得出爆破时振动速度沿不同方向传播衰减规律,如下所示:
其中,V为测点处的地表面质点振动速度,cm/s;Q为单段最大装药量为传播距离
在这一地下工程爆破中,利用最小二乘数法能够容易得出振动速度与装药量、传播距离之间的具体关系。并且利用式(9)对以后爆破进行指导都能起到较理想效果。故利用最小二乘数回归分析的爆破经验公式的方法是可行的
3 结语
本文通过对某地下工程监测数据利用最小二乘法对萨道夫斯基的经验公式进行回归分析得出其具体参数值,并且能够在具体工程中得到很好的利用。因此说明该方法能很好地解决爆破震动参数的回归分析进一步为工程提供爆破理论依据
摘要:针对爆破的不确定性与难估量性以及传播介质特性的复杂性,利用最小乘法原理对爆破震动进行回归分析,得出爆破震动速度与装药量、传播距离的关系,实践证明该方法能很好地解决爆破震动参数的回归分析。
关键词:爆破,最小二乘法,回归分析,爆破震动
参考文献
[1]姚勇,何川.董家山隧道小净距段爆破控制研究[J].公路,2005(11):212-216.
[2]秦根杰.某隧道扩建中的控制爆破[J].西部探矿工程,2005(sup):228-229.
[3]谭忠盛,杨小林,王梦恕.复线隧道施工爆破对既有隧道的影响分析[J].岩石力学与工程学报,2003,22(2):281-285.
[4]郑际汪,陈理真.爆破荷载作用下隧道围岩稳定性分析[J].矿山压力与顶板管理,2004(4):53-55.
分段最小二乘法 篇8
关键词:MIMO雷达,单调约束方向图设计,发射信号设计,凸优化
MIMO (multiple-input multiple-output )雷达是一种新体制雷达,它在信号检测、参数估计、空间分辨力等方面具有诸多优点[1,2,3]。根据收发阵元的间距可分为统计MIMO雷达和相干MIMO雷达[1]。波形优化设计是MIMO雷达的一个重要内容。目前MIMO雷达波形优化设计主要有两种思路:只给出发射信号的协方差矩阵;设计出具体的时间序列。在设计发射信号波形时可利用发射信号间的部分相关性,使雷达系统辐射的能量尽可能集中在感兴趣的区域内[4,5,6,7,8]。文献[4,5,6]提出一种半正定规划算法,主要包括:方向图匹配设计、最小旁瓣设计,可在多项式时间内求解。文献[7]提出发射方向图可以由一组基波束合成的思想,用此法综合出来的发射方向图具有较低的副瓣和较小的互相关。文献[8]采用梯度算法把方向图最优问题转化为一个四次多元多项式的最小值问题,极大降低了运算量,但不能保证收敛到全局最小值。以上研究仅仅只得到信号的协方差矩阵R,都没有对如何从R综合出具体的发射信号而进行深入研究,这类文献非常少。文献[9]提出了如何在给定的R的条件下设计恒模的波形的最优问题。
本文针对最小旁瓣方向图设计在主瓣范围内有时会出现较大的波动,提出主瓣单调约束方向图设计算法,然后利用凸优化算法CVX[10] 求解出最优的信号协方差矩阵R。同时基于将发射波形看作伪随机信号的思想将协方差矩阵R分解,提出了一种利用最小二乘法计算发射信号序列及其权值的算法。
1 发射信号模型
假设MIMO雷达系统由Mt个发射阵元、Mr个接收阵元组成的均匀线阵。xm(n)表示第m个阵元的发射信号,θ表示空间方位角,d表示阵元间距,λ表示载波波长。发射信号是窄带信号且不考虑衰减,则在远场处的辐射信号为[1]:
式(1)中,发射信号: x(n)=[x1(n) … xMt(n)]T。
方向导向矢量:
平均功率:p(θ)=aH(θ)Ra(θ) (2)
式(2)中,R=E{x(n)xH(n)}是发射信号的协方差矩阵。
由式(2)可知如何合理分配MIMO雷达发射能量的问题在数学上就变成如何设计信号x(n)的协方差矩阵R。实际中还必须根据矩阵R综合出具体的信号x(n)。
2 方向图设计
2.1 最小旁瓣方向图设计
通过选择合适的协方差矩阵R使得方向图的旁瓣最小,即,在预定的区域内最小化旁瓣电平并达到预定的3 dB主瓣宽度。
在数学上就是要解如下的最优问题[4]:
式(3)中,t :中间变量, θ0:主瓣波束指向 ,Ω:旁瓣范围 ,a(θ):方向导向矢量,θ1、θ2:3 dB主瓣宽度边界,Rmm:R的主对角线元素,R≥0:表示矩阵是半正定的。cp:信号功率。
2.2 主瓣单调约束方向图设计
式(3)在主瓣范围内有时会出现较大的波动,对目标的探测不利。这是因为原算法只对3 dB主瓣宽度边界点的幅度作了约束,主瓣内各点间的幅度没有相对约束。本文在此基础上进行改进,即方向图在半个3 dB主瓣宽度内受单调性约束,使得主瓣内无较大波动。
在式(3)的基础上增加以下约束条件:
aH(θ1)Ra(θ1)<aH(φl)Ra(φl), ∀φl∈Ω1 (4)
式(4)中,Ω1:3 dB主瓣宽度(Ω1=θ1-θ2)。
式(4)表示在3 dB主瓣宽度边界点的功率小于主瓣内任意一个方向的功率。
式(3)、式(4)的优化问题都是能收敛到全局最优解的,可以用凸优化算法(CVX)求解得到全局最优数值解R。
3 发射信号设计
3.1 协方差矩阵的伪随机信号分解
仅得到信号的协方差矩阵R是不够的,实际中需要知道具体的发射信号x(n)。本文讨论的发射波形是BPSK信号。这里把发射波形看作伪随机信号[9],这样就可以把R表示成伪随机序列的线性组合。当发射阵元个数为Mt时,si有2Mt种可能,即R可以表示成:
3.2 最小二乘法计算发射序列权值
因为R是实对称的且主对角线元素相等,所以从式(5)中最多可得到由(M
式(6)中,
利用最小二乘法求解线性方程组DP=c,就可解出P,即各序列的权值。得到与pi对应的序列si。单个阵元在一个脉冲周期内发射的信号与si的关系:
由式(7)即可得到x(n)。
4 仿真试验及结果分析
阵元间距为半波长的均匀线阵,Mt=12,3 dB主瓣范围为[-10°,10°],波束指向0°,信号功率:cp=10。
仿真1:最小旁瓣和主瓣单调约束方向图仿真
对原算法(式(3))和改进的算法(式(3)与式(4))利用CVX进行仿真:
a) Ω=(-90°,-20°)∪(20°,90°),如图1所示;
b) Ω=(-90°,-15°)∪(15°,90°),如图2所示。
由图1、图2可知,文献[4]中提出的最小旁瓣方向图设计方法在3 dB主瓣宽度内有时会出现较大波动。原算法在3 dB主瓣宽度内出现较大波动(峰峰值超过3 dB)时,本文的改进算法能有效抑制波动、峰值旁瓣仅高0.5 dB(图1);Ω=(-90°,-15°)∪(15°,90°)时,两者的方向图完全一样(图2), 即原算法在3dB主瓣宽度内无较大波动时,改进算法能保持其性能。
仿真2:最小旁瓣和发射信号方向图仿真
利用CVX根据式(3)与式(4)进行仿真,其中Ω=(-90°,-15°)∪(15°,90°),如图3中实线波形;针对用主瓣单调约束方向图设计法得到的协方差矩阵R,利用最小二乘法求解式(6)可得各序列的权值及相对应的序列,如表1所示;由综合出的伪随机发射信号仿真,如图3中虚线波形。
表1给出了各阵元在一个脉冲周期内的发射信号与si及其权值的关系。最后一列的和∑pi=0.93,不等于理论值1,这是由于由最小二乘法得到的解是近似解。图3中,发射信号方向图比最小旁瓣方向图的3 dB主瓣宽度展宽约0.3°,在主瓣波束指向上发射信号方向图比最小旁瓣方向图低约0.43 dB,发射信号方向图比最小旁瓣方向图的峰值旁瓣高约0.53 dB,两者比较接近,这说明本文提出的信号综合算法对BPSK信号是有效的。
5 结束语
本文针对最小旁瓣方向图设计法在主瓣范围内有时会出现较大的波动,提出主瓣单调约束方向图设计法,该法能有效抑制主瓣波动。同时提出了一种利用最小二乘法由协方差矩阵反推具体发射信号序列及其权值的信号综合算法,综合出的发射信号序列的方向图与由协方差矩阵得到的方向图比较接近,说明该算法对BPSK信号是有效的,为将这些波形优化算法应用到工程实际中提供了可能。
参考文献
[1] Fishler E,Haimovich A,Blum R,et al.MIMO radar:an idea whosetime has come.Proceedings of the IEEE Radar Conference,April2004:71—78
[2] Robey F C,Coutts S,Weikle D,et al.MIMO radar theory and expri-mental results.Proceedings of the 38th Asilomar Conference on Sig-nals,Systems and Computers,Nov.2004:300—304
[3] Fishler E,Haimovich A,Blum R,et al.Spatial Diversity in RadarsModels and Detection Performance.IEEE Transactions on Signal Pro-cessing,Mar.2006;54(3):823—838
[4] Jian Li,P.Stoica,Yao Xie.On Probing Signal Design for MIMO Ra-dar.IEEE Transactions on Signal Processing,2007;55(8):4151—4161
[5] Antonio G S,Fuhrmann D R.Beampattern Synthesis For WidebandMIMO Radar Systems.IEEE Conference on Radar,2005:105—108
[6]胡亮兵,刘宏伟,刘保昌,等.MIMO雷达发射方向图匹配和波形优化方法.西安电子科技大学学报(自然科学版),2009;36(6):1021—1026
[7]胡亮兵,刘宏伟,杨晓超,等.集中式MIMO雷达发射方向图快速设计方法.电子与信息学报,2010;32(2):481—484
[8] Aittomaki T,Koivunen V.Beampattern Optimiza-tion by Minimiza-tion of Quartic Polynomial.2009 IEEE/SP 15th Workshop on Statisti-cal Signal Processing:437—440
[9] Fuhrmann D R,Antonio G S.Transmit Beamforming for MIMO RadarSystems using SignalCross-Correlation.IEEE Transactions on Aero-space and Electronic Systems,2008;44(1):171—186
【分段最小二乘法】推荐阅读:
最小二乘法分析07-05
最小二乘法模型12-02
最小二乘法拟合直线05-21
偏最小二乘法模型05-26
递推最小二乘法06-04
线性回归最小二乘法06-19
最小二乘法曲线拟合07-11
偏最小二乘法回归模型07-05
非线性最小二乘法拟合06-24
改进的偏最小二乘法11-29