//=============================================================================
#include "movieIO.h"
#include "rgbImage.h"
#include "featureCorrelator.h"
#include "subWindow.h"
#include "percentBar.h"
#include "math.h"

extern "C" void svdcmp( float **a, int m, int n, float *w, float **v );
extern "C" void svbksb( float **u, float *w, float **v, int m, int n,
			float *b, float *x );
extern "C" void gaussj( float **a, int n, float **b, int m );

// C = A.B
#define MATMULT_3X3( A, B, C ) \
  (C)[0] = (A)[0]*(B)[0] + (A)[1]*(B)[3] + (A)[2]*(B)[6]; \
  (C)[1] = (A)[0]*(B)[1] + (A)[1]*(B)[4] + (A)[2]*(B)[7]; \
  (C)[2] = (A)[0]*(B)[2] + (A)[1]*(B)[5] + (A)[2]*(B)[8]; \
  (C)[3] = (A)[3]*(B)[0] + (A)[4]*(B)[3] + (A)[5]*(B)[6]; \
  (C)[4] = (A)[3]*(B)[1] + (A)[4]*(B)[4] + (A)[5]*(B)[7]; \
  (C)[5] = (A)[3]*(B)[2] + (A)[4]*(B)[5] + (A)[5]*(B)[8]; \
  (C)[6] = (A)[6]*(B)[0] + (A)[7]*(B)[3] + (A)[8]*(B)[6]; \
  (C)[7] = (A)[6]*(B)[1] + (A)[7]*(B)[4] + (A)[8]*(B)[7]; \
  (C)[8] = (A)[6]*(B)[2] + (A)[7]*(B)[5] + (A)[8]*(B)[8]
  
#define MATMULT_3X3b( A, B, c1, c2, c3 ) \
  (c1)[0] = (A)[0]*(B)[0] + (A)[1]*(B)[3] + (A)[2]*(B)[6]; \
  (c1)[1] = (A)[0]*(B)[1] + (A)[1]*(B)[4] + (A)[2]*(B)[7]; \
  (c1)[2] = (A)[0]*(B)[2] + (A)[1]*(B)[5] + (A)[2]*(B)[8]; \
  (c2)[0] = (A)[3]*(B)[0] + (A)[4]*(B)[3] + (A)[5]*(B)[6]; \
  (c2)[1] = (A)[3]*(B)[1] + (A)[4]*(B)[4] + (A)[5]*(B)[7]; \
  (c2)[2] = (A)[3]*(B)[2] + (A)[4]*(B)[5] + (A)[5]*(B)[8]; \
  (c3)[0] = (A)[6]*(B)[0] + (A)[7]*(B)[3] + (A)[8]*(B)[6]; \
  (c3)[1] = (A)[6]*(B)[1] + (A)[7]*(B)[4] + (A)[8]*(B)[7]; \
  (c3)[2] = (A)[6]*(B)[2] + (A)[7]*(B)[5] + (A)[8]*(B)[8]

// Determinant of a 3X3 matrix
#define DET_3X3(a1,a2,a3,b1,b2,b3,c1,c2,c3) \
  (a1)*((b2)*(c3)-(b3)*(c2)) \
 -(b1)*((a2)*(c3)-(a3)*(c2)) \
 +(c1)*((a2)*(b3)-(a3)*(b2))
    
#define DET_3X3(A,B,C) \
  (A)[0] * ((B)[1]*(C)[2]-(B)[2]*(C)[1]) \
 -(B)[0] * ((A)[1]*(C)[2]-(A)[2]*(C)[1]) \
 +(C)[0] * ((A)[1]*(B)[2]-(A)[2]*(B)[1])

// Determinant of a 4X4 matrix
#define DET_4X4(A,B,C,D) \
  (A)[0] * DET_3X3( &(B)[1], &(C)[1], &(D)[1] ) \
 -(B)[0] * DET_3X3( &(A)[1], &(C)[1], &(D)[1] ) \
 +(C)[0] * DET_3X3( &(A)[1], &(B)[1], &(D)[1] ) \
 -(D)[0] * DET_3X3( &(A)[1], &(B)[1], &(C)[1] )

#define DET_4X4(A) DET_4X4( &(A)[0], &(A)[4], &(A)[8], &(A)[12] )

#define 
void computeM( const float *F, const float *e1, const float *e2,
	       float *P1, float *P2 )
{
  int i;
	  
  // P1 is set to be I
  for ( i=1; i<12; i++ )
    P1[i]=0.0f;
  P1[0] = 1.0f;
  P1[4] = 1.0f;
  P1[8] = 1.0f;

  // P2 = [ [e2]XF | e2 ]
  float e2XMat[3][3];
  e2XMat[0][0]= 0.0f;   e2XMat[0][1]=-e2[2];  e2XMat[0][2]= e2[1];
  e2XMat[1][0]= e2[2];  e2XMat[1][1]= 0.0f;   e2XMat[1][2]=-e2[0];
  e2XMat[2][0]=-e2[1];  e2XMat[2][1]= e2[0];  e2XMat[2][2]= 0.0f;

  MATMULT_3X3b( e2XMat[0], F, &P2[0], &P2[4], &P2[8] );
  P2[ 3] = e2[0];
  P2[ 7] = e2[1];
  P2[11] = e2[2];
}
  
float findProjBasisT( float *pCorresp, float *T )
{
  int i, j, k;
  float aBuf[3][3];
  float *a[4];
  for ( i=0; i<3; i++ )
    a[i+1] = aBuf[i]-1;

  for ( i=0; i<3; i++ )
  {
    aBuf[0][i] = pCorresp[i*4];
    aBuf[1][i] = pCorresp[i*4+1];
    aBuf[2][i] = 1.0f;
  }
  
  float vBuf[3][3];
  float *v[4];
  for ( i=0; i<3; i++ )
    v[i+1]=vBuf[i]-1;

  float w[4];

  svdcmp( a, 3, 3, w, v );

  float x[4];
  float b[4];
  
  b[1] = pCorresp[12];
  b[2] = pCorresp[13];
  b[3] = 1.0f;

  svbksb( a, w, v, 3, 3, b, x );
  
  for ( i=0; i<3; i++ )
  {
    aBuf[0][i] = x[i+1]*pCorresp[i*4];
    aBuf[1][i] = x[i+1]*pCorresp[i*4+1];
    aBuf[2][i] = x[i+1];
  }

  gaussj( a, 3, 0, 0 );
  for ( i=0; i<12; i++ )
    T[i]=*(aBuf[0]+i);
  
//   svdcmp( a, 3, 3, w, v );
//   for ( i=1; i<4; i++ )
//     w[i] = 1.0f/w[i];
  
//   for ( i=0; i<3; i++ )
//     for ( j=0; j<3; j++ )
//   {
//     T[i*3+j] = 0.0f;
//     for ( k=0; k<3; k++ )
//       T[i*3+j] += vBuf[i][k]*w[k+1]*aBuf[j][k];
//   }
  
  float diff1 = ( T[0]*pCorresp[ 0]+T[1]*pCorresp[ 1]+T[2] - 1.0f );
  float diff2 = ( T[3]*pCorresp[ 0]+T[4]*pCorresp[ 1]+T[5] - 0.0f );
  float diff3 = ( T[6]*pCorresp[ 0]+T[7]*pCorresp[ 1]+T[8] - 0.0f );
  float diff4 = ( T[0]*pCorresp[ 4]+T[1]*pCorresp[ 5]+T[2] - 0.0f );
  float diff5 = ( T[3]*pCorresp[ 4]+T[4]*pCorresp[ 5]+T[5] - 1.0f );
  float diff6 = ( T[6]*pCorresp[ 4]+T[7]*pCorresp[ 5]+T[8] - 0.0f );
  float diff7 = ( T[0]*pCorresp[ 8]+T[1]*pCorresp[ 9]+T[9] - 0.0f );
  float diff8 = ( T[3]*pCorresp[ 8]+T[4]*pCorresp[ 9]+T[9] - 0.0f );
  float diff9 = ( T[6]*pCorresp[ 8]+T[7]*pCorresp[ 9]+T[9] - 1.0f );
  float diffA = ( T[0]*pCorresp[12]+T[1]*pCorresp[13]+T[2] - 1.0f )/10000.0f;
  float diffB = ( T[3]*pCorresp[12]+T[4]*pCorresp[13]+T[5] - 1.0f )/10000.0f;
  float diffC = ( T[6]*pCorresp[12]+T[7]*pCorresp[13]+T[8] - 1.0f )/10000.0f;

//   return
//     diff1*diff1 + diff2*diff2 + diff3*diff3 + diff4*diff4 +
//     diff5*diff5 + diff6*diff6 + diff7*diff7 + diff8*diff8 +
//     diff9*diff9 + diffA*diffA + diffB*diffB + diffC*diffC;
//   return
//     diff1 + diff2 + diff3 + diff4 +
//     diff5 + diff6 + diff7 + diff8 +
//     diff9 + diffA + diffB + diffC;
  return
    fabs(diff1) + fabs(diff2) + fabs(diff3) + fabs(diff4) +
    fabs(diff5) + fabs(diff6) + fabs(diff7) + fabs(diff8) +
    fabs(diff9) + fabs(diffA) + fabs(diffB) + fabs(diffC);
}

// void computeM( float *pCorresp, int numFeats,
// 	       float *e1, float *e2, float *M1, float *M2 )
// {
//   int i;
//   float T1[3][3], T2[3][3], tmpT1[3][3], tmpT2[3][3];

//   int c1, c2, c3, c4, cc1, cc2, cc3, cc4;
//   float tc[4][4];
//   float tmp;
//   float minProjBasisRes=FLT_MAX;
//   printf("Choosing minimum projective basis:\n");
// //   for ( c1=0; c1<numFeats-3; c1++ )
// //     for ( c2=c1+1; c2<numFeats-2; c2++ )
// //       for ( c3=c2+1; c3<numFeats-1; c3++ )
// // 	for ( c4=c3+1; c4<numFeats; c4++ )
// //   {
//   for ( c1=0; c1<numFeats; c1++ )
//     for ( c2=0; c2<numFeats; c2++ )
//       for ( c3=0; c3<numFeats; c3++ )
// 	for ( c4=0; c4<numFeats; c4++ )
//   {
//     if ( (c1==c2) || (c1==c3) || (c1==c4) ||
// 	 (c2==c3) || (c2==c4) || (c3==c4) ) continue;
//     for ( i=0; i<4; i++ )
//     {
//       tc[0][i] = pCorresp[c1*4+i];
//       tc[1][i] = pCorresp[c2*4+i];
//       tc[2][i] = pCorresp[c3*4+i];
//       tc[3][i] = pCorresp[c4*4+i];
//     }
//     float projBasisRes =
//       findProjBasisT( tc[0], tmpT1[0] ) +
//       findProjBasisT( tc[0]+2, tmpT2[0] );
//     float p1p[4][3];
//     for ( i=0; i<4; i++ )
//     {
//       p1p[i][0] = tmpT1[0][0]*pCorresp[i*4] + tmpT1[0][1]*pCorresp[i*4+1] + tmpT1[0][2];
//       p1p[i][1] = tmpT1[1][0]*pCorresp[i*4] + tmpT1[1][1]*pCorresp[i*4+1] + tmpT1[1][2];
//       p1p[i][2] = tmpT1[2][0]*pCorresp[i*4] + tmpT1[2][1]*pCorresp[i*4+1] + tmpT1[2][2];
//       printf("1:(%6.3f %6.3f %6.3f) ", p1p[i][0], p1p[i][1], p1p[i][2] );
//     }
//     printf("\n");
//     for ( i=0; i<4; i++ )
//     {
//       p1p[i][0] = tmpT2[0][0]*pCorresp[i*4+2] + tmpT2[0][1]*pCorresp[i*4+3] + tmpT2[0][2];
//       p1p[i][1] = tmpT2[1][0]*pCorresp[i*4+2] + tmpT2[1][1]*pCorresp[i*4+3] + tmpT2[1][2];
//       p1p[i][2] = tmpT2[2][0]*pCorresp[i*4+2] + tmpT2[2][1]*pCorresp[i*4+3] + tmpT2[2][2];
//       printf("2:(%6.3f %6.3f %6.3f) ", p1p[i][0], p1p[i][1], p1p[i][2] );
//     }
//     printf("\n");
  
//     if ( projBasisRes<minProjBasisRes )
//     {
//       minProjBasisRes=projBasisRes;
//       memcpy( T1[0], tmpT1[0], 9*sizeof(float) );
//       memcpy( T2[0], tmpT2[0], 9*sizeof(float) );
//       cc1 = c1;
//       cc2 = c2;
//       cc3 = c3;
//       cc4 = c4;
//     }
//   }
//   for ( i=0; i<4; i++ )
//   {
//     tmp = pCorresp[cc1*4+i];
//     pCorresp[cc1*4+i] = pCorresp[i];
//     pCorresp[i] = tmp;
//     tmp = pCorresp[cc2*4+i];
//     pCorresp[cc2*4+i] = pCorresp[4+i];
//     pCorresp[4+i] = tmp;
//     tmp = pCorresp[cc3*4+i];
//     pCorresp[cc3*4+i] = pCorresp[8+i];
//     pCorresp[8+i] = tmp;
//     tmp = pCorresp[cc4*4+i];
//     pCorresp[cc4*4+i] = pCorresp[12+i];
//     pCorresp[12+i] = tmp;
//   }
    
//   float p1p[5][3];
//   for ( i=0; i<5; i++ )
//   {
//     p1p[i][0] = T1[0][0]*pCorresp[i*4] + T1[0][1]*pCorresp[i*4+1] + T1[0][2];
//     p1p[i][1] = T1[1][0]*pCorresp[i*4] + T1[1][1]*pCorresp[i*4+1] + T1[1][2];
//     p1p[i][2] = T1[2][0]*pCorresp[i*4] + T1[2][1]*pCorresp[i*4+1] + T1[2][2];
//     printf("(%6.3f %6.3f %6.3f)\n", p1p[i][0], p1p[i][1], p1p[i][2] );
//   }
  
//   float e1p[3];
//   e1p[0] = T1[0][0]*e1[0] + T1[0][1]*e1[1] + T1[0][2]*e1[2];
//   e1p[0] = T1[1][0]*e1[0] + T1[1][1]*e1[1] + T1[1][2]*e1[2];
//   e1p[0] = T1[2][0]*e1[0] + T1[2][1]*e1[1] + T1[2][2]*e1[2];

//   float p2p[5][3];
//   for ( i=0; i<5; i++ )
//   {
//     p2p[i][0] = T2[0][0]*pCorresp[i*4+2] + T2[0][1]*pCorresp[i*4+3] + T2[0][2];
//     p2p[i][1] = T2[1][0]*pCorresp[i*4+2] + T2[1][1]*pCorresp[i*4+3] + T2[1][2];
//     p2p[i][2] = T2[2][0]*pCorresp[i*4+2] + T2[2][1]*pCorresp[i*4+3] + T2[2][2];
//     printf("(%6.3f %6.3f %6.3f)\n", p2p[i][0], p2p[i][1], p2p[i][2] );
//   }

//   float e2p[3];
//   e2p[0] = T2[0][0]*e2[0] + T2[0][1]*e2[1] + T2[0][2]*e2[2];
//   e2p[0] = T2[1][0]*e2[0] + T2[1][1]*e2[1] + T2[1][2]*e2[2];
//   e2p[0] = T2[2][0]*e2[0] + T2[2][1]*e2[1] + T2[2][2]*e2[2];

//   float p5cross[3];
//   p5cross[0] = p1p[4][1]*p2p[4][2] - p1p[4][2]*p2p[4][1];
//   p5cross[1] = p1p[4][2]*p2p[4][0] - p1p[4][0]*p2p[4][2];
//   p5cross[2] = p1p[4][0]*p2p[4][1] - p1p[4][1]*p2p[4][0];
  
//   float x1 = ( e2p[0]*p5cross[0] + e2p[1]*p5cross[1] + e2p[2]*p5cross[2] ) /
//              ( p1p[4][0]*e2p[0]*p5cross[0] +
// 	       p1p[4][1]*e2p[1]*p5cross[1] +
// 	       p1p[4][2]*e2p[2]*p5cross[2] );
//   float x2 = ( e1p[0]*p5cross[0] + e1p[1]*p5cross[1] + e1p[2]*p5cross[2] ) /
//              ( p2p[4][0]*e1p[0]*p5cross[0] +
// 	       p2p[4][1]*e1p[1]*p5cross[1] +
// 	       p2p[4][2]*e1p[2]*p5cross[2] );

//   for ( i=0; i<12; i++ )
//   {
//     M1[i]=0.0f;
//     M2[i]=0.0f;
//   }

//   M1[ 0] = p1p[4][0]*x1-1.0f;
//   M1[ 3] = 1.0f;
//   M1[ 5] = p1p[4][1]*x1-1.0f;
//   M1[ 7] = 1.0f;
//   M1[10] = p1p[4][2]*x1-1.0f;
//   M1[11] = 1.0f;
  
//   M2[ 0] = p2p[4][0]*x2-1.0f;
//   M2[ 3] = 1.0f;
//   M2[ 5] = p2p[4][1]*x2-1.0f;
//   M2[ 7] = 1.0f;
//   M2[10] = p2p[4][2]*x2-1.0f;
//   M2[11] = 1.0f;
// }

void computeFp( const float *F, float *Fp, float *e1, float *e2 )
{
  int i,j,k;
  
  float aBuf[3][3];
  float *a[4];
  for ( i=0; i<3; i++ )
  {
    a[i+1] = aBuf[i]-1;
    for ( j=0; j<3; j++ )
      aBuf[i][j] = F[i*3+j];
  }
  
  float vBuf[3][3];
  float *v[4];
  for ( i=0; i<3; i++ )
    v[i+1]=vBuf[i]-1;

  float w[4];

  svdcmp( a, 3, 3, w, v );

  w[3] = 0.0f;

  for ( i=0; i<3; i++ )
    for ( j=0; j<3; j++ )
  {
    Fp[i*3+j] = 0.0f;
    for ( k=0; k<3; k++ )
      Fp[i*3+j] += aBuf[i][k]*w[k+1]*vBuf[j][k];
  }

  for ( i=0; i<3; i++ )
    e1[i] = vBuf[2][i];
  
  for ( i=0; i<3; i++ )
    e2[i] = aBuf[i][2];
}

void computeF( const float *P1, const float *P2, float *F )
{
  F[0] =  DET_4X4( &P1[1], &P1[2], &P2[1], &P[2] );
  F[1] = -DET_4X4( &P1[0], &P1[2], &P2[1], &P[2] );
  F[2] =  DET_4X4( &P1[0], &P1[1], &P2[1], &P[2] );
  
  F[3] = -DET_4X4( &P1[1], &P1[2], &P2[0], &P[2] );
  F[4] =  DET_4X4( &P1[0], &P1[2], &P2[0], &P[2] );
  F[5] = -DET_4X4( &P1[0], &P1[1], &P2[0], &P[2] );
  
  F[6] =  DET_4X4( &P1[1], &P1[2], &P2[0], &P[1] );
  F[7] = -DET_4X4( &P1[0], &P1[2], &P2[0], &P[1] );
  F[8] =  DET_4X4( &P1[0], &P1[1], &P2[0], &P[1] );
}

void computeF( const vector<subWindowListAnim> &features,
	       int frame1, int frame2, int imgSizeX, int imgSizeY,
	       float *F )
{
  float x1, y1, x2, y2, compX, compY, compW;
  vector<subWindowListAnim>::const_iterator iter;
  int i, j;
  
  int numFeatures=0;
  for ( iter = features.begin(); iter!=features.end(); iter++ )
  {
    if ( (iter->includesFrame(frame1)) &&
	 (iter->includesFrame(frame2)) )
    { // feature is present in both images
      numFeatures++;
    }
  }
  printf("Using %i feature corespondences\n", numFeatures);

  float *pCBuf = (float*)malloc( numFeatures * 4 * sizeof(float) );
  float **pointCorresp = (float**)malloc( numFeatures * sizeof(float*) );
  for ( i=0; i<numFeatures; i++ )
    pointCorresp[i]=&pCBuf[i*4];
  for ( iter = features.begin(), i=0; iter!=features.end(); iter++ )
  {
    if ( (iter->includesFrame(frame1)) &&
	 (iter->includesFrame(frame2)) )
    { // feature is present in both images
      pointCorresp[i][0] = ((float)iter->getSubWindow(frame1).centerX()) /
	(float)imgSizeX * 2.0f - 1.0f;
      pointCorresp[i][1] = ((float)iter->getSubWindow(frame1).centerY()) /
	(float)imgSizeY * 2.0f - 1.0f;
      pointCorresp[i][2] = ((float)iter->getSubWindow(frame2).centerX()) /
	(float)imgSizeX * 2.0f - 1.0f;
      pointCorresp[i][3] = ((float)iter->getSubWindow(frame2).centerY()) /
	(float)imgSizeY * 2.0f - 1.0f;
      i++;
    }
  }
  
  float *aBuf = (float*)malloc( numFeatures * 9 * sizeof(float) );
  float **a = (float**)malloc( (numFeatures+1) * sizeof(float*) );
  for ( i=0; i<numFeatures; i++ )
    a[i+1] = &aBuf[i*9]-1;
  for ( i=0; i<numFeatures; i++ )
  {
    x1 = pointCorresp[i][0];
    y1 = pointCorresp[i][1];
    x2 = pointCorresp[i][2];
    y2 = pointCorresp[i][3];
    
    a[i+1][1] = x2*x1;
    a[i+1][2] = x2*y1;
    a[i+1][3] = x2;
    a[i+1][4] = y2*x1;
    a[i+1][5] = y2*y1;
    a[i+1][6] = y2;
    a[i+1][7] = x1;
    a[i+1][8] = y1;
    a[i+1][9] = 1.0f;
  }

  float vBuf[9][9];
  float *v[10];
  for ( i=0; i<9; i++ )
    v[i+1]=vBuf[i]-1;

  float w[10];

  svdcmp( a, numFeatures, 9, w, v );

  float f[3][3];

  for ( i=0; i<3; i++ )
    for ( j=0; j<3; j++ )
      F[i*3+j] = f[i][j]  = vBuf[8][i*3+j];

  printf("Fundamental matrix from SVD directly\n");
  for ( i=0; i<numFeatures; i++ )
  {
    x1 = pointCorresp[i][0];
    y1 = pointCorresp[i][1];
    x2 = pointCorresp[i][2];
    y2 = pointCorresp[i][3];
      
    compX = x1*f[0][0] + y1*f[1][0] + f[2][0];
    compY = x1*f[0][1] + y1*f[1][1] + f[2][1];
    compW = x1*f[0][2] + y1*f[1][2] + f[2][2];
    if (compW==0.0f) {printf("<><><>");compW=0.000001f;}
    compX *= x2;
    compY *= y2;
    compW *= 1.0f;

    printf( "diff %i: (%5.2f, %5.2f, %5.2f)\n",
	    i, compX, compY, compW );
  }

  free(aBuf);
  free(a);
  free(pCBuf);
  free(pointCorresp);
}

int main( int argc, char *argv[] )
{
  int width, height, length;
  vector<subWindowListAnim> features;
  if ( argc>=2 )
  {
    movieIO movie;
    printf( "Preparing \"%s\"...\n", argv[1] );
    if ( !movie.open( argv[1] ) )
    {
      printf( "Error in preparing %s\n", argv[1] );
      exit(0);
    }
    width=movie.getWidth();
    height=movie.getHeight();
    length=movie.getHeight();

    printf("\n");
    featureCorrelator corl( "Initializing correlator: ", movie, true );
    printf("\n");
    features = corl.buildTrackedFeatureList( "Tracking Features: ", 10, 5 );
    printf("\n");
  }
  else
  {
    printf("Loading previously tracked feature data...\n");
    #define FEATURE_LIST_NAME features
    #define WIDTH_VAR width
    #define HEIGHT_VAR height
    #define LENGTH_VAR length
    #include "track1.out"
    printf("Loaded %i tracked features\n", features.size() );
  }

  float ***F = malloc()
  float F[3][3], Fp[3][3], e1[3], e2[3], M1[3][4], M2[3][4];
  computeF( features, 17, 23, width, height,
	    F[0], Fp[0], e1, e2, M1[0], M2[0] );
  

//   percentBar bar(" Adding squares: " );
//   vector<subWindowListAnim>::iterator iter = features.begin();
//   unsigned char colour[3];
//   float curr=0.0f, end=features.size();
//   for ( ; iter!=features.end(); iter++ )
//   {
//     colour[0] = rand() % 256;
//     colour[1] = rand() % 256;
//     colour[2] = rand() % 256;
//     for ( int srcFrame=iter->getInterval().getEnd();
// 	  srcFrame>=iter->getInterval().getStart(); srcFrame-- )
//     {
//       rgbBuf = movie.getRGBFrame(srcFrame);
//       rgbBuf->addSquare( (*iter)[srcFrame], colour[0], colour[1], colour[2],
// 			 1, true );
//       movie.setRGBFrame( srcFrame, rgbBuf );
//       movie.releaseRGBFrame( rgbBuf );
//     }
//     bar.updatePercent( (curr+=1.0f)/end );
//   }
//   bar.endProgressBar();
  
//   printf("\n\n");
//   movie.play();
//   movie.close();
}