fnimage.cpp

在人脸检测的基础之上,对嘴部的运动表情进行分析,进行语音模拟.

// FnImage.cpp: implementation of the FnImage class.
//
//////////////////////////////////////////////////////////////////////

#include "stdafx.h"
#include "FnImage.h"
#include "Matrix.h"
#include "FnMatrix.h"
#include "FnMath.h"
#include <math.h>
#include <time.h>
#include <stdio.h>

#define M_PI 3.14159265
#define MAX_RAND        (pow(2.0,31.0)-1.0)     //(2**31)-1

static FnMath fMath;

static CWnd   fWnd;

//--------------------------------------------------------------------------
// Finds the luminance of a color image 
// all Red components first, Green components next 
// and Blue Components last.
//--------------------------------------------------------------------------
matrix FindLum( const matrix &image)
{
	ASSERT( image.ysize > 0 && image.xsize > 0);
	ASSERT( image.numbands == 3);
	
	double coef[ 3] = { 0.299, 0.587, 0.114};
	
	matrix lum( image.ysize, image.xsize, 1);
	
	for(int u = 0; u < 3; u++) 
	{
		int offset = u * image.ysize * image.xsize;
		for(int p = 0; p < image.ysize * image.xsize; p++) 
		{
			lum.elem[ p] += coef[ u] * image.elem[ p + offset];
		}
	}
	
	return lum;
	
}

//--------------------------------------------------------------------------
matrix  FindLumRGB(const matrix &image)  
//--------------------------------------------------------------------------
// Finds the luminance of a color image 
// RGB component together side by side 

{
	
	int           xsi,ysi,i,tot;
	double       *t1,*t2;
	
	xsi = image.xsize ;   
	ysi = image.ysize;
	matrix luminance(ysi,xsi/3);
	tot = xsi*ysi/3;      
	t1  = image.elem;
	t2  = luminance.elem;
	
	for(i = 0 ; i < tot ; ++i) 
    {
		*t2++ = 0.299 * (*t1) + 0.587 * (*(t1+1)) + 0.114 * (*(t1+2));
		t1   += 3;  
    }
	return(luminance);
	
}

//------------------------------------------------------------------------
matrix   rotate(matrix &a1, float angle ,int direction, char type)
//------------------------------------------------------------------------
// angle is given in degrees
// enter 'n' for nearest neighbourhood interpolation
// enter 'b' for bilinear interpolation during rotation
// direction is given as 1 for rotation in the clockwise direction
// direction is given as -1 for rotation in the counter clockwise direction.

{
	
	matrix  b,a;
	int     i,j,xshift,yshift,x1,y1;
	double  x,y;
	
	a = a1;
	/*
	int BLUR = 5;
	matrix  fil(1,BLUR);
	for (i = 1; i <= BLUR ; ++i)
	{
    fil(0,i-1) = 1.0 / BLUR;  // for uniform blur
	}
	a = (a)|(fil);
	a = (a)|(!fil);
	*/
	b = zeros(a.ysize,a.xsize);
	xshift = (int) ((a.xsize / 2 ));
	yshift = (int) ((a.ysize / 2 ));
	angle = (float) (angle * PI / 180);
	
	if (type == 'n')
    {
		for (j = 0; j < b.ysize; ++j)
		{
			for (i = 0; i < b.xsize; ++i)
			{
				x = cos(angle) * (i-xshift) + direction * sin(angle) * (j-yshift); 
				y = -direction * sin(angle) * (i-xshift) + cos(angle) * (j-yshift);
				x = x + xshift;
				y = y + yshift;
				x1 = (int)(x);
				y1 = (int)(y);
				if ((x1 >= 0) && (x1 < a.xsize) && (y1 >= 0) && (y1 < a.ysize))
				{
					b(j,i) = a(y1,x1);
				}
			}
		}
    }
	else if  (type == 'b')
    {
		double  f1,f2;
		for (j = 0; j < b.ysize; ++j)
		{
			for (i = 0; i < b.xsize; ++i)
			{
				x = cos(angle) * (i-xshift) + direction * sin(angle) * (j-yshift); 
				y = -direction * sin(angle) * (i-xshift) + cos(angle) * (j-yshift);
				x = x + xshift;
				y = y + yshift;
				x1 = (int)(x);
				y1 = (int)(y);
				if ((x1 > 0) && (x1 < (a.xsize-1)) && (y1 > 0) && (y1 < (a.ysize-1)))
				{
					f1 = x - floor(x);
					f2 = y - floor(y);
					b(j,i)  = a(y1-1,x1-1) * (1-f1) * (1-f2);
					b(j,i) += a(y1-1,x1+1) * (f1)   * (1-f2);
					b(j,i) += a(y1+1,x1-1) * (1-f1) * (f2);
					b(j,i) += a(y1+1,x1+1) * f1     * f2    ;
				}
			}
		}
    }
	
	return b;
}
//-----------------------------------------------------------
matrix  RotateImage(matrix &image, double ang)
//-----------------------------------------------------------
{
	// rotate image frm1 about its center by angle ang (in degrees)
	// we assume  -pi/2 <= ang <= pi/2
	// from the book by George Woolberg    
	int     x, y, wmax, newh, neww, h, w,flag=0;
	double  sine, tangent, offst, cosine; 
	double *dst, *src;
	
	matrix frm1 = image; 
	
	if ((ang >= -90) && (ang <0))
    {
		frm1 = flip_vertical(frm1);
		ang = fabs(ang); 
		flag =1; 
    }
	
	h = frm1.ysize; w = frm1.xsize;
	
	
	/* the dimensions of the rotated image as it is processed are
	(h)(w) -> (h)(wmax) -> (newh)(wmax) -> (newh)(neww)
	+1 will be added to dimensions due to last fractional pixel
	Temporary image TMP is used to hold the intermediate image
	*/
	sine = sin(ang*PI/180.0);
	cosine = cos(ang*PI/180.0);
	tangent = tan(ang*PI/(180.0*2.0));
	
	wmax = (int)(w + h*tangent + 1);       /* width of intermediate image */
	newh = (int)(w*sine + h*cosine + 1);   /* final image height*/
	neww = (int)(h*sine + w*cosine + 1);   /* final image width*/
	matrix frm2(h,wmax); 
	
	/* 1st pass: skew x(horizontal scanlines) */
	
	for (y = 0; y < h; y++)                        /* visit each row in frm1 */
    {
		src = frm1.rowptr[y];                  /* input scanline pointer */ 
		dst = frm2.rowptr[y];                  /* output scanline pointer */
		skew(src, w, wmax, y*tangent,1, dst); /* skew row */
    }
	/* 2nd pass:skew y(vertical scanlines). Use TMP for intermediate image*/
	
	matrix  TMP(newh,wmax);
	offst = (w-1)*sine;               /* offset from top of image */
	for (x = 0; x < wmax; x++)        /* visit each column in frm2 */
    {
		src = &frm2.elem[x];          /* input scanline pointer */
		dst = &TMP.elem[x];           /* output scanline pointer */ 
		skew(src, h, newh, offst-x*sine, wmax, dst);  /* skew column */
    }
	/* 3rd pass : skew x(horizontal scanlines)*/ 
	
	frm2.resize(newh,neww);
	for (y = 0; y < newh; y++)       /* visit each row in TMP*/
    {
		src = TMP.rowptr[y];
		dst = frm2.rowptr[y]; 
		skew(src, wmax, neww, (y-offst)*tangent, 1, dst);
    }
	if (flag)
    {
		frm2 =  flip_vertical(frm2);
    }
	return frm2;  
}
// translates  the image with the parameters (d1, d2) in the x and y directions
//-----------------------------------------------------------------------
matrix          shift_integer(matrix &a,  int d1, int d2)
//-----------------------------------------------------------------------
{
	matrix b;
	
	//b = zeros(a.ysize, a.xsize, a.numbands);
	b = a;
	int i, j, k;
	
	for (i = 0; i < b.xsize ; ++i)
	{
		for (j = 0; j < b.ysize ; ++j)
		{
			if (((i-d1) >=0) && ((i-d1) < a.xsize) && ((j-d2) >=0) && ((j-d2) < a.ysize))
			{
				for ( k = 0; k < b.numbands; k++)
				{
					
					b(j, i, k) = a(j - d2, i - d1, k);
				}
			}
		}
	}
	return b;
}
//------------------------------------------------------------------
void   skew(double *src,int len, int nlen, double strt, int offst,double *dst)
//------------------------------------------------------------------
// skew scanline in src(length len) into dst(length nlen)
// starting at position strt. The offset between each 
// scanline pixel is offst, offst = 1 for rows ; 
// offst = width for columns
{
	int     i, istrt, lim;
	double  f, g, w1,w2;
	
	/* process left end of output: either prepare for clipping or add padding */
	
	istrt = (int) strt;                   /* integer index                           */
	if (istrt < 0) src -= (offst*istrt);  /* advance input ptr for clipping          */
	lim = MIN(len+istrt, nlen);           /* find index for right edge               */
	for (i=0; i < istrt; i++)             /* visit all null ouput pixels at left edge*/
    {
		*dst = 0;                         /* pad with zero                           */
		dst += offst;                     /* advance output pointer                  */
    }
	f = fabs(strt - istrt);               /* weight for right straddle               */
	g = 1.-f;                             /* weight for left straddle                */
	
	if (f==0.)                            /* simple integer shift: no interpolation  */
    {
		for(;i < lim; i++)                /* visit all pixels in valid range         */
		{
			*dst = *src;                  /* copy input to output                    */
			src += offst;                 /* advance input pointer                   */
			dst += offst;                 /* advance output pointer                  */ 
		}
    }
	else                                  /* fractional shift:interpolate            */
    {
		if (strt > 0.)
		{
			w1 = f;                       /* weight for left pixel                   */ 
			w2 = g;                       /* weight for right pixel                  */ 
			*dst = g*src[0];              /* first pixel                             */
			dst += offst;                 /* advance ouput pointer                   */
			i++;
		}                                 /* increment index                         */ 
		else 
		{
			w1 = g;                       /* weight for left pixel                   */ 
			w2 = f;                       /* weight for right pixel                  */ 
			if (lim < nlen) lim--;
		}
		for (; i < lim; i++)              /* visit all pixels in valid range         */
		{
			/* src[0] is left (top) pixel, and src[offst] is right(bottom) pixel     */
			*dst = w1 * src[0] + w2 * src[offst];  /* linear interpolation           */ 
			dst += offst;                 /* advance output pointer                  */  
			src += offst;                 /* advance input pointer                   */
		}
		if (i < nlen)
		{
			*dst = w1 * src[0];          /* src[0] is last pixel                     */
			dst += offst;                /* advance output pointer                   */  
			i++;
		}
    }
	
	for (; i <nlen; i++)                 /* visit all remaining pixels at right edge */
    {
		*dst = 0;                        /* pad with zero                            */
		dst += offst;                    /* advance output pointer                   */ 
    }
	
}
//-----------------------------------------------------------
matrix flip_vertical( matrix &a)
//-----------------------------------------------------------
{
	int      i,j;
	matrix   temp(a.ysize, a.xsize);
	
	for (j=0; j < temp.ysize; ++j)
		for (i=0; i < temp.xsize; ++i)
			temp(j,i) = a(temp.ysize - 1 - j,i); 
		return(temp); 
		
}


//-----------------------------------------------------------
matrix flip_horizontal( matrix &a)
//-----------------------------------------------------------
{
	int      i,j;
	matrix   temp(a.ysize, a.xsize);
	
	for (i=0; i < temp.xsize; ++i)
		for (j=0; j < temp.ysize; ++j)
			temp(j,i) = a(j, temp.ysize -1 - i); 
		return(temp); 
}


//-------------------------------------------------------------------------
matrix kucult(matrix& img, int xs, int ys)
//-------------------------------------------------------------------------
{
	int i, j, at_x, at_y;
	matrix out_image(ys,xs);
	
	at_x = (img.xsize-xs)/2;
	at_y = (img.ysize-ys)/2; 
	
	for (j = 0; j < ys; ++j)
		for (i = 0; i < xs; ++i)
			out_image(j,i) = img(j + at_y, i + at_x);
		
	return (out_image);
}

//--------------------------------------------------------------
matrix ScaleImage(matrix &im, double xs, double ys)
//--------------------------------------------------------------
// uses bilinear interpolation 
{
	// xs : the scale factor in the x direction
	// ys : the scale factor in the y direction
	int     i, j, ci;
	double  c, cd, f; 
	matrix  temp(im.ysize, im.xsize),temp2(im.ysize, im.xsize);
	ASSERT ((xs > 0 ) && (ys > 0));
	
	for (i = 0; i < temp.xsize; ++i) //scale in x-direction
    {
		c = (double)i/xs;
		ci = (int) c;
		cd = c - ci;
		f = 1-cd;
		for (j = 0; j < temp.ysize; ++j)
		{
			temp(j,i) = ((ci+1) < im.xsize)?(im(j,ci)*f + im(j,ci+1)*cd): 0.0;
		}
    }
	for (j = 0; j < temp.ysize; ++j) //scale in y-direction
    {
		c = (double)j/ys;
		ci = (int) c;
		cd = c - ci;
		f = 1-cd;
		for (i = 0; i < temp.xsize; ++i)
		{
			temp2(j,i) = ((ci+1) < temp.ysize)?(temp(ci,i)*f + temp(ci+1,i)*cd):0.0;
		}
    }
	return(temp2); 
}
//-------------------------------------------------------------
matrix ZoomImage( matrix &image, int x_factor, int y_factor)
//------------------------------------------------------------
{
	// Enlarges the image by replicating rows and columns

	ASSERT( (image.xsize > 0) && (image.ysize >0));
	ASSERT( (x_factor > 1) && (y_factor > 1) );
	ASSERT( (x_factor < 5) && (y_factor < 5));

	matrix temp = image, temp2;
	int y, x, k;

	temp = UpsampleRow( temp, y_factor);
	for ( y = 0; y < temp.ysize; y++)
	{
		for ( x = 0; x < temp.xsize; x++)
		{
			for ( k = 0; k < image.numbands; k++ )
			{
				temp(y, x, k) = image( y / y_factor, x, k); 
			}
		}
	}

	temp2 = UpsampleColumn( temp, x_factor);
	for ( y = 0; y < temp2.ysize; y++)
	{
		for ( x = 0; x < temp2.xsize; x++)
		{
			for ( k = 0; k < image.numbands; k++ )
			{
				temp2(y, x, k) = temp( y, x / x_factor, k); 
			}
		}
	}
	return temp2;
}

//------------------------------------------------------------------------
double maximum(matrix &a , matrix& b)
//------------------------------------------------------------------------

{	// returns the value of the maximum element 
	// and returns the location of the maximum 
	// element in b 
	ASSERT( (b.xsize == 2 && b.ysize == 1) || (b.ysize == 2 && b.xsize == 1) ); 
	
	int     i,j;
	double  max;
	
	max = a(0,0);
	b(0) = 0;
	b(1) = 0; 
	
	for ( j = 0; j < a.ysize ; ++j)
	{
		for (i = 0; i < a.xsize ; ++i)
		{
			if ( a(j,i)  > max)
			{
				max = a(j,i);
				b(0) = i;
				b(1) = j;
			}
		}
	}
	return max;
}

//------------------------------------------------------------------------
double maximum3D(matrix &a , matrix& b)
//------------------------------------------------------------------------

{	// returns the value of the maximum element of the 3D matrix a
	// and returns the location of the maximum 
	// element in b 
	ASSERT( (b.xsize == 3 && b.ysize == 1) || (b.ysize == 3 && b.xsize == 1) ); 
	
	int     i,j,k;
	double  max;
	
	max = a(0,0,0);
	b(0) = 0;
	b(1) = 0; 
	b(2) = 0;
	
	for ( j = 0; j < a.ysize ; ++j)
	{
		for (i = 0; i < a.xsize ; ++i)
		{
			for (k = 0; k < a.numbands; k++)
			{
				if ( a(j,i,k)  > max)
				{
					max = a(j,i,k);
					b(0) = j;
					b(1) = i;
					b(2) = k;
				}
			}
		}
	}
	return max;
}
//------------------------------------------------------------------------
double maximum(matrix &a )
//------------------------------------------------------------------------

{	// returns the value of the maximum element 
	// and returns the location of the maximum 
	// element in b 
	ASSERT( a.numbands == 1);
	int     i,j;
	double  max;
	
	max = a(0,0);
	
	for ( j = 0; j < a.ysize ; ++j)
	{
		for (i = 0; i < a.xsize ; ++i)
		{
			if ( a(j,i)  > max)
			{
				max = a(j,i);
			}
		}
	}
	return max;
}

//------------------------------------------------------------------------
double minimum(matrix &a ,  matrix& b)
//------------------------------------------------------------------------

{
	ASSERT( (b.xsize == 2 && b.ysize == 1) || (b.ysize == 2 && b.xsize == 1) ); 
	int     i,j;
	double  min;
	
	min = a(0,0);
	b(0) = 0;
	b(1) = 0;