Ask Your Question

Revision history [back]

click to hide/show revision 1
initial version

Kalman filter code

Hi everyone, this is my first time on a forum so if I do something wrong please forgive me. I'm trying to build my own Kalman Filter, not because I don't like OpenCV's but because I think that when you want to really understand how something works you should get your hands dirty on it. So I made my researches and then modified a Kalman filter tutorial I found online, in order to see the different behaviour of my KF and OpenCV's one. The results were not really optimal: my KF seems to have a perfect and fast tracking but with way too much noise...OpenCV's one is a bit slower but much more stable...and for what I'm going to do with KF I prefer stability to velocity. So the question is: what am I doing wrong? Any help will be really appreciated. I don't think there's something wrong in the algorithm itself, I found plenty of papers online dealing with the algorithm...so I think the problem it somwhere in the coding part, but I really can't figure out where. Main things I added to the original: the GBKalman structure, the functions GBKalman_init, GBKalman_predict, GBKalman_correct, GBKalman_update. Here's the code:

#include "opencv2/highgui/highgui.hpp"
#include "opencv2/video/tracking.hpp"
#include <Windows.h>
#include <iostream>

#define drawCross( center, color, d )                                 \
line( img, Point( center.x - d, center.y - d ), Point( center.x + d, center.y + d ), color, 2, CV_AA, 0); \
line( img, Point( center.x + d, center.y - d ), Point( center.x - d, center.y + d ), color, 2, CV_AA, 0 )

using namespace cv;
using namespace std;

typedef struct GBKalman
{
    int n;
    int m;
    Mat A;
    Mat H;
    Mat Q;
    Mat R;
    Mat x;
    Mat x_priori;
    Mat x_previous;
    Mat P;
    Mat P_priori;
    Mat P_previous;
    Mat Kk;
    Mat Htransp;
    Mat Atransp;
    Mat I;
    Mat temp;
};
void GBKalman_init(GBKalman* G, int a, int c, Mat D, Mat F, Mat I, Mat J, Mat X, Mat Y)
{
    G->n = a;
    G->m = c;
    G->A = D;
    G->H = F;
    G->Q = I;
    G->R = J;
    G->x_previous = X;
    G->x_priori = Mat_<float>(a, 1);
    G->x = Mat_<float>(a, 1);

    G->P_previous = Y;
    G->P_priori = Mat_<float>(a, a);
    G->P = Mat_<float>(a, a);

    G->Kk = Mat_<float>(a, c);

    G->Htransp = Mat_<float>(a, c); G->Htransp = G->H.t();
    G->Atransp = Mat_<float>(a, a); G->Atransp = G->A.t();
    G->I.eye(a, a, CV_32F);
    G->temp = Mat(c, c, CV_32F);
}
void GBKalman_predict(GBKalman* G)
{
    G->x_priori = G->A * G->x_previous;
    G->P_priori = G->A * G->P_previous * G->Atransp + G->Q;
}
void GBKalman_correct(GBKalman* G, Mat zk)
{
    G->temp = G->R + G->H * G->P_priori * G->Htransp;
    G->Kk = G->P_priori * G->Htransp * G->temp.inv();
    G->x = G->x_priori + G->Kk * (zk - G->H * G->x_priori);
    G->P = (G->I - G->Kk * G->H) * G->P_priori;
}
void GBKalman_update(GBKalman *G)
{
    G->x_previous = G->x;
    G->P_previous = G->P;
}

int main()
{

    KalmanFilter KF(4, 2, 0);
    POINT mousePos;
    GetCursorPos(&mousePos);

    // intialization of KF...
    KF.transitionMatrix = *(Mat_<float>(4, 4) << 1, 0, 1, 0, 0, 1, 0, 1, 0, 0, 1, 0, 0, 0, 0, 1);
    Mat_<float> measurement(2, 1); measurement.setTo(Scalar(0));

    KF.statePre.at<float>(0) = mousePos.x;
    KF.statePre.at<float>(1) = mousePos.y;
    KF.statePre.at<float>(2) = 0;
    KF.statePre.at<float>(3) = 0;
    setIdentity(KF.measurementMatrix);
    setIdentity(KF.processNoiseCov, Scalar::all(1e-4));
    setIdentity(KF.measurementNoiseCov, Scalar::all(1));
    setIdentity(KF.errorCovPost, Scalar::all(.1));

    // initialization of my Kalman filter
    GBKalman GBK;
    Mat A = (Mat_<float>(4, 4) << 1, 0, 1, 0, 0, 1, 0, 1, 0, 0, 1, 0, 0, 0, 0, 1);
    Mat H = (Mat_<float>(2, 4) << 1, 0, 0, 0, 0, 1, 0, 0);
    Mat Q = (Mat_<float>(4, 4) << 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1);
    Mat R = (Mat_<float>(2, 2) << 3, 0, 0, 3);
    Mat x0 = (Mat_<float>(4, 1) << mousePos.x, mousePos.y, 0, 0);
    Mat P0 = (Mat_<float>(4, 4) << 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1);
    GBKalman_init(&GBK, 4, 2, A, H, Q, R, x0, P0);

    // Image to show mouse tracking
    Mat img(600, 800, CV_8UC3);
    vector<Point> mousev, kalmanv;
    mousev.clear();
    kalmanv.clear();

    Mat_<float> v(2, 1); setIdentity(v, Scalar::all(1e-1));

    while (1)
    {
        // First predict, to update the internal statePre variable
        Mat prediction = KF.predict();
        Point predictPt(prediction.at<float>(0), prediction.at<float>(1));

        //gaussian measurement noise
        randn(v, 0, 5);

        // Get mouse point without noise. Noise will be added when calling the update functions
        GetCursorPos(&mousePos);
        measurement(0) = mousePos.x;
        measurement(1) = mousePos.y;

        // my kalman routine
        GBKalman_predict(&GBK);
        GBKalman_correct(&GBK, measurement+v);
        GBKalman_update(&GBK);

        // The update phase
        Mat estimated = KF.correct(measurement+v);

        Point statePt(estimated.at<float>(0), estimated.at<float>(1));  // opencv estimated point
        Point measPt(measurement(0), measurement(1));                   // true point without noise

        Point truepoint(measurement(0)+v(0), measurement(1)+v(1));      // true mesured point (with noise)
        Point myPt(GBK.x.at<float>(0),GBK.x.at<float>(1));              // my KF estimated point

        // plot points
        imshow("mouse kalman", img);
        img = Scalar::all(0);

        drawCross(myPt, Scalar(255, 0, 0), 5); // blue : estimated position by my KF
        drawCross(statePt, Scalar(255, 255, 255), 5); // white : estimated by OpenCV's KF
        //drawCross(measPt, Scalar(0, 0, 255), 5); // red: real position without noise
        //drawCross(truepoint, Scalar(0, 255, 0), 5); // green: real position + noise

        waitKey(10);
    }
    return 0;
}

thank everybody in advance

Kalman filter code

Hi everyone, this is my first time on a forum so if I do something wrong please forgive me. I'm trying to build my own Kalman Filter, not because I don't like OpenCV's but because I think that when you want to really understand how something works you should get your hands dirty on it. So I made my researches and then modified a Kalman filter tutorial I found online, in order to see the different behaviour of my KF and OpenCV's one. The results were not really optimal: my KF seems to have a perfect and fast tracking but with way too much noise...OpenCV's one is a bit slower but much more stable...and for what I'm going to do with KF I prefer stability to velocity. So the question is: what am I doing wrong? Any help will be really appreciated. I don't think there's something wrong in the algorithm itself, I found plenty of papers online dealing with the algorithm...so I think the problem it somwhere in the coding part, but I really can't figure out where. Main things I added to the original: the GBKalman structure, the functions GBKalman_init, GBKalman_predict, GBKalman_correct, GBKalman_update. Here's the code:

#include "opencv2/highgui/highgui.hpp"
#include "opencv2/video/tracking.hpp"
#include <Windows.h>
#include <iostream>

#define drawCross( center, color, d )                                 \
line( img, Point( center.x - d, center.y - d ), Point( center.x + d, center.y + d ), color, 2, CV_AA, 0); \
line( img, Point( center.x + d, center.y - d ), Point( center.x - d, center.y + d ), color, 2, CV_AA, 0 )

using namespace cv;
using namespace std;

typedef struct GBKalman
{
    int n;
    int m;
    Mat A;
    Mat H;
    Mat Q;
    Mat R;
    Mat x;
    Mat x_priori;
    Mat x_previous;
    Mat P;
    Mat P_priori;
    Mat P_previous;
    Mat Kk;
    Mat Htransp;
    Mat Atransp;
    Mat I;
    Mat temp;
};
void GBKalman_init(GBKalman* G, int a, int c, Mat D, Mat F, Mat I, Mat J, Mat X, Mat Y)
{
    G->n = a;
    G->m = c;
    G->A = D;
    G->H = F;
    G->Q = I;
    G->R = J;
    G->x_previous = X;
    G->x_priori = Mat_<float>(a, 1);
    G->x = Mat_<float>(a, 1);

    G->P_previous = Y;
    G->P_priori = Mat_<float>(a, a);
    G->P = Mat_<float>(a, a);

    G->Kk = Mat_<float>(a, c);

    G->Htransp = Mat_<float>(a, c); G->Htransp = G->H.t();
    G->Atransp = Mat_<float>(a, a); G->Atransp = G->A.t();
    G->I.eye(a, a, CV_32F);
    G->temp = Mat(c, c, CV_32F);
}
void GBKalman_predict(GBKalman* G)
{
    G->x_priori = G->A * G->x_previous;
    G->P_priori = G->A * G->P_previous * G->Atransp + G->Q;
}
void GBKalman_correct(GBKalman* G, Mat zk)
{
    G->temp = G->R + G->H * G->P_priori * G->Htransp;
    G->Kk = G->P_priori * G->Htransp * G->temp.inv();
    G->x = G->x_priori + G->Kk * (zk - G->H * G->x_priori);
    G->P = (G->I - G->Kk * G->H) * G->P_priori;
}
void GBKalman_update(GBKalman *G)
{
    G->x_previous = G->x;
    G->P_previous = G->P;
}

int main()
{

    KalmanFilter KF(4, 2, 0);
    POINT mousePos;
    GetCursorPos(&mousePos);

    // intialization of KF...
    KF.transitionMatrix = *(Mat_<float>(4, 4) << 1, 0, 1, 0, 0, 1, 0, 1, 0, 0, 1, 0, 0, 0, 0, 1);
    Mat_<float> measurement(2, 1); measurement.setTo(Scalar(0));

    KF.statePre.at<float>(0) = mousePos.x;
    KF.statePre.at<float>(1) = mousePos.y;
    KF.statePre.at<float>(2) = 0;
    KF.statePre.at<float>(3) = 0;
    setIdentity(KF.measurementMatrix);
    setIdentity(KF.processNoiseCov, Scalar::all(1e-4));
    setIdentity(KF.measurementNoiseCov, Scalar::all(1));
    setIdentity(KF.errorCovPost, Scalar::all(.1));

    // initialization of my Kalman filter
    GBKalman GBK;
    Mat A = (Mat_<float>(4, 4) << 1, 0, 1, 0, 0, 1, 0, 1, 0, 0, 1, 0, 0, 0, 0, 1);
    Mat H = (Mat_<float>(2, 4) << 1, 0, 0, 0, 0, 1, 0, 0);
    Mat Q = (Mat_<float>(4, 4) << 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1);
    Mat R = (Mat_<float>(2, 2) << 3, 0, 0, 3);
    Mat x0 = (Mat_<float>(4, 1) << mousePos.x, mousePos.y, 0, 0);
    Mat P0 = (Mat_<float>(4, 4) << 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1);
    GBKalman_init(&GBK, 4, 2, A, H, Q, R, x0, P0);

    // Image to show mouse tracking
    Mat img(600, 800, CV_8UC3);
    vector<Point> mousev, kalmanv;
    mousev.clear();
    kalmanv.clear();

    Mat_<float> v(2, 1); setIdentity(v, Scalar::all(1e-1));

    while (1)
    {
        // First predict, to update the internal statePre variable
        Mat prediction = KF.predict();
        Point predictPt(prediction.at<float>(0), prediction.at<float>(1));

        //gaussian measurement noise
        randn(v, 0, 5);

        // Get mouse point without noise. Noise will be added when calling the update functions
        GetCursorPos(&mousePos);
        measurement(0) = mousePos.x;
        measurement(1) = mousePos.y;

        // my kalman routine
        GBKalman_predict(&GBK);
        GBKalman_correct(&GBK, measurement+v);
        GBKalman_update(&GBK);

        // The update phase
        Mat estimated = KF.correct(measurement+v);

        Point statePt(estimated.at<float>(0), estimated.at<float>(1));  // opencv estimated point
        Point measPt(measurement(0), measurement(1));                   // true point without noise

        Point truepoint(measurement(0)+v(0), measurement(1)+v(1));      // true mesured point (with noise)
        Point myPt(GBK.x.at<float>(0),GBK.x.at<float>(1));              // my KF estimated point

        // plot points
        imshow("mouse kalman", img);
        img = Scalar::all(0);

        drawCross(myPt, Scalar(255, 0, 0), 5); // blue : estimated position by my KF
        drawCross(statePt, Scalar(255, 255, 255), 5); // white : estimated by OpenCV's KF
        //drawCross(measPt, Scalar(0, 0, 255), 5); // red: real position without noise
        //drawCross(truepoint, Scalar(0, 255, 0), 5); // green: real position + noise

        waitKey(10);
    }
    return 0;
}

thank everybody in advance