加速度监测数据校正(1)

加速度检测仪数据校正

摘要

本文针对声屏障检测仪内部加速度检测器的数据校正，利用题中所给的加速度数据，在校正的过程中，我们结合物理运动规律、离散型随机变量等知识，建立了加速度-速度、加速度-位移仿真模型，运用组合辛普森算法、卡尔曼滤波器、正负补偿法等方法消除了随机误差和系统误差的影响，使物体运动规律符合实际情况，并将模型推广运用到其他领域。

针对问题一中声屏障的速度、位移的仿真计算以及误差分析，我们基于s

-,基

a-

a

v

本物理公式建立计算声屏障运动速度和位移的仿真模型，采用高精度的组合辛普森积分公式，通过仿真计算分别绘制出三种情形下的速度、位移-时间关系图，将仿真图与理想情况相对比，我们发现存在系统误差和随机误差。因此我们分别从系统误差和随机误差2个角度对数据进行定性和定量的误差分析。

问题二中，利用问题一中速度和位移的数值积分计算模型和误差分析结果，以尽量消除系统误差与随机误差，使得速度和位移的计算结果基本符合物体运动事实为目标，对加速度数据进行校正。我们利用卡尔曼滤波对加速度数据进行降噪处理，然后利用正负补偿法消除系统误差，有效的校正了数据，将校正后的数据代入问题一所建模型中建议，发现物体运动规律符合实际情况，即最终速度为0，位移为一稳定值。

问题三中，改进后的加速度数据校正模型可以推广应用到生产生活中，如电梯的加速度测量、油井示功图位移测量技术、惯性导航系统、胎儿心率检测仪等等。

本文最大的特色在于利用精确度较高的组合辛普森算法，并且综合卡尔曼滤波法和正负补偿法分别减小了随机误差和系统误差的影响，具有一定的可靠性。

关键词：组合辛普森算法、卡尔曼滤波、正负补偿法、数据校正

一、问题提出

声屏障是一种控制铁路、公路、高速铁路等各种道路行车对周围环境的噪声污染有效措施之一，随着列车的大幅度加速，脉动风交替出现在列车两侧，从而引起对声屏障的拉压作用，声屏障发生摆动。正常状态下，声屏障的摆动应当在一定的范围内，当超过正常范围则需要对其进行加固维修。由于声屏障维修或重建费用高昂，故需声屏障检测仪对声屏障的工作状态进行检测，有针对性的对声屏障进行维修。

声屏障检测仪的工作原理是：通过内部的加速度传感器来记录车辆经过时声屏障振动而产生的加速度数值(密集采样)。将加速度数据通过数值积分，按照加速度-位移的物理公式将加速度数据转化为震动的位移，并通过震动位移对声屏障状态进行判断。

在试验中，传感器测得的数据通常会存在误差，误差包括系统误差、随机误差。其中系统误差，又称为固有误差，一般其存在是具有一定的规律性，是可以被分析掌握的；随机误差，又称为测量误差，一般它的出现是不具有规律并且不可避免的。由于误差的存在，在使用数值积分方法计算振动位移的过程中，就会累积较多的干扰，故而在测得数据后，需要经过系统误差校正、随机误差数据滤波等对数据进行校正。现在请建立数学模型解决如下问题：

1.建立适当的数学模型，基于加速度-速度和加速度-位移物理公式，通过数值积分的方法计算声屏障的速度、位移，并基于给定数据对模型进行仿真计算，判断声屏障检测仪是否存在明显误差，从随机误差、系统误差2个角度对数据进行误差分析；

2.基于速度和位移的数值积分计算模型和误差分析结果，建立数学模型来对加速度数据进行校正，要求能尽量消除系统误差与随机误差，使得速度和位移的计算结果基本符合物体运动事实；

3.对你所建立的数据处理方法和模型进行推广，所改进过的加速度检测仪除了可以用于声屏障监测以外，还可以应用于哪些场景，请结合改进方案阐述理由。

二、基本假设

1、假设：在三组数据采集过程中，加速度检测仪所处外界环境基本相同；

2、假设：在加速度检测仪工作的初始时刻,声屏障的加速度、速度、位移均为0；

3、假设：不考虑环境因素如温度，气候对加速度监测仪的影响；

4、假设：声屏障发生的形变总是弹性形变，侧向外力消除后可恢复原状并回复原位；

5、假设：试验过程中，仅考虑声屏障在水平方向的运动，即仅存在水平方向的加速度。

三、符号说明

四、问题分析

本文主要解决的是加速度检测仪数据校正的问题。

对于问题一，基于v a -和s a -物理公式，以及组合辛普森积分公式，通过matlab 程序给出基于给定数据的模型仿真。将实际值与理论值进行对比，判断声屏障检测仪是否存在明显误差。根据随机误差的正向分布，采用贝塞尔公式求出三组数据的偏差，判断随机误差的大小。其次，通过对三组加速度数据积分所得速度-时间图像、位移-时间

图像对比，总结出系统误差规律。

对于问题二，我们的研究思路如下：

对于问题三，总结上述对加速度检测仪所采集数据的校正方法，将其合理推广于加速度变化情况与之相似的场景中，并作合理阐述。

五、模型的建立与求解

5.1 问题一模型建立与求解

5.1.1 问题一的分析

首先，基于物理学基本公式建立计算声屏障运动速度和位移的数值积分模型；其次采用高精度的组合辛普森求积公式，通过仿真计算分别绘制出三种情形下的速度、位移-时间关系图；最后将仿真图与理想情况相对比，从系统误差和随机误差2个角度对数据进行误差分析

5.1.2.1数值积分模型的建立

用加速度仪的采样周期T 作为时间间隔，将声屏障的振动过程分成n 个大小相同的区间；

1Step ：分析每个区间上速度的增量（或者位移的增量），

第k 个区间],)1-[(kT T k 上速度的增量)(k v ?和位移的增量)(k s ?为

)()()()()1()1(??

?

??=?=???--kT

T k kT T k dt

t v k s dt

t a k v 2Step ：依据加速度-速度和加速度-位移物理公式得出数值积分模型，

t 时刻声屏障运动的速度)(t v 和位移)(t s 为

n k dt dt t a dt t v t s dt t a t v T t k kT T k T t k kT T k T t k kT

T k T t k kT

T k 1,2,... ))(()()()()(/1)1(/1)1(/1)1(/1

)1(=???

????

===∑?∑?∑?∑?=-=-=-=-

对于单方向运动从A-B 的过程1395=n ； 从C-D ，再从D-C 的双向运动过程3134=n ；

从E 点到F 点，再由F 到E ，并再重复一次的过程2397=n ； 5.1.2.2引用组合辛普森公式求解数值积分模型

由于加速度传感器采样频率较高，且采集次数较多，因此在不增加计算量并且保证

精确度足够高的条件下，可引用组合辛普森积分公式近似求解，过程如下： 设区间[]kT T k ,)1-(被等距节点M h t 2,,1,0, ==ααα，分为宽度为M

T

h 2=的M 2个子区间[]1,+ααt t ,M 2个子区间的组合辛普森公式表示)(t v 的定积分为：

[]∑∑∑?=-==--+++==T

t M M

T

T k t a h t a h T a a h dt t a t v /1

k 111122k )1()(34)(32)()0(3)()(αααα 同理可得)(t s 的定积分表达式为：

[]∑∑∑?-==-=-+++==111122/1

k )1()(34)(32)()0(3)()(M M T

t kT T k t v h t v h T v v h dt t v t s αααα 组合辛普森公式计算误差为 )()(180

),(E 4)

4(4

h O c f h T h f s =-

=，其中f 为)()(t v t a 或.

由于每个时间区间与整个运动全程相比是非常小的，因此我们假设)()(t v t a 、在区间[]kT T k ,)1-(范围内是呈线性变化的，这样即可求出)()(ααt v t a 、的函数表达式。 5.1.2.3模型的参数计算

采样周期s Hz

T 001.010001

==

设M =4，即将区间[]kT T k ,)1-(分为8个子区间， 则子区间宽度s 800012==

M T h , 节点8000

ααα==h t 5.1.2.4数值积分仿真结果

利用上述组合辛普森数值积分算法编程，在matlab 中进行仿真计算，得到三组数据下速度-时间关系图和位移-时间关系图：

（1）单方向B →A 运动的加速度-时间关系图

（2）单方向B →A 运动的速度-时间关系图

（3）单方向B →A 运动的位移-时间关系图

图5-1 单方向B →A 仿真结果

（1）从C D →→C 运动的加速度-时间关系图

（2）从C D →→C 运动的速度-时间关系图

（3）从C D →→C 运动的位移-时间关系图

图5-2 单方向C D →→C 仿真结果

(1)从E F E F →→→→E 运动的加速度-时间关系图

(2)从E F E F →→→→E 运动的速度-时间关系图

(3)从E F E F →→→→E 运动的速度-时间关系图

图5-3 从E F E F →→→→E 运动的仿真结果

5.1.2.4数值积分仿真结果的分析

由积分运算的性质知，原始数据中的微小误差会随着积分运算的过程不断累加形成

巨大误差，如图321-5、、

中速度-时间关系图和位移-时间关系图都与实际情况存在较大出入。

我们以单方向B →A 运动为例分析：

图5-4 单方向B →A 运动速度、位移-时间关系图

理想情况下，当速度-时间关系图中a 点之前b 点之后的区域加速度为0，速度也为0，曲线图应靠近虚线部分；对应的位移-时间关系图中c 点之后的位移达到最大值，位移曲线应当靠近虚线，而实际仿真的波形图显示结果并非如此，因此认为声屏障检测仪存在明显误差。下面我们就从系统误差、随机误差2个角度对误差进行定量和定性分析。

5.1.3误差分析

加速度测量值a '与理论值a 之间存在测量误差。测量误差一般是由随机误差δ与系统误差ε组成，即有：a a εδ'=++。

5.1.3.1随机误差的分析

随机误差也成为偶然误差，它的出现从表面上看是毫无规律的。在任何测量中，随机误差都是不可避免的，而且在同一条件下重复进行的各次测量中，随机误差的出现或大或小，或正或负，其大小和方向均不固定。随着测量次数的增加，正负误差可以相互补偿，误差的平均值将趋向于0。但就其总体来说,却具有某些内在的共性，利用概率论的一些理论和统计学的一些方法，可以掌握随机误差的若干规律，并对误差大小进行估算：正负误差出现的概率相等，小误差出现的概率大，大误差出现的概率小，总体按照误差大小服从正态分布。

由于密集采样过程中，样本容量均远大于30，故认为随机误差δ服从正态分布，亦称高斯分布，它满足的概率密度分布函数为：

()22

2

f δσδ-

=

在实际测量中，由于样本容量不能无限大，故此时的算术平均值不是真值，对标准误差σ的实际处理只能进行估算。利用数理统计理论，可以得到对偏差进行估计的贝塞尔公式：

∑=--=n

i i a a n 1

2)(11?σ

n 为样本容量，σ

?是σ的最佳估计值。借助excel 软件，可得到本次实验中的三组数据的标准偏差分别为：024779517.1?,706158082.0?,582910513.0?321===σ

σσ 由概率密度分布函数可知，测量值的随机误差出现在),3,2,1)(3,3(为标准误差?=??-j j j 范围内的概率为99.7%，即认为所测得的全部数据中，将有99.7%的数据其随机误差落

入置信区间。

5.1.3.2系统误差的分析

系统误差的出现一般是有规律的，在同一测量条件下，多次重复测量同一值时，观测值总往一个方向偏差，测量误差的绝对值和正负号在重复多次测量中几乎相同，都保持基本不变。产生系统误差的原因可能各不相同，但是它们具有共同的特点，即确定的变化规律。

由于加速度散点图图像密集不易区分判断，故这里采用组合辛普森积分公式对加速度积分后，得到的速度-时间图像进行分析。此外，用组合辛普森积分产生的误差可以表示为)(O 4h ，其中h 为一个极小的常数，又因为1000f Hz =，所以积分过程对速度以及位移的影响极小，故在此忽略不计。

图5-5 三组数据速度-时间关系比较图

根据图5-5速度-时间分布规律，不难看出，每经过一次单方向运动后的速度值总会减小，亦即加速度散点图中的波峰积分值小于波谷积分值。从物理学角度讲，即认为声屏障检测仪所测得的反向加速度值普遍高于正向运动时的加速度值，认为此类误差为定值系统误差。

5.2 问题二模型建立与求解

5.2.1 问题二的分析

此问利用问题一中速度和位移的数值积分计算模型和误差分析结果，以尽量消除系统误差与随机误差，使得速度和位移的计算结果基本符合物体运动事实为目标，对加速度数据进行校正。为消除数据中存在的误差，我们首先利用卡尔曼滤波法对加速度数据进行降噪处理，有效的消除随机误差，然后进行系统误差研究，对于本题中这种定值系统误差，因测量误差的绝对值和正负号在重复多次测量中几乎相同的特性，我们可以利用反向补偿法来消除系统误差，即改变原加速度方向是误差符号相反，从而抵消不变得系统误差。

5.2.2 问题二模型的建立

结合问题分析，我们须引用以下算法：

卡尔曼滤波法：一个动态系统的状态由于受到外部的干扰时，可以分为两个组成部分，一部分是由已知的运动方程正确的预测出来，另一部分可以看作是随机变量，对于线性离散时间系统可以分为两个部分，即使用n 维状态方程和m 维测量方程描述： ()(1)()

(x t A x t B u t w t =-++ (2.1)

()(1)()y t C x t v t =-+ (2.2) 其中：()x t 是t 时刻的系统状态，()u t 是t 时刻对系统的控制量，A 和B 是系统参数对于多模型系统，它们为矩阵，()y t 为t 时刻的测量值，C 为测量系统的参数，对于多测量系统，其为矩阵。()()w t v t 和分别表示过程噪声和测量噪声，它们被假设成高斯白噪声，他们的协方差分别是R Q ,。

首先我们要利用系统的过程模型，来预测下一状态的系统。假设现在的系统状态是k ，根据系统的模型，可以基于系统的上一状态而预测出现在状态： (|1)(1|1)(

x t t A x t t

B u k -=--+ (2.3) 式中，(|1)x t t -是利用上一状态预测的结果，(1|1)x t t --是上一状态最优的结果，()u k 为现在状态的控制量，如果没有控制量，它可以为0。

我们用p 表示协方差：

(|1)(1|1)p t t A p t t A Q '-=--+ (2.4) 式中，(|1)p t t -是(|1)x t t -对应的协方差， A '表示A 的转置矩阵，Q 是系统过程的协方差。

现在我们有了现在状态的预测结果，然后我们再收集现在状态的测量值。结合预测值和测量值，我们可以得到现在状态t 的最优化估算值(|)x t t ：

(|)(|1)()(()(|x t t x t t k g t z t C x t t =-+-- (2.5) 其中kg 为卡尔曼增益：

(|1)()(|1)p t t C kg t Cp t t C R

'

-=

'-+

为了要另卡尔曼滤波器不断的运行直到系统过程结束，我们还要更新t 状态下(|)x t t 的协方差：

(|)(|1)(1())p t t p t t Ckg t =--

利用卡尔曼滤波法，在MA TLAB 环境下进行降噪处理，下图为经过多次滤波后的散点图：

A-B 数据原始图

滤波后数据

C-D 数据原始图

E-F 数据原始图

滤波后数据

在数据进行有效消除随机误差之后我们利用正负补偿法对其进行系统误差处理，在声屏障的每次运动当中，根据加速度的变化规律知道，声屏障在速度增大和减小的过程中存在对应的加速度点，其值互为相反数。我们将其近似处理为弹簧振子的震动，如图：

假设N 点的加速度理论值为(),0t t a a >，实际值为t a '，误差为1a ?，则有：

1t t a a a '=+?。与其对称的点'N 加速度理论值为t a -，实际值为t a ''，误差为2a ?，则有：2t t a a a ''=-+?。此时N 与'N 实际值的绝对值之差为12a a ?+?。

现对N 与'N 的实际值取相反数，分别为：t a '-、t a ''-，则t a '-为反向加速度，t a ''

-

为正向加速度。将新得到的正向与反向加速度分别与原来的正、反向加速度相加后取平均值：

()()1222

t t t t a a a a a a '''+?--+?-=

()()()1

2

2

2

t t t

t

a a a a a a '''

----+?+-+?=

则得到新的正反向加速度分别为：122t a a a ?-?+

，21

2

t a a a ?-?-+。 由于

1221

12,22

a a a a a a ?-??-?

5.3问题三的求解

声屏障检测仪的工作原理是：通过内部的加速度传感器来记录车辆经过时声屏障振动而产生的加速度数值(密集采样)。将加速度数据通过数值积分，按照加速度-位移的物理公式将加速度数据转化为震动的位移，并通过震动位移对声屏障状态进行判断。

在现实生活中，油井示功图位移测量技术可以了解油井实时动态及抽油装置是否合理，其核心硬件为加速度信号调理电路。加速度传感器负责采集油杆上下运动的加速度信号，并对原始数据进行奇异值处理、相应滤波技术处理以最大限度的将噪声信号剔除，通过加速度双重积分得到示功图测试仪的位移和冲程。因此我们可以将改进后的加速度检测仪应用于此，从而消除加速度信号由于电源纹波和信号干扰的影响所引起波形的微小畸变。

改进后的加速度检测仪还可以应用于惯性导航系统。惯性导航系统是通过对加速度测量值进行积分而完成定位的。因此，虽然惯性传感器测量值中的误差、噪声可能非常小，系统在一小段时间内的相对定位精度也可能很高，但是这些测量误差、噪声在积分计算过程中会随着时间的推移而积累成越来越大的定位误差，再不采取任何措施的情况下定位误差最终可达无穷大，因此我们可以将改进后的加速度检测仪应用于此，以消除系统的测量误差和噪声。

该加速度检测仪还可以用于测量电梯的加速度。电梯上升和下降的过程加速度变化情况也与声屏障类似，服从一个正弦变换。速度的变化也是先快速增大，缓慢达到极值，速度下降由慢到快，直至速度降为0，且最终的位移也是一个稳定值。而对测量出的电梯的加速度积分求出其位移，可以用来判断电梯是否有故障。

在这三个场景中，均需通过加速度来求位移，而正反方向的加速度应该对称，因此用我们的处理数据的方法同样可以减小误差，最后使位移趋于一个稳定值而不是一个变化的值。

六、模型的评价与推广

6.1 模型的评价

对于问题一，在根据原加速度数据求解速度和位移的过程中我们采用了组合辛普森积分算法，考虑到加速度监测仪采样频率为1000HZ求速度和位移会产生很小的误差，可以忽略不计，该算法不仅保证了模型的精确度，而且避免了繁琐的计算过程。在误差分析的过程中，我们从定性分析和定量计算两个方面分析了该加速度检测仪的系统误差

和随机误差，并且运用了卡尔曼滤波算法和3次7点平滑滤波法尽可能地减小了随机误差和系统误差。经我们改进后的模型和数据处理方法具有普适性和一般性，在大样本情况下依然能够快速地对数据进行有效处理，可以推广应用到许多生产和生活中，具有很强的说服力和可靠性。

但是我们建立的模型也有一些缺点，比如说建模方法较单一，对于同一个问题没有建立多个模型进行比较，就第一问来说，我们可以利用精确度更高而计算量等同的龙贝格算法。再比如说减少了系统误差和随机误差之后所得到的一些数据还是会产生一定的不可避免的误差，而且我们也没有较全面地分析系统误差在有无脉动风时的差异性，这些都需要我们对模型做进一步的改进。

6.2 模型的推广

6.2.1推广领域

基于加速度传感器的油井示功图位移测量技术研究

6.2.2原理介绍

基于加速度传感器的示功图测试仪的硬件框图如图1所示，主要分为四部分：电源管理系统、载荷信号调理电路、加速度信号调理电路、主控及辅助设备电路。其中，加速度传感器负责采集油杆上下运动的加速度信号，通过积分算法而得到位移和冲程，是示功图测试仪硬件的设计核心。

示功图测试仪的位移和冲程是通过加速度双重积分得到的，考虑到油杆上下周期性运动的特点，将一个周期加速度测量值减去其平均值，令边界条件为零，对修正后的加速度积分得到速度；令边界条件为零，对修正后的速度积分得到冲程。示功仪加速度在两次去除边界积分后并不能得到准确的冲程，往往对于同一口井会得出两个差异很大的冲程。因为加速度量的电压信号很小，3V供电系统造成加速度与电压的比例系数很小(0.56V/g)，MCU采集加速度电压信号受干扰严重，所以必须对采集的加速度信号进行合适的滤波后再双重积分得到冲程。首先，将采集到一个周期的加速度的数据存放在RAM 中，对加速度数据进行奇异值的滤除；然后对加速度进行3次7点平滑窗滤波，最大限度地将噪声信号滤除；最后，应用周期去边界的双重积分得到各点的位移值，对应各点载荷量，在液晶屏上画出示功图，并将油井信息和示功图信息储存到外部EEPROM中。

七、参考文献

[1] [美]Charles K. Chui，[中] Guanrong Chen等, 《卡尔曼滤波及其实时应用》（第四版），清华大学出版社，2014年4月

[2] 施浒立,赵彦,《误差设计新理念与方法》，北京：科学出版社，2007年：4

[3] 赵黎，王丰，《大学物理实验》，上海：复旦大学出版社，2013年：7-

[4]陈为真，汪秉文，胡晓娅，《基于时域积分的加速度信号处理》，《华中科技大学学报(自然科学版》，第38卷，第1期，2010年1月

[5] 吴畅，《浅析系统误差消除法》，《计量技术》，2001年

[6] 司守奎等，《数学建模算法及应用》，国防工业出版社，2014 年 6 月

[7] 卓金武等，《MATLAB在数学建模中的应用》，北京航空航天大学出版社，2011年4月

八、附录

8.1 附录清单

附录1：求解问题一的LINGO程序

附录2：求解问题二的Mathematica程序

附录3：求解问题一的中间数据

附录4：问题二的完整数据结果

8.2 附录正文

附录1：求解问题一的matlab程序

附录2：求解问题二的matlab程序

format long g

[A_value]=textread('A.txt','%n',1395);

subplot(1,2,1),plot(A_value,'k+-');hold on

grid on

title('\fontsize{20}A组数据原始图');

set(gca,'xtick',[0:200:1500]); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%对A组数据进行二次卡尔曼滤波消除随机误差

N=1395;

A_filter_1=zeros(1,N);

A_filter_1(1)=A_value(1);%为A_filter_1(k)赋初值

p(1)=1;%A_filter_1(1)对应的协方差

Q=cov(randn(1,N));%过程噪声的协方差

R=cov(randn(1,N));%测量噪声的协方差

for k=2:N%循环里面是卡尔曼滤波的具体过程

A_filter_1(k)=A_filter_1(k-1);

p(k)=p(k-1)+Q;

Kg(k)=p(k)/(p(k)+R);%Kg为Kalman Gain,卡尔曼增益

A_filter_1(k)=A_filter_1(k)+Kg(k)*(A_value(k)-A_filter_1(k));

p(k)=(1-Kg(k))*p(k);

end

plot(A_filter_1,'+k');hold on

A_filter_2=zeros(1,N);

A_filter_2(1)=A_filter_1(1);%为x(k)赋初值

p(1)=1;%x(1)对应的协方差

Q=cov(randn(1,N));%过程噪声的协方差

R=cov(randn(1,N));%测量噪声的协方差

for k=2:N%循环里面是卡尔曼滤波的具体过程

A_filter_2(k)=A_filter_2(k-1);

p(k)=p(k-1)+Q;

Kg(k)=p(k)/(p(k)+R);%Kg为Kalman Gain,卡尔曼增益

A_filter_2(k)=A_filter_2(k)+Kg(k)*(A_filter_1(k)-A_filter_2(k)); p(k)=(1-Kg(k))*p(k);

end

A_std_2=std(A_filter_2)

subplot(1,2,2),plot(A_filter_2,'+k');hold on

title('\fontsize{20}二次滤波后');

grid on

axis([0 1400 -3 4]);

附录3：求解问题一的中间数据

附录4：问题二的完整数据结果