aoptflowhs.cpp

微软的基于HMM的人脸识别原代码, 非常经典的说

/*M///////////////////////////////////////////////////////////////////////////////////////
//
//  IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
//
//  By downloading, copying, installing or using the software you agree to this license.
//  If you do not agree to this license, do not download, install,
//  copy or use the software.
//
//
//                        Intel License Agreement
//                For Open Source Computer Vision Library
//
// Copyright (C) 2000, Intel Corporation, all rights reserved.
// Third party copyrights are property of their respective owners.
//
// Redistribution and use in source and binary forms, with or without modification,
// are permitted provided that the following conditions are met:
//
//   * Redistribution's of source code must retain the above copyright notice,
//     this list of conditions and the following disclaimer.
//
//   * Redistribution's in binary form must reproduce the above copyright notice,
//     this list of conditions and the following disclaimer in the documentation
//     and/or other materials provided with the distribution.
//
//   * The name of Intel Corporation may not be used to endorse or promote products
//     derived from this software without specific prior written permission.
//
// This software is provided by the copyright holders and contributors "as is" and
// any express or implied warranties, including, but not limited to, the implied
// warranties of merchantability and fitness for a particular purpose are disclaimed.
// In no event shall the Intel Corporation or contributors be liable for any direct,
// indirect, incidental, special, exemplary, or consequential damages
// (including, but not limited to, procurement of substitute goods or services;
// loss of use, data, or profits; or business interruption) however caused
// and on any theory of liability, whether in contract, strict liability,
// or tort (including negligence or otherwise) arising in any way out of
// the use of this software, even if advised of the possibility of such damage.
//
//M*/

#include "CvTest.h"

/* Testing parameters */
static char FuncName[] = "cvCalcOpticalFlowHS";
static char TestName[] = "Optical flow (Horn & Schunck)";
static char TestClass[] = "Algorithm";

static long lImageWidth;
static long lImageHeight;
static long  lNumIterations;
static float flEpsilon;
static float lambda;

#define EPSILON 0.0001f

static int fmaCalcOpticalFlowHS( void )
{
    /* Some Variables */
    int            i,j,k;
    
    uchar*         roiA;
    uchar*         roiB;

    float*         VelocityX;
    float*         VelocityY;

    float*         auxVelocityX;
    float*         auxVelocityY;

    float*         DerX;
    float*         DerY;
    float*         DerT;

    long            lErrors = 0;

    CvTermCriteria criteria;

    int usePrevious;

    int Stop = 0;
    int iteration = 0;
    float epsilon = 0;

    static int  read_param = 0;

    /* Initialization global parameters */
    if( !read_param )
    {
        read_param = 1;
        /* Reading test-parameters */
        trslRead( &lImageHeight, "300", "Image height" );
        trslRead( &lImageWidth, "300", "Image width" );
        trssRead( &lambda, "20", "lambda" );
    }

    /* initialization - for warning disable */
    criteria.epsilon = 0;
    criteria.maxIter = 0;
    criteria.type  = 1;

    /* Allocating memory for all frames */
    IplImage* imgA = cvCreateImage( cvSize(lImageWidth,lImageHeight), IPL_DEPTH_8U, 1 );
    IplImage* imgB = cvCreateImage( cvSize(lImageWidth,lImageHeight), IPL_DEPTH_8U, 1 );

    IplImage* testVelocityX = cvCreateImage( cvSize(lImageWidth,lImageHeight), IPL_DEPTH_32F, 1 );
    IplImage* testVelocityY = cvCreateImage( cvSize(lImageWidth,lImageHeight), IPL_DEPTH_32F, 1 );

    VelocityX = (float*)icvAlloc(  lImageWidth*lImageHeight * sizeof(float) );
    VelocityY = (float*)icvAlloc(  lImageWidth*lImageHeight * sizeof(float) );

    auxVelocityX = (float*)icvAlloc(  lImageWidth*lImageHeight * sizeof(float) );
    auxVelocityY = (float*)icvAlloc(  lImageWidth*lImageHeight * sizeof(float) );

    DerX = (float*)icvAlloc( lImageWidth*lImageHeight * sizeof(float) );
    DerY = (float*)icvAlloc( lImageWidth*lImageHeight * sizeof(float) );
    DerT = (float*)icvAlloc( lImageWidth*lImageHeight * sizeof(float) );

    /* Filling images */
    ats1bInitRandom( 0, 255, (uchar*)imgA->imageData, lImageWidth * lImageHeight );
    ats1bInitRandom( 0, 255, (uchar*)imgB->imageData, lImageWidth * lImageHeight );

    /* set ROI of images */
    roiA = (uchar*)imgA->imageData;
    roiB = (uchar*)imgB->imageData;

    /* example of 3*3 ROI*/
    /*roiA[0] = 0;
    roiA[1] = 1;
    roiA[2] = 2;
    roiA[lImageWidth] = 0;
    roiA[lImageWidth+1] = 1;
    roiA[lImageWidth+2] = 2;
    roiA[2*lImageWidth] = 0;
    roiA[2*lImageWidth+1] = 1;
    roiA[2*lImageWidth+2] = 2;

    roiB[0] = 1;
    roiB[1] = 2;
    roiB[2] = 3;
    roiB[lImageWidth] = 1;
    roiB[lImageWidth+1] = 2;
    roiB[lImageWidth+2] = 3;
    roiB[2*lImageWidth] = 1;
    roiB[2*lImageWidth+1] = 2;
    roiB[2*lImageWidth+2] = 3;*/
/****************************************************************************************\
*                  Calculate derivatives                                                 *
\****************************************************************************************/
    for (i=0; i<lImageHeight; i++)
    {
        for(j=0; j<lImageWidth; j++)
        {
            int jr,jl,it,ib;

            if ( j==lImageWidth-1 )
                jr = lImageWidth-1;
            else jr = j + 1;

            if ( j==0 )
                jl = 0;
            else jl = j - 1;

            if ( i==(lImageHeight - 1) )
                ib = lImageHeight - 1;
            else ib = i + 1;

            if ( i==0 )
                it = 0;
            else it = i - 1;

            DerX[ i*lImageWidth + j ] = (float)
                (roiA[ (it)*imgA->widthStep + jr ]
                - roiA[ (it)*imgA->widthStep + jl ]
                + 2*roiA[ (i)*imgA->widthStep + jr ]
                - 2*roiA[ (i)*imgA->widthStep + jl ]
                + roiA[ (ib)*imgA->widthStep + jr ]
                - roiA[ (ib)*imgA->widthStep + jl ])/8 ;

            DerY[ i*lImageWidth + j ] = (float)
                ( roiA[ (ib)*imgA->widthStep + jl ]
                + 2*roiA[ (ib)*imgA->widthStep + j  ]
                + roiA[ (ib)*imgA->widthStep + jr ]
                - roiA[ (it)*imgA->widthStep + jl ]
                - 2*roiA[ (it)*imgA->widthStep + j  ]
                - roiA[ (it)*imgA->widthStep + jr ])/8  ;

            DerT[ i*lImageWidth + j ] = (float)
                (roiB[i*imgB->widthStep + j] - roiA[i*imgA->widthStep + j]);
        }
    }
for( usePrevious = 0; usePrevious < 2; usePrevious++ )
{
/****************************************************************************************\
*                    Cases                                                               *
\****************************************************************************************/
    for ( k = 0; k < 4; k++ )
    {
        switch (k)
        {
        case 0:
            {
                criteria.type = CV_TERMCRIT_ITER;
                criteria.maxIter = 3;

                trsWrite( ATS_LST|ATS_CON,
                         "usePrevious = %d, criteria = ITER, maxIter = %d\n",
                         usePrevious, criteria.maxIter);

                break;
            }
        case 1:
            {
                criteria.type = CV_TERMCRIT_EPS;
                criteria.epsilon = 0.001f;