6

我为平滑曲线编写了这段代码。它在一个点旁边取 5 个点并将它们相加并取平均值。

/* Smoothing */
void smoothing(vector<Point2D> &a)
{
    //How many neighbours to smooth
    int NO_OF_NEIGHBOURS=10;
    vector<Point2D> tmp=a;
    for(int i=0;i<a.size();i++)
    {

        if(i+NO_OF_NEIGHBOURS+1<a.size())
        {
            for(int j=1;j<NO_OF_NEIGHBOURS;j++)
            {
                a.at(i).x+=a.at(i+j).x;
                a.at(i).y+=a.at(i+j).y;
            }
            a.at(i).x/=NO_OF_NEIGHBOURS;
            a.at(i).y/=NO_OF_NEIGHBOURS;

        }
        else
        {
            for(int j=1;j<NO_OF_NEIGHBOURS;j++)
            {
                a.at(i).x+=tmp.at(i-j).x;
                a.at(i).y+=tmp.at(i-j).y;
            }
            a.at(i).x/=NO_OF_NEIGHBOURS;
            a.at(i).y/=NO_OF_NEIGHBOURS;
        }

    }

}

但是我为每个点得到了非常高的值,而不是与前一点相似的值。形状最大化了很多,这个算法出了什么问题?

4

6 回答 6

10

你在这里看到的是一个有限脉冲响应 (FIR) 滤波器的低音实现,它实现了一个boxcar 窗口函数。考虑 DSP 方面的问题,您需要vector使用NO_OF_NEIGHBOURS相等的 FIR 系数过滤输入,每个系数的值都为1/NO_OF_NEIGHBOURS. 通常最好使用已建立的算法,而不是重新发明轮子。

这是一个相当邋遢的实现,我快速敲定了过滤器加倍。您可以轻松地对其进行修改以过滤您的数据类型。该演示显示了对上升锯函数 (0,.25,.5,1) 的几个周期的过滤,仅用于演示目的。它可以编译,因此您可以使用它。

#include <iostream>
#include <vector>

using namespace std;

class boxFIR
{
    int numCoeffs; //MUST be > 0
    vector<double> b; //Filter coefficients
    vector<double> m; //Filter memories

public:
    boxFIR(int _numCoeffs) :
    numCoeffs(_numCoeffs)
    {
        if (numCoeffs<1)
            numCoeffs = 1; //Must be > 0 or bad stuff happens

        double val = 1./numCoeffs;
        for (int ii=0; ii<numCoeffs; ++ii) {
            b.push_back(val);
            m.push_back(0.);
        }
    }    

    void filter(vector<double> &a)
    {
        double output;

        for (int nn=0; nn<a.size(); ++nn)
        {
            //Apply smoothing filter to signal
            output = 0;
            m[0] = a[nn];
            for (int ii=0; ii<numCoeffs; ++ii) {
                output+=b[ii]*m[ii];
            }

            //Reshuffle memories
            for (int ii = numCoeffs-1; ii!=0; --ii) {
                m[ii] = m[ii-1];
            }                        
            a[nn] = output;
        }
    }


};

int main(int argc, const char * argv[])
{
    boxFIR box(1); //If this is 1, then no filtering happens, use bigger ints for more smoothing

    //Make a rising saw function for demo
    vector<double> a;
    a.push_back(0.); a.push_back(0.25); a.push_back(0.5); a.push_back(0.75); a.push_back(1.);
    a.push_back(0.); a.push_back(0.25); a.push_back(0.5); a.push_back(0.75); a.push_back(1.);
    a.push_back(0.); a.push_back(0.25); a.push_back(0.5); a.push_back(0.75); a.push_back(1.);
    a.push_back(0.); a.push_back(0.25); a.push_back(0.5); a.push_back(0.75); a.push_back(1.);

    box.filter(a);

    for (int nn=0; nn<a.size(); ++nn)
    {
        cout << a[nn] << endl;
    }
}

使用这条线增加滤波器系数的数量,以查看逐渐平滑的输出。只有 1 个滤波器系数,没有平滑。

boxFIR box(1);

该代码足够灵活,您甚至可以根据需要更改窗口形状。通过修改构造函数中定义的系数来做到这一点。

注意:这将为您的实现提供稍微不同的输出,因为这是一个因果过滤器(仅取决于当前样本和以前的样本)。您的实现不是因果关系,因为它会在未来的样本中提前进行平均,这就是为什么您需要条件语句来处理接近向量末尾的情况。如果您想要输出类似于您尝试使用此算法对您的过滤器执行的操作,请通过此算法反向运行您的向量(只要窗口函数是对称的,这工作正常)。这样,您可以获得类似的输出,而无需算法的讨厌的条件部分。

于 2012-10-19T12:01:13.383 回答
3

过滤有利于“记忆”平滑。这是learnvst答案的反向传递,以防止相位失真

for (int i = a.size(); i > 0; --i)
{
    // Apply smoothing filter to signal
    output = 0;
    m[m.size() - 1] = a[i - 1];

    for (int j = numCoeffs; j > 0; --j) 
        output += b[j - 1] * m[j - 1];

    // Reshuffle memories
    for (int j = 0; j != numCoeffs; ++j) 
        m[j] = m[j + 1];

    a[i - 1] = output;
}

更多关于 MATLAB 中的零相位失真 FIR 滤波器:http: //www.mathworks.com/help/signal/ref/filtfilt.html

于 2015-04-07T14:20:59.870 回答
3

在以下块中:

            for(int j=0;j<NO_OF_NEIGHBOURS;j++)
            {
                a.at(i).x=a.at(i).x+a.at(i+j).x;
                a.at(i).y=a.at(i).y+a.at(i+j).y;
            }

对于每个邻居,您分别将 a.at(i)'sx 和 y 添加到邻居值。

我理解正确,应该是这样的。

            for(int j=0;j<NO_OF_NEIGHBOURS;j++)
            {
                a.at(i).x += a.at(i+j+1).x
                a.at(i).y += a.at(i+j+1).y
            }
于 2012-10-19T11:34:19.780 回答
1

You make addition with point itself when you need to take neighbor points - just offset index by 1:

for(int j=0;j<NO_OF_NEIGHBOURS;j++)
 {
    a.at(i).x += a.at(i+j+1).x
    a.at(i).y += a.at(i+j+1).y
 }
于 2012-10-19T12:05:06.953 回答
1

该点的当前值被使用了两次:一次是因为你使用+=,一次是如果y==0。因此,您正在构建例如 6 个点的总和,但仅除以 5。这个问题在 IF 和 ELSE 情况下都存在。另外:您应该检查向量是否足够长,否则您的 ELSE 案例将在负索引处读取。

以下本身不是问题,而只是一个想法:您是否考虑过使用只触及每个点两次的算法?:您可以存储一个临时 xy 值(初始化为与第一个点相同),然后当您访问每个点您只需添加新点并减去最旧的点,如果它比您的NEIGHBOURS背部更远。您为每个点更新这个“运行总和”,并将这个值除以NEIGHBOURS-number 存储到新点中。

于 2012-10-19T11:42:22.467 回答
0

这对我来说很好:

for (i = 0; i < lenInput; i++)
{
    float x = 0;
    for (int j = -neighbours; j <= neighbours; j++)
    {
        x += input[(i + j <= 0) || (i + j >= lenInput) ? i : i + j];
    }
    output[i] = x / (neighbours * 2 + 1);
}
于 2019-06-09T14:24:32.983 回答