--- /dev/null
+#include <stdio.h>
+#include <string.h>
+#include <stdlib.h>
+#include <unistd.h>
+#include <math.h>
+#include "../src/img.h"
+#include <opencv2/core/core_c.h>
+#include <opencv2/imgproc/imgproc_c.h>
+#include <opencv2/video/tracking.hpp>
+#include<opencv2/highgui/highgui.hpp>
+
+/* test single flow without smoothing */
+//#define TEST_SINGLE
+
+/* test input image */
+//#define TEST_INPUT
+
+int flow_window = 30, border = 0;
+double limit_frames = 0;
+
+static inline void smooth_pixle(double *array, int width, int height, int i, int j)
+{
+ double val = 0;
+ int num = 0;
+ /* left-top */
+ if (i > 0 && j > 0) {
+ val += array[i*width+j-1-width];
+ num++;
+ }
+ /* top */
+ if (i > 0) {
+ val += array[i*width+j-width];
+ num++;
+ }
+ /* right-top */
+ if (i > 0 && j < width-1) {
+ val += array[i*width+j+1-width];
+ num++;
+ }
+ /* left */
+ if ( j > 0) {
+ val += array[i*width+j-1];
+ num++;
+ }
+ /* right */
+ if ( j < width-1) {
+ val += array[i*width+j+1];
+ num++;
+ }
+ /* left-bottom */
+ if (i < height-1 && j > 0) {
+ val += array[i*width+j-1+width];
+ num++;
+ }
+ /* bottom */
+ if (i < height-1) {
+ val += array[i*width+j+width];
+ num++;
+ }
+ /* right-bottom */
+ if (i < height-1 && j < width-1) {
+ val += array[i*width+j+1+width];
+ num++;
+ }
+ if (num)
+ array[i*width+j] = val/num;
+}
+
+/* interpolate missing pixles in array */
+static void smooth(double *array, double *mask, int width, int height)
+{
+ int i, j, iter;
+
+ for (iter = 0; iter < 3; iter++) {
+ for (i = 0; i < height; i++) {
+ for (j = 0; j < width; j++) {
+ if (*mask++ == 0.0)
+ continue;
+ smooth_pixle(array, width, height, i, j);
+ }
+ }
+ for (i = height-1; i >= 0; i--) {
+ for (j = width-1; j >= 0; j--) {
+ if (*--mask == 0.0)
+ continue;
+ smooth_pixle(array, width, height, i, j);
+ }
+ }
+ }
+}
+
+/* use optical flow to interpolate two images with given offset */
+static int interpolate(unsigned short *img1, unsigned short *img2, unsigned short *res, int width, int height, double offset, int frame)
+{
+ static int lastframe = -1;
+ static CvMat *flow1 = NULL, *flow2 = NULL;
+ static IplImage *image1, *image2;
+ double *array1, *array2;
+ CvSize size;
+ int i, j, x, y;
+ unsigned short *i1, *i2;
+ unsigned int grey;
+ double alpha;
+ const float *f;
+ const float pyrScale = 0.5;
+ const float levels = 3;
+ const float winsize = flow_window;
+ const float iterations = 8;
+ const float polyN = 5;
+ const float polySigma = 1.2;
+ const int flags = 0;
+
+ if (lastframe != frame && flow1) {
+ cvReleaseMat(&flow1);
+ flow1 = NULL;
+ }
+ if (lastframe != frame && flow2) {
+ cvReleaseMat(&flow2);
+ flow2 = NULL;
+ }
+ lastframe = frame;
+
+ array1 = calloc(width*height*5*sizeof(*array1), 1);
+ array2 = calloc(width*height*5*sizeof(*array2), 1);
+ if (!array1 || !array2) {
+ printf("Failed to create image array\n");
+ return -1;
+ }
+
+ /* if we need to make flow maps, we do it */
+ if (!flow1 || !flow2) {
+ printf(" -> calculating flow %d<->%d\n", frame, frame+1);
+ flow1 = cvCreateMat(height, width, CV_32FC2);
+ flow2 = cvCreateMat(height, width, CV_32FC2);
+ if (!flow1 || !flow2) {
+ printf("Failed to create flow map\n");
+ return -1;
+ }
+
+ size.width = width;
+ size.height = height;
+// image1 = cvCreateImage(size, IPL_DEPTH_32F, 1);
+// image2 = cvCreateImage(size, IPL_DEPTH_32F, 1);
+ image1 = cvCreateImage(size, IPL_DEPTH_16U, 1);
+ image2 = cvCreateImage(size, IPL_DEPTH_16U, 1);
+ if (!image1 || !image2) {
+ printf("Failed to create image for optical flow\n");
+ return -1;
+ }
+ i1 = img1;
+ i2 = img2;
+ for (i = 0; i < height; i++) {
+ for (j = 0; j < width; j++) {
+ grey = *i1++;
+ grey += *i1++;
+ grey += *i1++;
+// CV_IMAGE_ELEM(image1, float, i, j) = grey / 76;
+ CV_IMAGE_ELEM(image1, ushort, i, j) = grey / 3;
+ grey = *i2++;
+ grey += *i2++;
+ grey += *i2++;
+// CV_IMAGE_ELEM(image2, float, i, j) = grey / 76;
+ CV_IMAGE_ELEM(image2, ushort, i, j) = grey / 3;
+ }
+ }
+#if 0
+cvNamedWindow( "Image1", 0 );// Create a window for display.
+cvShowImage( "Image1", image1); // Show our image inside it.
+cvNamedWindow( "Image2", 0 );// Create a window for display.
+cvShowImage( "Image2", image2); // Show our image inside it.
+#endif
+
+
+ cvCalcOpticalFlowFarneback(image1, image2, flow1, pyrScale, levels, winsize, iterations, polyN, polySigma, flags);
+ cvCalcOpticalFlowFarneback(image2, image1, flow2, pyrScale, levels, winsize, iterations, polyN, polySigma, flags);
+
+#ifndef TEST_INPUT
+ cvReleaseImage(&image1);
+ cvReleaseImage(&image2);
+#endif
+ }
+
+ /* apply image to array using flow map */
+ for (i = 0; i < height; i++) {
+ f = (const float *)(flow1->data.ptr + flow1->step * i);
+ for (j = 0; j < width; j++) {
+ x = j + round(*f++ * offset);
+ y = i + round(*f++ * offset);
+ if (j < border || i < border || j >= width-border || i >= height-border)
+ continue;
+ if (x < 0 || x >= width || y < 0 || y >= height)
+ continue;
+ array1[y*width+x] += (double)img1[(i*width+j)*3+0] / 65535;
+ array1[y*width+x + width*height] += (double)img1[(i*width+j)*3+1] / 65535;
+ array1[y*width+x + width*height*2] += (double)img1[(i*width+j)*3+2] / 65535;
+ array1[y*width+x + width*height*3] += 1.0;
+ }
+ }
+ for (i = 0; i < height; i++) {
+ f = (const float *)(flow2->data.ptr + flow2->step * i);
+ for (j = 0; j < width; j++) {
+ x = j + round(*f++ * (1-offset));
+ y = i + round(*f++ * (1-offset));
+ if (j < border || i < border || j >= width-border || i >= height-border)
+ continue;
+ if (x < 0 || x >= width || y < 0 || y >= height)
+ continue;
+ array2[y*width+x] += (double)img2[(i*width+j)*3+0] / 65535;
+ array2[y*width+x + width*height] += (double)img2[(i*width+j)*3+1] / 65535;
+ array2[y*width+x + width*height*2] += (double)img2[(i*width+j)*3+2] / 65535;
+ array2[y*width+x + width*height*3] += 1.0;
+ }
+ }
+
+ /* mask unset pixles */
+ for (i = 0; i < height; i++) {
+ for (j = 0; j < width; j++) {
+ alpha = array1[i*width+j + width*height*3];
+ if (alpha == 0.0) {
+ array1[i*width+j] = 0.0;
+ array1[i*width+j + width*height] = 0.0;
+ array1[i*width+j + width*height*2] = 0.0;
+ array1[i*width+j + width*height*3] = 0.0;
+ /* mask */
+ array1[i*width+j + width*height*4] = 1.0;
+ }
+ alpha = array2[i*width+j + width*height*3];
+ if (alpha == 0.0) {
+ array2[i*width+j] = 0.0;
+ array2[i*width+j + width*height] = 0.0;
+ array2[i*width+j + width*height*2] = 0.0;
+ array2[i*width+j + width*height*3] = 0.0;
+ /* mask */
+ array2[i*width+j + width*height*4] = 1.0;
+ }
+ }
+ }
+
+#ifndef TEST_SINGLE
+ /* smooth image by interpolating missing (masked) pixles */
+ smooth(array1, array1+width*height*4, width, height);
+ smooth(array1+width*height, array1+width*height*4, width, height);
+ smooth(array1+width*height*2, array1+width*height*4, width, height);
+ smooth(array1+width*height*3, array1+width*height*4, width, height);
+ smooth(array2, array2+width*height*4, width, height);
+ smooth(array2+width*height, array2+width*height*4, width, height);
+ smooth(array2+width*height*2, array2+width*height*4, width, height);
+ smooth(array2+width*height*3, array2+width*height*4, width, height);
+#endif
+
+ /* normalize smoothed image and compose */
+ for (i = 0; i < height; i++) {
+ for (j = 0; j < width; j++) {
+ alpha = array1[i*width+j + width*height*3];
+ if (alpha > 0.0) {
+ array1[i*width+j] /= alpha;
+ array1[i*width+j + width*height] /= alpha;
+ array1[i*width+j + width*height*2] /= alpha;
+ array1[i*width+j + width*height*3] = 1.0;
+ }
+ alpha = array2[i*width+j + width*height*3];
+ if (alpha > 0.0) {
+ array2[i*width+j] /= alpha;
+ array2[i*width+j + width*height] /= alpha;
+ array2[i*width+j + width*height*2] /= alpha;
+ array2[i*width+j + width*height*3] = 1.0;
+ }
+#ifndef TEST_SINGLE
+ array1[i*width+j] = (1.0-offset)*array1[i*width+j] + offset*array2[i*width+j];
+ array1[i*width+j + width*height] = (1.0-offset)*array1[i*width+j + width*height] + offset*array2[i*width+j + width*height];
+ array1[i*width+j + width*height*2] = (1.0-offset)*array1[i*width+j + width*height*2] + offset*array2[i*width+j + width*height*2];
+#endif
+#ifdef TEST_INPUT
+ if (i < (height>>1)) {
+ array1[i*width+j] = (double)CV_IMAGE_ELEM(image2, ushort, i, j) / 65535;
+ array1[i*width+j + width*height] = (double)CV_IMAGE_ELEM(image1, ushort, i, j) / 65535;
+ array1[i*width+j + width*height*2] = 0;
+ }
+#endif
+ }
+ }
+
+ array2img_short(array1, width, height, res, width, height, 0);
+
+ free(array1);
+ free(array2);
+
+ return 0;
+}
+
+int main(int argc, char *argv[])
+{
+ unsigned short *img1 = NULL, *img2 = NULL, *res = NULL;
+ int width, height, width2, height2;
+ const char *inputname, *outputname;
+ double oldfps, newfps, step, inputframe, baseframe, offset;
+ int outputframe;
+ int _argc = argc - 1;
+ int _arg = 1;
+ int rc;
+
+ if (_argc > 1 && !strcmp(argv[_arg], "--depth")) {
+ save_depth = atoi(argv[_arg+1]);
+ _arg += 2;
+ _argc -= 2;
+ }
+
+ if (_argc > 1 && !strcmp(argv[_arg], "--border")) {
+ border = atoi(argv[_arg+1]);
+ _arg += 2;
+ _argc -= 2;
+ }
+
+ if (_argc > 1 && !strcmp(argv[_arg], "--window")) {
+ flow_window = atoi(argv[_arg+1]);
+ _arg += 2;
+ _argc -= 2;
+ }
+
+ if (_argc > 1 && !strcmp(argv[_arg], "--frames")) {
+ limit_frames = atoi(argv[_arg+1]);
+ _arg += 2;
+ _argc -= 2;
+ }
+
+ if (_argc != 6) {
+ printf("This tool will change frame rate using optical flow and interpolation.\n");
+ printf("Usage: %s [--depth x] [--border x] [--window x] <input> <first> <output> <first> <oldfps> <newfps>\n", argv[0]);
+ printf(" --border <pixles> Remove border pixles, if darker than image (try 1)\n");
+ printf(" --depth <bits> Change color depth from default %d to given bits\n", save_depth);
+ printf(" --window <size> Change optical flow window size (default %d)\n", flow_window);
+ printf(" --frames <number> Limit number of input frames (default infinite)\n");
+ printf(" input Input image sequence (e.g. input_%%05d.ppm)\n");
+ printf(" output Output image sequence (e.g. output_%%05d.ppm)\n");
+ printf(" first Index of first input and output image\n");
+ printf(" oldfps Old frame rate\n");
+ printf(" newfps New frame rate\n");
+ return 0;
+ }
+
+ inputname = argv[_arg];
+ inputframe = atoi(argv[_arg+1]);
+ outputname = argv[_arg+2];
+ outputframe = atoi(argv[_arg+3]);
+ oldfps = strtod(argv[_arg+4], NULL);
+ newfps = strtod(argv[_arg+5], NULL);
+
+ step = oldfps / newfps;
+ printf("Time Scale Factor: %.4f\n", step);
+ printf("Optical Flow window: %d\n", flow_window);
+ if (limit_frames)
+ printf("Limit number of input frames: %.0f\n", limit_frames);
+
+next_frame:
+ baseframe = floor(inputframe);
+ offset = inputframe - baseframe;
+ if (offset < 0.001)
+ offset = 0.0;
+ else if (offset > 0.999) {
+ offset = 0.0;
+ baseframe++;
+ }
+ if (offset == 0.0) {
+ printf("Copy input frame %d, output frame %d\n", (int)baseframe, outputframe);
+ res = load_img(&width, &height, inputname, baseframe);
+ if (!res)
+ goto end;
+ rc = save_img(res, width, height, 0, outputname, outputframe);
+ if (rc < 0)
+ goto out;
+ } else {
+ printf("Processing input frame %.4f, output frame %d\n", inputframe, outputframe);
+ img1 = load_img(&width, &height, inputname, baseframe);
+ if (!img1)
+ goto end;
+ img2 = load_img(&width2, &height2, inputname, baseframe+1);
+ if (!img2)
+ goto end;
+ if (width != width2) {
+ printf("Frame %d and %d have different width (%d != %d)\n", (int)baseframe, (int)baseframe+1, width, width2);
+ goto out;
+ }
+ if (height != height2) {
+ printf("Frame %d and %d have different height (%d != %d)\n", (int)baseframe, (int)baseframe+1, height, height2);
+ goto out;
+ }
+ res = malloc(width * height * 3 * sizeof(unsigned short));
+ if (!res)
+ goto out;
+ rc = interpolate(img1, img2, res, width, height, offset, (int)baseframe);
+ if (rc < 0)
+ goto out;
+ rc = save_img(res, width, height, 0, outputname, outputframe);
+ if (rc < 0)
+ goto out;
+ }
+
+ free(img1);
+ img1 = NULL;
+ free(img2);
+ img2 = NULL;
+ free(res);
+ res = NULL;
+
+ outputframe++;
+ inputframe += step;
+ if (limit_frames) {
+ limit_frames -= step;
+ if (limit_frames <= 0)
+ goto out;
+ }
+ goto next_frame;
+
+end:
+ printf("We end here, since there is no more image to process\n");
+
+out:
+ free(img1);
+ free(img2);
+ free(res);
+ return 0;
+}