/*****************************************************************************
   This routine must ensure that the psImageFitPolynomial() and
   psImageEvalPolynomial() rotines work correctly.
 
   XXX: After you debug, set CHEBY_X_DIM != CHEBY_Y_DIM.
   XXX: use non-square images.
*****************************************************************************/
#include <stdio.h>
#include "pslib.h"
#include "psTest.h"
#include "psMemory.h"
#include "psImage.h"
#include "psImageStats.h"
#include <float.h>
#include <math.h>
#define IMAGE_SIZE 64
#define CHEBY_X_DIM 8
#define CHEBY_Y_DIM 8
#define THRESHOLD 10
psS32 main()
{
    psStats * myStats = NULL;
    psImage *tmpImage = NULL;
    psImage *outImage = NULL;
    psImage *nullImage = NULL;
    psPolynomial2D *my2DPoly = NULL;
    psPolynomial2D *null2DPoly = NULL;
    psS32 testStatus = true;
    psS32 memLeaks = 0;
    psS32 i = 0;
    psS32 j = 0;
    psS32 currentId = 0;

    currentId = psMemGetId();
    /*************************************************************************/
    /*  Allocate and initialize data structures                      */
    /*************************************************************************/
    tmpImage = psImageAlloc( IMAGE_SIZE, IMAGE_SIZE, PS_TYPE_F32 );
    outImage = psImageAlloc( IMAGE_SIZE, IMAGE_SIZE, PS_TYPE_F32 );
    for ( i = 0;i < IMAGE_SIZE;i++ ) {
        for ( j = 0;j < IMAGE_SIZE;j++ ) {
            tmpImage->data.F32[ i ][ j ] = 4.0;
            tmpImage->data.F32[ i ][ j ] = ( float ) ( i + j );
            tmpImage->data.F32[ i ][ j ] = ( float ) ( i + j ) + ( 4.0 * ( float ) i );
            outImage->data.F32[ i ][ j ] = 0.0;
        }
    }
    my2DPoly = psPolynomial2DAlloc( CHEBY_X_DIM, CHEBY_Y_DIM, PS_POLYNOMIAL_CHEB );
    /*************************************************************************/
    /*  Calculate Chebyshev Polynomials, no mask                     */
    /*************************************************************************/
    printPositiveTestHeader( stdout,
                             "psImageStats functions",
                             "Calculate Chebyshev Polynomials, no mask" );

    my2DPoly = psImageFitPolynomial( my2DPoly, tmpImage );
    for ( i = 0;i < CHEBY_X_DIM;i++ ) {
        for ( j = 0;j < CHEBY_Y_DIM;j++ ) {
            printf( "Cheby Polynomial (%d, %d) coefficient is %.2f\n", i, j, 0.0001f+my2DPoly->coeff[ i ][ j ] );
        }
    }

    psLogMsg(__func__,PS_LOG_INFO,"Following should be an error message for NULL coeffs argument.");
    null2DPoly = psImageFitPolynomial( NULL, tmpImage );
    if ( null2DPoly != NULL ) {
        printf("ERROR: psImageFitPolynomial did not return NULL with invalid coeffs argument.");
    }

    psLogMsg(__func__,PS_LOG_INFO,"Following should be a error message for NULL input argument.");
    null2DPoly = psImageFitPolynomial( null2DPoly, NULL );
    if ( null2DPoly != NULL ) {
        printf("ERROR: psImageFitPolynomial did not return NULL with null input image.");
    }

    psMemCheckCorruption( 1 );
    printFooter( stdout,
                 "psImageStats functions",
                 "Calculate Chebyshev Polynomials, no mask",
                 testStatus );
    /*************************************************************************/
    /*  Evaluate Chebyshev Polynomials, no mask                      */
    /*************************************************************************/
    printPositiveTestHeader( stdout,
                             "psImageStats functions",
                             "Calculate Chebyshev Polynomials, no mask" );

    outImage = psImageEvalPolynomial( outImage, my2DPoly );
    for ( i = 0;i < IMAGE_SIZE;i++ ) {
        for ( j = 0;j < IMAGE_SIZE;j++ ) {

            //             printf("pixel[%d][%d] is (%f, %f)\n", i, j,
            //                     tmpImage->data.F32[i][j],
            //                     outImage->data.F32[i][j]);
            if ( fabs( outImage->data.F32[ i ][ j ] - tmpImage->data.F32[ i ][ j ] ) > THRESHOLD ) {
                printf( "Pixel (%d, %d) is %.2f, should be %.2f\n", i, j,
                        outImage->data.F32[ i ][ j ],
                        tmpImage->data.F32[ i ][ j ] );
            }

        }
    }


    psLogMsg(__func__,PS_LOG_INFO,"Following should be an error message for NULL coeffs argument.");
    nullImage = psImageEvalPolynomial( outImage, NULL );
    if ( nullImage != NULL ) {
        printf("ERROR: psImageEvalPolynomial did not return NULL with invalid coeffs argument.");
    }

    psLogMsg(__func__,PS_LOG_INFO,"Following should be a error message for NULL input argument.");
    nullImage = psImageEvalPolynomial( NULL, my2DPoly );
    if ( nullImage != NULL ) {
        printf("ERROR: psImageEvalPolynomial did not return NULL with null input image.");
    }


    psMemCheckCorruption( 1 );
    printFooter( stdout,
                 "psImageStats functions",
                 "Calculate Chebyshev Polynomials, no mask",
                 testStatus );

    /*************************************************************************/
    /*  Deallocate data structures                                   */
    /*************************************************************************/
    printPositiveTestHeader( stdout,
                             "psImageStats functions",
                             "Deallocate the psStats/psImage structure." );
    psFree( myStats );
    psFree( tmpImage );
    psFree( outImage );
    psFree( my2DPoly );

    psMemCheckCorruption( 1 );
    memLeaks = psMemCheckLeaks(currentId,NULL,stderr,false);
    if ( 0 != memLeaks ) {
        psAbort( __func__, "Memory Leaks! (%d leaks)", memLeaks );
    }

    printFooter( stdout,
                 "psImageStats functions",
                 "Deallocate the psStats/psImage structure.",
                 testStatus );

    return ( !testStatus );
}
