最近在写一些数据处理的程序。经常需要对数据进行平滑处理。直接用FIR滤波器或IIR滤波器都有一个启动问题,滤波完成后总要对数据掐头去尾。因此去找了些简单的数据平滑处理的方法。

在一本老版本的《数学手册》中找到了几个基于最小二乘法的数据平滑算法。将其写成了C 代码,测试了一下,效果还可以。这里简单的记录一下,算是给自己做个笔记。

算法的原理很简单,以五点三次平滑为例。取相邻的5个数据点,可以拟合出一条3次曲线来,然后用3次曲线上相应的位置的数据值作为滤波后结果。简单的说就是 Savitzky-Golay 滤波器 。只不过Savitzky-Golay 滤波器并不特殊考虑边界的几个数据点,而这个算法还特意把边上的几个点的数据拟合结果给推导了出来。

不多说了,下面贴代码。首先是线性拟合平滑处理的代码. 分别为三点线性平滑、五点线性平滑和七点线性平滑。

void linearSmooth3 ( double in[], double out[], int N )
{
int i;
if ( N < 3 )
{
for ( i = 0; i <= N - 1; i++ )
{
out[i] = in[i];
}
}
else
{
out[0] = ( 5.0 * in[0] + 2.0 * in[1] - in[2] ) / 6.0;
for ( i = 1; i <= N - 2; i++ )
{
out[i] = ( in[i - 1] + in[i] + in[i + 1] ) / 3.0;
}
out[N - 1] = ( 5.0 * in[N - 1] + 2.0 * in[N - 2] - in[N - 3] ) / 6.0;
}
}
void linearSmooth5 ( double in[], double out[], int N )
{
int i;
if ( N < 5 )
{
for ( i = 0; i <= N - 1; i++ )
{
out[i] = in[i];
}
}
else
{
out[0] = ( 3.0 * in[0] + 2.0 * in[1] + in[2] - in[4] ) / 5.0;
out[1] = ( 4.0 * in[0] + 3.0 * in[1] + 2 * in[2] + in[3] ) / 10.0;
for ( i = 2; i <= N - 3; i++ )
{
out[i] = ( in[i - 2] + in[i - 1] + in[i] + in[i + 1] + in[i + 2] ) / 5.0;
}
out[N - 2] = ( 4.0 * in[N - 1] + 3.0 * in[N - 2] + 2 * in[N - 3] + in[N - 4] ) / 10.0;
out[N - 1] = ( 3.0 * in[N - 1] + 2.0 * in[N - 2] + in[N - 3] - in[N - 5] ) / 5.0;
}
}
void linearSmooth7 ( double in[], double out[], int N )
{
int i;
if ( N < 7 )
{
for ( i = 0; i <= N - 1; i++ )
{
out[i] = in[i];
}
}
else
{
out[0] = ( 13.0 * in[0] + 10.0 * in[1] + 7.0 * in[2] + 4.0 * in[3] +
in[4] - 2.0 * in[5] - 5.0 * in[6] ) / 28.0;
out[1] = ( 5.0 * in[0] + 4.0 * in[1] + 3 * in[2] + 2 * in[3] +
in[4] - in[6] ) / 14.0;
out[2] = ( 7.0 * in[0] + 6.0 * in [1] + 5.0 * in[2] + 4.0 * in[3] +
3.0 * in[4] + 2.0 * in[5] + in[6] ) / 28.0;
for ( i = 3; i <= N - 4; i++ )
{
out[i] = ( in[i - 3] + in[i - 2] + in[i - 1] + in[i] + in[i + 1] + in[i + 2] + in[i + 3] ) / 7.0;
}
out[N - 3] = ( 7.0 * in[N - 1] + 6.0 * in [N - 2] + 5.0 * in[N - 3] +
4.0 * in[N - 4] + 3.0 * in[N - 5] + 2.0 * in[N - 6] + in[N - 7] ) / 28.0;
out[N - 2] = ( 5.0 * in[N - 1] + 4.0 * in[N - 2] + 3.0 * in[N - 3] +
2.0 * in[N - 4] + in[N - 5] - in[N - 7] ) / 14.0;
out[N - 1] = ( 13.0 * in[N - 1] + 10.0 * in[N - 2] + 7.0 * in[N - 3] +
4 * in[N - 4] + in[N - 5] - 2 * in[N - 6] - 5 * in[N - 7] ) / 28.0;
}
}

然后是利用二次函数拟合平滑。

void quadraticSmooth5(double in[], double out[], int N)
{
int i;
if ( N < 5 )
{
for ( i = 0; i <= N - 1; i++ )
{
out[i] = in[i];
}
}
else
{
out[0] = ( 31.0 * in[0] + 9.0 * in[1] - 3.0 * in[2] - 5.0 * in[3] + 3.0 * in[4] ) / 35.0;
out[1] = ( 9.0 * in[0] + 13.0 * in[1] + 12 * in[2] + 6.0 * in[3] - 5.0 *in[4]) / 35.0;
for ( i = 2; i <= N - 3; i++ )
{
out[i] = ( - 3.0 * (in[i - 2] + in[i + 2]) +
12.0 * (in[i - 1] + in[i + 1]) + 17 * in[i] ) / 35.0;
}
out[N - 2] = ( 9.0 * in[N - 1] + 13.0 * in[N - 2] + 12.0 * in[N - 3] + 6.0 * in[N - 4] - 5.0 * in[N - 5] ) / 35.0;
out[N - 1] = ( 31.0 * in[N - 1] + 9.0 * in[N - 2] - 3.0 * in[N - 3] - 5.0 * in[N - 4] + 3.0 * in[N - 5]) / 35.0;
}
}
void quadraticSmooth7(double in[], double out[], int N)
{
int i;
if ( N < 7 )
{
for ( i = 0; i <= N - 1; i++ )
{
out[i] = in[i];
}
}
else
{
out[0] = ( 32.0 * in[0] + 15.0 * in[1] + 3.0 * in[2] - 4.0 * in[3] -
6.0 * in[4] - 3.0 * in[5] + 5.0 * in[6] ) / 42.0;
out[1] = ( 5.0 * in[0] + 4.0 * in[1] + 3.0 * in[2] + 2.0 * in[3] +
in[4] - in[6] ) / 14.0;
out[2] = ( 1.0 * in[0] + 3.0 * in [1] + 4.0 * in[2] + 4.0 * in[3] +
3.0 * in[4] + 1.0 * in[5] - 2.0 * in[6] ) / 14.0;
for ( i = 3; i <= N - 4; i++ )
{
out[i] = ( -2.0 * (in[i - 3] + in[i + 3]) +
3.0 * (in[i - 2] + in[i + 2]) +
6.0 * (in[i - 1] + in[i + 1]) + 7.0 * in[i] ) / 21.0;
}
out[N - 3] = ( 1.0 * in[N - 1] + 3.0 * in [N - 2] + 4.0 * in[N - 3] +
4.0 * in[N - 4] + 3.0 * in[N - 5] + 1.0 * in[N - 6] - 2.0 * in[N - 7] ) / 14.0;
out[N - 2] = ( 5.0 * in[N - 1] + 4.0 * in[N - 2] + 3.0 * in[N - 3] +
2.0 * in[N - 4] + in[N - 5] - in[N - 7] ) / 14.0;
out[N - 1] = ( 32.0 * in[N - 1] + 15.0 * in[N - 2] + 3.0 * in[N - 3] -
4.0 * in[N - 4] - 6.0 * in[N - 5] - 3.0 * in[N - 6] + 5.0 * in[N - 7] ) / 42.0;
}
}

最后是三次函数拟合平滑。

/**
* 五点三次平滑
*
*/
void cubicSmooth5 ( double in[], double out[], int N )
{
int i;
if ( N < 5 )
{
for ( i = 0; i <= N - 1; i++ )
out[i] = in[i];
}
else
{
out[0] = (69.0 * in[0] + 4.0 * in[1] - 6.0 * in[2] + 4.0 * in[3] - in[4]) / 70.0;
out[1] = (2.0 * in[0] + 27.0 * in[1] + 12.0 * in[2] - 8.0 * in[3] + 2.0 * in[4]) / 35.0;
for ( i = 2; i <= N - 3; i++ )
{
out[i] = (-3.0 * (in[i - 2] + in[i + 2])+ 12.0 * (in[i - 1] + in[i + 1]) + 17.0 * in[i] ) / 35.0;
}
out[N - 2] = (2.0 * in[N - 5] - 8.0 * in[N - 4] + 12.0 * in[N - 3] + 27.0 * in[N - 2] + 2.0 * in[N - 1]) / 35.0;
out[N - 1] = (- in[N - 5] + 4.0 * in[N - 4] - 6.0 * in[N - 3] + 4.0 * in[N - 2] + 69.0 * in[N - 1]) / 70.0;
}
return;
}
void cubicSmooth7(double in[], double out[], int N)
{
int i;
if ( N < 7 )
{
for ( i = 0; i <= N - 1; i++ )
{
out[i] = in[i];
}
}
else
{
out[0] = ( 39.0 * in[0] + 8.0 * in[1] - 4.0 * in[2] - 4.0 * in[3] +
1.0 * in[4] + 4.0 * in[5] - 2.0 * in[6] ) / 42.0;
out[1] = ( 8.0 * in[0] + 19.0 * in[1] + 16.0 * in[2] + 6.0 * in[3] -
4.0 * in[4] - 7.0* in[5] + 4.0 * in[6] ) / 42.0;
out[2] = ( -4.0 * in[0] + 16.0 * in [1] + 19.0 * in[2] + 12.0 * in[3] +
2.0 * in[4] - 4.0 * in[5] + 1.0 * in[6] ) / 42.0;
for ( i = 3; i <= N - 4; i++ )
{
out[i] = ( -2.0 * (in[i - 3] + in[i + 3]) +
3.0 * (in[i - 2] + in[i + 2]) +
6.0 * (in[i - 1] + in[i + 1]) + 7.0 * in[i] ) / 21.0;
}
out[N - 3] = ( -4.0 * in[N - 1] + 16.0 * in [N - 2] + 19.0 * in[N - 3] +
12.0 * in[N - 4] + 2.0 * in[N - 5] - 4.0 * in[N - 6] + 1.0 * in[N - 7] ) / 42.0;
out[N - 2] = ( 8.0 * in[N - 1] + 19.0 * in[N - 2] + 16.0 * in[N - 3] +
6.0 * in[N - 4] - 4.0 * in[N - 5] - 7.0 * in[N - 6] + 4.0 * in[N - 7] ) / 42.0;
out[N - 1] = ( 39.0 * in[N - 1] + 8.0 * in[N - 2] - 4.0 * in[N - 3] -
4.0 * in[N - 4] + 1.0 * in[N - 5] + 4.0 * in[N - 6] - 2.0 * in[N - 7] ) / 42.0;
}
}

上面的代码经过了简单的测试,可以放心使用。

最近还有人向我要9点线性平滑和11点线性平滑的系数。我简单算了算,把系数列在了后面。不过我觉得这两组系数意义不大,与其这样做还不如设计一个FIR 滤波器或者 SG 滤波器。

9点线性平滑

yy[0]=(17*y[0] + 14*y[1] + 11*y[2] + 8*y[3] + 5*y[4] + 2*y[5] - y[6] - 4*y[7] - 7*y[8])/45;
yy[1]=(56*y[0] + 47*y[1] + 38*y[2] + 29*y[3] + 20*y[4] + 11*y[5] + 2*y[6] - 7*y[7] - 16*y[8])/180;
yy[2]=(22*y[0] + 19*y[1] + 16*y[2] + 13*y[3] + 10*y[4] + 7*y[5] + 4*y[6] + y[7] - 2*y[8])/90;
yy[3]=(32*y[0] + 29*y[1] + 26*y[2] + 23*y[3] + 20*y[4] + 17*y[5] + 14*y[6] + 11*y[7] + 8*y[8])/ 180;

yy[4]=(y[0] + y[1] +  y[2] + y[3] + y[4] + y[5] + y[6] + y[7] + y[8])/9;

yy[5]=(8*y[0] + 11*y[1] + 14*y[2] + 17*y[3] + 20*y[4] + 23*y[5] + 26*y[6] + 29*y[7] + 32*y[8])/ 180;
yy[6]=(-2*y[0] + y[1] + 4*y[2] + 7*y[3] + 10*y[4] + 13*y[5] + 16*y[6] + 19*y[7] + 22*y[8])/90;
yy[7]=(-16*y[0] - 7*y[1] +  2*y[2] + 11*y[3] + 20*y[4] + 29*y[5] + 38*y[6] + 47*y[7] + 56*y[8])/ 180;
yy[8]=(-7*y[0] - 4*y[1] - y[2] + 2*y[3] + 5*y[4] + 8*y[5] + 11*y[6] + 14*y[7] + 17*y[8])/45;

11点线性平滑

yy[0]=0.3181818181818182*
      y[0] + 0.2727272727272727*y[1] + 0.22727272727272727*y[2] +
    0.18181818181818182*y[3] + 0.13636363636363635*y[4] +
    0.09090909090909091*y[5] + 0.045454545454545456*y[6] +
    1.6152909428804953*^-17*y[7] - 0.045454545454545456*y[8] -
    0.09090909090909091*y[9] - 0.13636363636363635*y[10]

yy[1]=0.27272727272727276*
      y[0] + 0.23636363636363633*y[1] + 0.19999999999999998*y[2] +
    0.16363636363636364*y[3] + 0.12727272727272726*y[4] +
    0.09090909090909091*y[5] + 0.05454545454545455*y[6] +
    0.018181818181818202*y[7] - 0.01818181818181818*y[8] -
    0.054545454545454536*y[9] - 0.09090909090909088*y[10]

yy[2]=0.2272727272727273*
      y[0] + 0.19999999999999996*y[1] + 0.17272727272727267*y[2] +
    0.14545454545454542*y[3] + 0.11818181818181815*y[4] +
    0.09090909090909091*y[5] + 0.06363636363636363*y[6] +
    0.03636363636363638*y[7] + 0.009090909090909094*y[8] -
    0.01818181818181816*y[9] - 0.0454545454545454*y[10]

yy[3]=0.18181818181818188*
      y[0] + 0.16363636363636355*y[1] + 0.1454545454545454*y[2] +
    0.12727272727272723*y[3] + 0.10909090909090904*y[4] +
    0.0909090909090909*y[5] + 0.07272727272727272*y[6] +
    0.054545454545454564*y[7] + 0.036363636363636376*y[8] +
    0.01818181818181823*y[9] + 8.326672684688674*^-17*y[10]

yy[4]=0.13636363636363644*
      y[0] + 0.12727272727272718*y[1] + 0.11818181818181811*y[2] +
    0.10909090909090903*y[3] + 0.09999999999999995*y[4] +
    0.0909090909090909*y[5] + 0.08181818181818182*y[6] +
    0.07272727272727275*y[7] + 0.06363636363636364*y[8] +
    0.05454545454545459*y[9] + 0.04545454545454555*y[10]

yy[5]=0.090909090909091*
      y[0] + 0.09090909090909077*y[1] + 0.09090909090909083*y[2] +
    0.09090909090909083*y[3] + 0.09090909090909084*y[4] +
    0.0909090909090909*y[5] + 0.09090909090909091*y[6] +
    0.09090909090909094*y[7] + 0.09090909090909093*y[8] +
    0.090909090909091*y[9] + 0.09090909090909102*y[10]

yy[6]=0.04545454545454558*
      y[0] + 0.0545454545454544*y[1] + 0.06363636363636352*y[2] +
    0.07272727272727264*y[3] + 0.08181818181818173*y[4] +
    0.0909090909090909*y[5] + 0.1*y[6] + 0.10909090909090911*y[7] +
    0.11818181818181821*y[8] + 0.12727272727272737*y[9] +
    0.13636363636363652*y[10]

yy[7]=1.1102230246251565*^-\
16*y[0] + 0.01818181818181802*y[1] + 0.03636363636363624*y[2] +
    0.05454545454545445*y[3] + 0.07272727272727264*y[4] +
    0.0909090909090909*y[5] + 0.10909090909090909*y[6] +
    0.12727272727272732*y[7] + 0.1454545454545455*y[8] +
    0.16363636363636372*y[9] + 0.181818181818182*y[10]

yy[8]=-0.0454545454545453*
      y[0] - 0.018181818181818354*y[1] + 0.009090909090908955*y[2] +
    0.03636363636363624*y[3] + 0.06363636363636355*y[4] +
    0.09090909090909088*y[5] + 0.11818181818181818*y[6] +
    0.1454545454545455*y[7] + 0.17272727272727273*y[8] +
    0.2000000000000001*y[9] + 0.22727272727272746*y[10]

yy[9]=-0.09090909090909072*
      y[0] - 0.05454545454545473*y[1] - 0.018181818181818327*y[2] +
    0.01818181818181805*y[3] + 0.05454545454545444*y[4] +
    0.09090909090909088*y[5] + 0.12727272727272726*y[6] +
    0.16363636363636366*y[7] + 0.2*y[8] + 0.23636363636363647*y[9] +
    0.27272727272727293*y[10]

yy[10]=-0.1363636363636362*
      y[0] - 0.09090909090909116*y[1] - 0.04545454545454561*y[2] -
    1.6653345369377348*^-16*y[3] + 0.04545454545454533*y[4] +
    0.09090909090909088*y[5] + 0.13636363636363635*y[6] +
    0.18181818181818188*y[7] + 0.2272727272727273*y[8] +
    0.27272727272727293*y[9] + 0.3181818181818184*y[10]

几个简单的数据点平滑处理算法相关推荐

  1. 等宽分箱_数据分析师-数据挖掘如何分箱以及对箱子中的数据进行平滑处理

    题干: 假定用于分析的数据包含属性age.数据元组中age的值如下(按递增序):13,15,16,16,19,20,20,21,22,22, 25,25,25,30,33,33,35,35,36,40 ...

  2. 用python做一个数据查询软件_Python实现功能简单的数据查询及可视化系统

    欢迎点击右上角关注小编,除了分享技术文章之外还有很多福利,私信学习资料可以领取包括不限于Python实战演练.PDF电子文档.面试集锦.学习资料等. image.png 前言 数据时代,数据的多源集成 ...

  3. 【项目实战课】人人免费可学!基于Pytorch的图像分类简单任务数据增强实战

    欢迎大家来到我们的项目实战课,本期内容是<基于Pytorch的图像分类简单任务数据增强实战>.所谓项目实战课,就是以简单的原理回顾+详细的项目实战的模式,针对具体的某一个主题,进行代码级的 ...

  4. [云炬python3玩转机器学习笔记] 3-12 数据加载和简单的数据探索

    读取数据和简单的数据探索 In [1]: import numpy as np In [2]: import matplotlib as mpl import matplotlib.pyplot as ...

  5. 你真的会玩SQL吗?简单的数据修改

    你真的会玩SQL吗?系列目录 你真的会玩SQL吗?之逻辑查询处理阶段 你真的会玩SQL吗?和平大使 内连接.外连接 你真的会玩SQL吗?三范式.数据完整性 你真的会玩SQL吗?查询指定节点及其所有父节 ...

  6. 【MATLAB统计分析与应用100例】案例004:matlab调用smoothts函数对金融时间序列数据进行平滑处理

    时间序列数据进行平滑处理效果预览: Excel原数据部分预览: 1. 读取金融时间序列数据,绘制波形图 x = xlsread('examp03_02.xls'); % 从文件examp03_02.x ...

  7. Python之数据分析(案例:净额成交量、矢量化、数据曲线平滑)

    文章目录 一.净额成交量 二.标量函数矢量化 三.数据曲线平滑 一.净额成交量 1.案例分析 符号数组a: [10 -20 30 40 -50] 用法:numpy.sign(a) 结果是[1 -1 1 ...

  8. AJAX for Java简单表数据查询实例

    AJAX for Java简单表数据查询实例<?XML:NAMESPACE PREFIX = O /> AJAX WebShop 3对Java开发具有良好的支持,同时也提供了各种层次的后台 ...

  9. python实现excel计算_用python实现简单EXCEL数据统计的实例

    任务: 用python时间简单的统计任务-统计男性和女性分别有多少人. 用到的物料:xlrd 它的作用-读取excel表数据 代码:import xlrd workbook = xlrd.open_w ...

最新文章

  1. ARM研发进展与企业清单
  2. 【Android 安全】DEX 加密 ( Application 替换 | 创建用户自定义 Application | 替换 ContextImpl 对象的 mOuterContext 成员 )
  3. “约见”面试官系列之常见面试题之第七十篇之==和===(建议收藏)
  4. 音频开发基本知识总结
  5. 如何实现分布式锁?赶紧收藏备战金三银四!
  6. 以后再也不用看“教程”!概括【配置环境】的原理,小白也能举一反三的python配置环境过程!
  7. AtCoder Beginner Contest 171 B - Mix Juice
  8. Python爬取图片时,urllib提示没有属性urlretrieve的问题
  9. 改善程序与设计的55个具体做法 day5
  10. day04-selenium滚动和常见反爬
  11. python学习笔记——类与对象、常用函数
  12. mysql压缩修复数据库_压缩修复Access数据库
  13. 2019年面试大概问题
  14. Rockchip DRM RK628 Porting Guide
  15. pd调节规律_PD 控制规律的特点是具有超前控制功能。( )
  16. 32微型计算机中的32室,32 位微型计算机中的32是指(   )_学小易找答案
  17. 基于python的个人博客系统的设计开题报告_基于JavaSSM框架的个人博客系统设计与实现开题报告...
  18. spring @Autowired用法
  19. jdbc连接teradata仓库_teradata的jdbc连接,
  20. 我是如何通过系统架构设计师考试的(2017年软考)

热门文章

  1. office 2010很强大很好用
  2. 2222222222222
  3. 因影片里白人太多,Twitch平台撤下BLM视频
  4. 深度学习论文精读05——用于乳腺癌组织病理学图像分类的深度特征
  5. servlet和jsp简介
  6. 亲测有效,解决MAC无法联网的问题
  7. (业务)自动生成邀请码
  8. mysql安装流程以及各类问题解决
  9. Python自动化工具(自动化操作)
  10. 用星坐标(Star Coordinates)表示高维数据