Ask Your Question

Revision history [back]

click to hide/show revision 1
initial version

Video Stabilization using goodFeaturesToTrack and calcOpticalFlowPyrLK

Hi all,

I'm currently trying to achieve video stabilization using goodFeaturesToTrack and calcOpticalFlowPyrLK.

I get the good features and calculate the good optical flow. I then get the image trajectory and smooth out the image trajectory and apply it to the video frames.

This all seems to be working correct, but the problem I'm facing is that the processed frames contain large black portions. How can I get the frames centered in the video without seeing the large black areas? Any help or suggestions will really be appreciated.

Below is the code I'm currently using:

#include <opencv2/opencv.hpp>
#include <iostream>
#include <cassert>
#include <cmath>
#include <fstream>
#include <ctime>

using namespace std;
using namespace cv;

const int SMOOTHING_RADIUS = 30; // In frames. The larger the more stable the video, but less reactive to sudden panning
const int HORIZONTAL_BORDER_CROP = 60; // In pixels. Crops the border to reduce the black borders from stabilisation being too noticeable.

struct TransformParam {
    TransformParam() {}
    TransformParam(double _dx, double _dy, double _da)
    {
        dx = _dx;
        dy = _dy;
        da = _da;
    }

    double dx;
    double dy;
    double da; // angle
};

struct Trajectory {
    Trajectory() {}
    Trajectory(double _x, double _y, double _a)
    {
        x = _x;
        y = _y;
        a = _a;
    }

    double x;
    double y;
    double a; // angle
};

void stabilize(const vector<Mat> &images, vector<Mat> &resultImages) {
  vector<TransformParam> prev_to_cur_transform;
  Mat last_T;

  for (int i = 1; i < images.size(); i++) {
    Mat prev = images[i-1];
    Mat cur = images[i];

    // vector from prev to cur
    vector<Point2f> prevCorner, curCorner;
    vector<Point2f> prevCorner2, curCorner2;
    vector<uchar> status;
    vector<float> error;

    Mat curGrey;
    Mat prevGrey;

    cvtColor(cur, curGrey, COLOR_BGR2GRAY);
    cvtColor(prev, prevGrey, COLOR_BGR2GRAY);

    size_t totalEntries = curGrey.cols * curGrey.rows;

    goodFeaturesToTrack(prevGrey, prevCorner, totalEntries, 0.01, 1);
    calcOpticalFlowPyrLK(prevGrey, curGrey, prevCorner, curCorner, status, error);

    // weed out bad matches
    for (int i = 0; i < status.size(); i++) {
      if (status[i]) {
        prevCorner2.push_back(prevCorner[i]);
        curCorner2.push_back(curCorner[i]);
      }
    }

    // translation + rotation only
    // false = rigidTransform, no scaling/shearing
    Mat T = estimateRigidTransform(prevCorner2, curCorner2, false);

    // in rare cases no valid transform is found, use the last known good transform
    if (T.data == NULL) {
      last_T.copyTo(T);
    }

    T.copyTo(last_T);

    // decompose T
    double dx = T.at<double>(0, 2);
    double dy = T.at<double>(1, 2);
    double da = atan2(T.at<double>(1, 0), T.at<double>(0, 0));

    prev_to_cur_transform.push_back(TransformParam(dx, dy, da));
    cout << "i: " << i << " " << "dx: " << dx << " " << "dy: " << dy << " " << "da: " << da << endl;
    cout << endl;

    curGrey.copyTo(prevGrey);
    cout << "Image: " << i << "/" << images.size() << " - good optical flow: " << prevCorner2.size() << endl;
  }

  // Step 2 - Accumulate transformations to get image trajectory
  // Accumulated frame to frame transformation
  double a = 0;
  double x = 0;
  double y = 0;

  vector<Trajectory> trajectory; // trajectory at all frames

  for (int i = 0; i < prev_to_cur_transform.size(); i++) {
    x += prev_to_cur_transform[i].dx;
    y += prev_to_cur_transform[i].dy;
    a += prev_to_cur_transform[i].da;

    trajectory.push_back(Trajectory(x, y, a));
    cout << "Trajectory " << (i+1) << ":" << " " << "x: " << x << " " << "y: " << y << " " << "a: " << a << endl;
  }

  // Step 3 - Smooth out trajectory using an averaging window
  vector<Trajectory> smoothed_trajectory; // trajectory at all frames

  for (int i = 0; i < trajectory.size(); i++) {
    double sum_x = 0;
    double sum_y = 0;
    double sum_a = 0;
    int count = 0;

    for (int j = -SMOOTHING_RADIUS; j <= SMOOTHING_RADIUS; j++) {
      if (i+j >= 0 && i+j < trajectory.size()) {
        sum_x += trajectory[i+j].x;
        sum_y += trajectory[i+j].y;
        sum_a += trajectory[i+j].a;

        count++;
      }
    }

    double avg_a = sum_a / count;
    double avg_x = sum_x / count;
    double avg_y = sum_y / count;

    smoothed_trajectory.push_back(Trajectory(avg_x, avg_y, avg_a));
    cout << "Smoothed Trajectory " << (i+1) << ":" << " " << "avg_x: " << avg_x << " " << "avg_y: " << avg_y << " " << "avg_a: " << avg_a << endl;
  }

  // Step 4 - Generate new set of previous to current transform, such that trajectory ends up being the same as smoothed_trajectory
  vector<TransformParam> new_prev_to_cur_transform;

  // Accumulated frame to frame transformation
  a = 0;
  x = 0;
  y = 0;

  for (int i = 0; i < prev_to_cur_transform.size(); i++) {
    x += prev_to_cur_transform[i].dx;
    y += prev_to_cur_transform[i].dy;
    a += prev_to_cur_transform[i].da;

    // target - current
    double diff_x = smoothed_trajectory[i].x - x;
    double diff_y = smoothed_trajectory[i].y - y;
    double diff_a = smoothed_trajectory[i].a - a;

    double dx = prev_to_cur_transform[i].dx + diff_x;
    double dy = prev_to_cur_transform[i].dy + diff_y;
    double da = prev_to_cur_transform[i].da + diff_a;

    new_prev_to_cur_transform.push_back(TransformParam(dx, dy, da));
    cout << "New Transform " << (i+1) << ":" << " " << "dx: " << dx << " " << "dy: " << dy << " " << "da: " << da << endl;
  }

  // Step 5 - Apply new transformation
  Mat T(2, 3, CV_64F);

  // get correct aspect ratio
  int vert_border = HORIZONTAL_BORDER_CROP * images[0].rows / images[0].cols;

  for (int i = 0; i < images.size(); i++) {
    Mat cur = images[i];
    T.at<double>(0, 0) = cos(new_prev_to_cur_transform[i].da);
    T.at<double>(0, 1) = -sin(new_prev_to_cur_transform[i].da);
    T.at<double>(1, 0) = sin(new_prev_to_cur_transform[i].da);
    T.at<double>(1, 1) = cos(new_prev_to_cur_transform[i].da);

    T.at<double>(0, 2) = new_prev_to_cur_transform[i].dx;
    T.at<double>(1, 2) = new_prev_to_cur_transform[i].dy;

    Mat cur2;
    warpAffine(cur, cur2, T, cur.size());

    cur2 = cur2(Range(vert_border, cur2.rows - vert_border), Range(HORIZONTAL_BORDER_CROP, cur2.cols - HORIZONTAL_BORDER_CROP));
    resize(cur2, cur2, cur.size());

    resultImages.push_back(cur2);
  }
  cout << "resultImages size: " << resultImages.size() << endl;
}

int main() {
  clock_t tStart = clock();

  string input = "vietnam.avi";
  string output = "vietnam.mov";

  vector<Mat> images;
  vector<Mat> resultImages;

  int k = 180;

  for (int i = 0; i < k; i++) {
    Mat image;
    image = imread(cv::format("/Users/martinj/Downloads/opencv-master/samples/cpp/images_b/%d.jpg", i), IMREAD_COLOR);

    if (image.empty()) {
      cout << "Image was empty" << endl;
      return 1;
    }
    images.push_back(image);
    cout << "Reading image " << i << "/" << k << endl;
  }
  cout << "images size: " << images.size() << endl;

  stabilize(images, resultImages);

  VideoWriter videoWriter(output, 0x6765706a, 25, cv::Size(images[0].size().width, images[0].size().height), true);

  for (int i = 0; i < resultImages.size(); i++) {
    // imwrite(cv::format("/Users/martinj/Downloads/opencv-master/samples/cpp/images_processed/%d.jpg", i), resultImages[i]);
    videoWriter.write(resultImages[i]);
    cout << "Writing image " << i << "/" << resultImages.size() << endl;
  }

  printf("Time taken: %.2fs\n", (double)(clock() - tStart) / CLOCKS_PER_SEC);
}