雷达物位计软件算法说明

1、1. 法兰距离计算1.1. 流程图1.2. 信号加窗信号加窗用于减小频谱泄露，可选择三种窗函数之一：矩形窗、汉宁窗、布莱克曼窗。假设中频信号电压采样数组为vN，采样点数为N（N=1199）；加窗实际上是构建一个N点的数组wN，将vN和wN进行点乘运算；信号加窗后的输出数组xN可表示为：xn=vn*wn （0n<N） (1-1)1) 矩形窗矩形窗的wN是一个全1的数组，相当于对vN不进行加窗处理。2) 汉宁窗汉宁窗的wN表示如下：wn=0.5-0.5cos2nN （0n<N） (1-2)3) 布莱克曼窗布莱克曼窗的wN表示如下：wn=0.42-0.5cos2nN+0.08cos4nN

2、 （0n<N） (1-3)1.3. 信号补零将N（N=1199）点的数组xN扩展成M（M=4096）点的数组yM：yn= xn （0n<N） 0 （Nn<M） (1-4)1.4. FFT运算对M点的数组yM进行实数快速傅里叶变换并求模，具体算法请参考数字信号处理相关书籍；TMS320F28335有现成的库函数RFFT_f32(RFFT_F32_STRUCT *)、RFFT_f32_mag(RFFT_F32_STRUCT *)供调用， FFT运算的输出为数组YM/2。1.5. 频谱峰值探测根据MATLAB信号处理工具箱中Peak Finder的算法思路：处于两个波谷间且大于波谷

3、一定阈值的点即为波峰；对应的流程图为：频谱峰值探测的输出为波峰索引数组。1.6. 回波筛选有效回波必须满足一定的位置条件和幅值条件。系统的距离分辨单元为：D=C2B*11994096 (1-5)式中，C 真空中光速 B 扫频带宽则回波的位置和幅度可以表示为： D=D*i-TCLA=Yi (1-6)式中，i 波峰索引 TCL TCL长度1) 位置条件有效回波位置必须处于盲区和罐底之间，将处于该范围之外的回波剔除。2) 幅值条件有两种幅值条件：统一阈值、ATP阈值。统一阈值是将峰值小于某阈值的回波剔除；ATP阈值是由位置阈值构成的一条折线，将峰值处于折线下方的回波剔除。回波筛选的输出为回波索引数组

4、。1.7. 谱估计对回波索引数组中的每一个回波D，Yi,根据该回波前后各1个点D-D，Yi-1、D+D，Yi+1的值，利用二次曲线拟合法估算回波的真实位置为：D0=D+Yi-1-Yi+12Yi-1+Yi+1-2Yi*D (1-7)回波的幅度为：A0=D0-DD0-D-DYi-1-2D0-D-DD0-D+DYi+D0-DD0-D+DYi+12D2 (1-8)1.8. 回波跟踪如果选择首回波算法，则回波跟踪过程可以用如下状态机表示：如果选择回波搜索算法，则回波跟踪过程可以用如下状态机表示：回波跟踪的输出为液面状态、位置和幅度。1.9. 滤波处理有三种滤波方法：跳跃滤波、阻尼、最小二乘法滤波；其中阻

5、尼不能和最小二乘法滤波同时使用，但跳跃滤波能和阻尼或最小二乘法滤波同时使用。1) 跳跃滤波跳跃滤波用于减少液面的大幅度快速抖动，跳跃滤波用到两个参数：跳跃滤波阈值JFilt_Thresh、跳跃滤波延迟JFilt_Delay，其算法思路如下：i. 如果当前液位与前一时刻液位输出值相距大于JFilt_Thresh，跳跃滤波计数器累加；否则跳跃滤波计数器清零, 用当前液位值作为当前时刻的输出。ii. 如果跳跃滤波计数器大于JFilt_Delay，用当前液位值作为当前时刻的输出；否则用前一时刻液位输出值作为当前时刻的输出。2) 阻尼阻尼用距离滤波因子来度量，施加阻尼后系统的输出为：Dn'=Dn

6、-1'+Dn-Dn-1'*DistFiltFactor (1-9)式中，Dn' 当前时间输出值 Dn-1' 前一时刻的输出值 Dn 当前时刻的测量值 DistFiltFactor 距离滤波因子距离滤波因子和阻尼时间常数之间可用如下关系式表示：DistFiltFactor=1-e-TT (1-10)式中，T 时间常数，单位为s T 液位计算的时间间隔(0.1s)3) 最小二乘法滤波最小二乘法滤波用于快速跟踪液位的变化，其算法步骤如下：iii. 利用前N（N=10）个时刻的输出值D'0、D'T、D'2T、D'N-1T，根据最小二乘法，

7、拟合一条输出随时间变化的直线：D''t=a+b*t (1-11)iv. 利用拟合直线，计算当前时刻的预测值：D''NT=a+b*NT (1-12)v. 根据当前时刻的预测值和测量值，计算当前时刻的输出值： D'NT=D''NT+DNT-D''NT*DistFiltFactor (1-13)式中，DNT 当前时刻的测量值。注：最小二乘法线性拟合见附录A。2. 液位计算以用户定义的液位零点为参考点，计算出的液位值为：Level=TankHeight_R-Flange_D-OffsetDist_G+CalibrationLeve

8、l+LevelOffset (2-1)式中，TankHeight_R 储罐高度 Flange_D 上一步计算出来的法兰距离 OffsetDist_G 距离偏移量 CalibrationLevel 液位校准值 LevelOffset 液位修正量3. 距离计算以用户定义的距离零点为参考点，计算出的距离值为：Ullage=Flange_D+OffsetDist_G-CalibrationLevel+UllageOffset (3-1)式中， UllageOffset 距离修正量4. 体积计算有4种体积计算方法：插值法、理想球罐法、理想垂直圆柱罐法、理想水平圆柱罐法。1) 插值法使用拉格朗日插值法进行

9、体积计算；假设体积插值表可表示为L0，V0、L1，V1、LN-1，VN-1，并且满足L0<L1<<LN-1、V0<V1<VN-1；首先将液位转换为以插值表参考零点为基准：Level_V=TankHeight_R-Flange_D-OffsetDist_G+CalibrationLevel+LevelOffs (4-1)式中，LevelOffs 插值表液位偏移量 线性插值寻找点Li，Vi、Li+1，Vi+1，使得LiLevel_V<Li+1，使用下式计算体积值：V'=Level_V-Li+1Li-Li+1*Vi+Level_V-LiLi+1-LiVi+

10、1 (4-2)如果Level_V<L0，用L0，V0、L1，V1进行计算；如果Level_VLN-1，用LN-2，VN-2、LN-1，VN-1计算。 二次曲线插值法寻找点Li，Vi、Li+1，Vi+1，使得LiLevel_V<Li+1，使用下式计算体积值：V'=Level_V-Li+1Level_V-Li+2Li-Li+1Li-Li+2*Vi +Level_V-LiLevel_V-Li+2Li+1-LiLi+1-Li+2*Vi+1 +Level_V-LiLevel_V-Li+1Li+2-LiLi+2-Li+1*Vi+2 (4-3)如果Level_V<L0，用L0，V0

11、、L1，V1、L2，V2计算；如果Level_VLN-2，用LN-3，VN-3、LN-2，VN-2、LN-1，VN-1计算。插值法的体积输出值为：Volume=V'+VolumeOffs (4-4)式中，VolumeOffs 体积修正量2) 理想球罐法Level_V=TankHeight_R-Flange_D-OffsetDist_G+CalibrationLevel+ZeroLevelToBottomVolume=163L1-2Level_V*Level_V2+VolumeOffs (4-5)式中，ZeroLevelToBottom 液位零点到罐底的距离 L1 球罐直径3) 理想垂直

12、圆柱罐法Level_V=TankHeight_R-Flange_D-OffsetDist_G+CalibrationLevel+ZeroLevelToBottomVolume=14L12*Level_V+VolumeOffs (4-6)式中，L1 圆柱体直径4) 理想水平圆柱罐法Level_V=TankHeight_R-Flange_D-OffsetDist_G+CalibrationLevel+ZeroLevelToBottomVolume=14atan2L12-L1-2Level_V2，L1-2Level_V*L12 - 14L1-2Level_V*L12-L1-2Level_V2+ Vo

13、lumeOffs (4-7)式中，L1 圆柱体直径 L2 圆柱体高度5. 液位速率计算液位速率由下式计算：LevelRate=Level_V-Level_V'T (5-1)式中，Level_V' 上一时刻的液位值 T 液位计算的时间间隔(0.1s)附录A 最小二乘法线性拟合最小二乘法线性拟合就是将一组符合Y=a+bX关系的测量数据，用计算的方法求出最佳的a和b。假设直线方程的表达式为：要根据测量数据求出最佳的a和b。对满足线性关系的一组等精度测量数据（xi，yi），假定自变量xi的误差可以忽略，则在同一xi下，测量点yi和直线上的点a+bxi的偏差di如下： 显然最好测量点都在直线上（即d1=d2=dn=0），求出的a和b是最理想的，但测量点不可能都在直线上，这样只有考虑d1、d2、dn为最小，也就是考虑d1+d2+dn为最小，但因d1、d2、dn有正有负，加起来可能相互抵消，因此不可取；而|d1|+ |d2|+ |dn|又不好解方程，因而不可行。现在采取一种等效方法：当d12+d22+dn2对a和b为最小时，d1、d2、dn也为最小。取（d12+d22+dn2）为最小值，求a和b