Unstability/misunderstanding of savgolfilt

[Dalecanka]

Hi.

I downloaded your code and did some test runs. However, I encountered some issues.

I need to use the code within the bigger project. I have a .dll library where I need to apply s-g filtering. Although the main file is .cu, I don't use any Cuda code now. This library is then called from MFC application to run some operations on data.
I am working with Visual Studio 2019, Eigen 3.3, on Windows 10.

Issue 1

This code is pretty much copied from example.cpp. I run the function "savgolfilt" with a given set of values 30 times and sometimes, I get different result from "Filter". I use F, x, x values from example.

My code:

// initialize values
    int F = 9;      
    int k = 5; 
    float x_min = 100, x_max = 1000;

// create vector of x values
    auto x = VectorXf::LinSpaced(F, x_min, x_max);

// calculate the same Filter 30 times
    for (int i = 0; i < 30; i++)
    {
        auto Filter = savgolfilt(x, k, F);
        log_file << "Filter =  " << Filter << "\n"; // save to a log_file, which is later saved as .txt
    }

Results are as follows: middle value is different sometimes.

Filter = 662.501 775 887.5 1000 550 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 550 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 1.96313e+06 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 957289 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 957289 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 550 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 1.96313e+06 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 1.96313e+06 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 550 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 1.96313e+06 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 957289 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 957289 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 1.96313e+06 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 957289 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 957289 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 957289 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 550 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 1.96313e+06 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 1.96313e+06 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 957289 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 1.96313e+06 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 550 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 1.96313e+06 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 550 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 550 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 550 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 550 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 1.96313e+06 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 1.96313e+06 99.9998 212.5 325 437.5 Filter = 662.501 775 887.5 1000 957289 99.9998 212.5 325 437.5

I run a couple of tests to see where the error happens. It seems that the line y << y_off.transpose().eval(), y_ss, y_on.transpose().eval(); in savgolfilt function is problematic - the vectors before concatenating into y seem just fine. Do you know what could be a problem?

Issue 2

I apply savgolfilt on my data. Typically, F = 129, k = 6, x is signal of length 400 (in other applications the length can be up to 1600) with float values in range for example (-0.0038, 0.0266). Here, I show example for F = 9 and k = 5 (in order not to show so many values).

        // initialize values
        int F = 9;
        int k = 5;

        // create vector of x values
        auto x = my_data;

        // calculate
        for (int i = 0; i < 10; i++)
        {
            auto Filter = savgolfilt(x, k, F);
            log_file << "Filter =  " << Filter << "\n";
        }

Result

Filter = 0.0066643 0.00634263 0.00642954 0.00900544 -nan 0.00869701 0.00574144 0.00707652 0.00823701 Filter = 0.0066643 0.00634263 0.00642954 0.00900544 -nan 0.00869701 0.00574144 0.00707652 0.00823701 Filter = 0.0066643 0.00634263 0.00642954 0.00900544 -nan 0.00869701 0.00574144 0.00707652 0.00823701 Filter = 0.0066643 0.00634263 0.00642954 0.00900544 -nan 0.00869701 0.00574144 0.00707652 0.00823701 Filter = 0.0066643 0.00634263 0.00642954 0.00900544 -nan 0.00869701 0.00574144 0.00707652 0.00823701 Filter = 0.0066643 0.00634263 0.00642954 0.00900544 -nan 0.00869701 0.00574144 0.00707652 0.00823701 Filter = 0.0066643 0.00634263 0.00642954 0.00900544 -nan 0.00869701 0.00574144 0.00707652 0.00823701 Filter = 0.0066643 0.00634263 0.00642954 0.00900544 -nan 0.00869701 0.00574144 0.00707652 0.00823701 Filter = 0.0066643 0.00634263 0.00642954 0.00900544 -nan 0.00869701 0.00574144 0.00707652 0.00823701 Filter = 0.0066643 0.00634263 0.00642954 0.00900544 -nan 0.00869701 0.00574144 0.00707652 0.00823701

There is nan in the middle. Is there any condition for data? Like they should be positive, in some range, or? Or, is there just problem with << operator, which displays values of Eigen variables?

Issue/question 3

I am sorry, but I got confused, even more after reading other issues here. How exactly should I use the code to get smoothed signal? I expect the smoothed signal will have the same length as input signal x, but it is not the case of savgolfilt. In Readme.md there is "savgolfilt(x, x_on, k, F);", but I don't see the four-argument input of this function nor in savgol.cpp, nor in its header.

Thank you for your answers. I am sorry if I ask something trivial, but I didn't find answer anywhere.

[journeytosilius]

I'm also having issues ... trying to feed savgolfilt with a vector with data and I get an exception:

error: invalid initialization of reference of type const VectorXf&’ {aka ‘const Eigen::Matrix<float, -1, 1>&’}

How are you supposed to load your own data to the filter ? Do I need to convert a vector of floats to An Eigen matrix ?

Thanks

[Dalecanka]

Hi. The definition of savgolfilt is VectorXd savgolfilt(VectorXd const & x, int k, int F), so try to feed it with vector of doubles (VectorXd) instead of floats. I used it this way.

[journeytosilius]

@Dalecanka the definition is a VectorXf : RowVectorXf savgolfilt(VectorXf const & x, int k, int F)

Can you share your code please ?

[Dalecanka]

Aah, I am sorry. It has been a few months and I forgot what changes did I make. I think I compared that with another SG filter in Matlab, perhaps combined, and made calculation in double. I send you the code; however I copied it in a .txt file since this system doesn't take .cpp and .h as input. Let me know if it helps.

savgol_cpp.txt savgol_header.txt