cxmathfuncs.cpp

opencv库在TI DM6437上的移植,目前包括两个库cv.lib和cxcore.lib的工程

        if( !CV_IS_MAT(ymat))
            CV_CALL( ymat = cvGetMat( ymat, &ystub, &coi2 ));        

        if( !CV_ARE_TYPES_EQ( angle, ymat ) )
            CV_ERROR_FROM_CODE( CV_StsUnmatchedFormats );
        
        if( !CV_ARE_SIZES_EQ( angle, ymat ) )
            CV_ERROR_FROM_CODE( CV_StsUnmatchedSizes );
        
        cont_flag &= ymat->type;
    }

    if( coi1 != 0 || coi2 != 0 || coi3 != 0 || coi4 != 0 )
        CV_ERROR( CV_BadCOI, "" );

    size = cvGetMatSize(angle);
    size.width *= CV_MAT_CN(angle->type);

    if( CV_IS_MAT_CONT( cont_flag ))
    {
        size.width *= size.height;
        size.height = 1;
    }

    block_size = MIN( size.width, ICV_MATH_BLOCK_SIZE );
    x_buffer = (float*)cvStackAlloc( block_size*sizeof(float));
    y_buffer = (float*)cvStackAlloc( block_size*sizeof(float));

    if( depth == CV_32F )
    {
        for( y = 0; y < size.height; y++ )
        {
            float* x_data = (float*)(xmat ? xmat->data.ptr + xmat->step*y : 0);
            float* y_data = (float*)(ymat ? ymat->data.ptr + ymat->step*y : 0);
            float* mag_data = (float*)(mag ? mag->data.ptr + mag->step*y : 0);
            float* angle_data = (float*)(angle->data.ptr + angle->step*y);

            for( x = 0; x < size.width; x += block_size )
            {
                int i, len = MIN( size.width - x, block_size );

                icvSinCos_32f( angle_data+x, y_buffer, x_buffer, len, angle_in_degrees );

                for( i = 0; i < len; i++ )
                {
                    float tx = x_buffer[i];
                    float ty = y_buffer[i];

                    if( mag_data )
                    {
                        float magval = mag_data[x+i];
                        tx *= magval;
                        ty *= magval;
                    }

                    if( xmat )
                        x_data[x+i] = tx;
                    if( ymat )
                        y_data[x+i] = ty;
                }
            }
        }
    }
    else
    {
        for( y = 0; y < size.height; y++ )
        {
            double* x_data = (double*)(xmat ? xmat->data.ptr + xmat->step*y : 0);
            double* y_data = (double*)(ymat ? ymat->data.ptr + ymat->step*y : 0);
            double* mag_data = (double*)(mag ? mag->data.ptr + mag->step*y : 0);
            double* angle_data = (double*)(angle->data.ptr + angle->step*y);
            double C = angle_in_degrees ? CV_PI/180. : 1;

            for( x = 0; x < size.width; x++ )
            {
                double phi = angle_data[x]*C;
                double magval = mag_data ? mag_data[x] : 1.;
                if( xmat )
                    x_data[x] = cos(phi)*magval;
                if( ymat )
                    y_data[x] = sin(phi)*magval;
            }
        }
    }

    __END__;
}

/****************************************************************************************\
*                                          E X P                                         *
\****************************************************************************************/

typedef union
{
    struct {
#if ( defined( WORDS_BIGENDIAN ) && !defined( OPENCV_UNIVERSAL_BUILD ) ) || defined( __BIG_ENDIAN__ )
        int hi;
        int lo;
#else
        int lo;
        int hi;
#endif
    } i;
    double d;
}
DBLINT;

#define EXPTAB_SCALE 6
#define EXPTAB_MASK  (1 << EXPTAB_SCALE) - 1

#define EXPPOLY_32F_A0 .9670371139572337719125840413672004409288e-2

static const double icvExpTab[] = {
    1.0 * EXPPOLY_32F_A0,
    1.0108892860517004600204097905619 * EXPPOLY_32F_A0,
    1.0218971486541166782344801347833 * EXPPOLY_32F_A0,
    1.0330248790212284225001082839705 * EXPPOLY_32F_A0,
    1.0442737824274138403219664787399 * EXPPOLY_32F_A0,
    1.0556451783605571588083413251529 * EXPPOLY_32F_A0,
    1.0671404006768236181695211209928 * EXPPOLY_32F_A0,
    1.0787607977571197937406800374385 * EXPPOLY_32F_A0,
    1.0905077326652576592070106557607 * EXPPOLY_32F_A0,
    1.1023825833078409435564142094256 * EXPPOLY_32F_A0,
    1.1143867425958925363088129569196 * EXPPOLY_32F_A0,
    1.126521618608241899794798643787 * EXPPOLY_32F_A0,
    1.1387886347566916537038302838415 * EXPPOLY_32F_A0,
    1.151189229952982705817759635202 * EXPPOLY_32F_A0,
    1.1637248587775775138135735990922 * EXPPOLY_32F_A0,
    1.1763969916502812762846457284838 * EXPPOLY_32F_A0,
    1.1892071150027210667174999705605 * EXPPOLY_32F_A0,
    1.2021567314527031420963969574978 * EXPPOLY_32F_A0,
    1.2152473599804688781165202513388 * EXPPOLY_32F_A0,
    1.2284805361068700056940089577928 * EXPPOLY_32F_A0,
    1.2418578120734840485936774687266 * EXPPOLY_32F_A0,
    1.2553807570246910895793906574423 * EXPPOLY_32F_A0,
    1.2690509571917332225544190810323 * EXPPOLY_32F_A0,
    1.2828700160787782807266697810215 * EXPPOLY_32F_A0,
    1.2968395546510096659337541177925 * EXPPOLY_32F_A0,
    1.3109612115247643419229917863308 * EXPPOLY_32F_A0,
    1.3252366431597412946295370954987 * EXPPOLY_32F_A0,
    1.3396675240533030053600306697244 * EXPPOLY_32F_A0,
    1.3542555469368927282980147401407 * EXPPOLY_32F_A0,
    1.3690024229745906119296011329822 * EXPPOLY_32F_A0,
    1.3839098819638319548726595272652 * EXPPOLY_32F_A0,
    1.3989796725383111402095281367152 * EXPPOLY_32F_A0,
    1.4142135623730950488016887242097 * EXPPOLY_32F_A0,
    1.4296133383919700112350657782751 * EXPPOLY_32F_A0,
    1.4451808069770466200370062414717 * EXPPOLY_32F_A0,
    1.4609177941806469886513028903106 * EXPPOLY_32F_A0,
    1.476826145939499311386907480374 * EXPPOLY_32F_A0,
    1.4929077282912648492006435314867 * EXPPOLY_32F_A0,
    1.5091644275934227397660195510332 * EXPPOLY_32F_A0,
    1.5255981507445383068512536895169 * EXPPOLY_32F_A0,
    1.5422108254079408236122918620907 * EXPPOLY_32F_A0,
    1.5590044002378369670337280894749 * EXPPOLY_32F_A0,
    1.5759808451078864864552701601819 * EXPPOLY_32F_A0,
    1.5931421513422668979372486431191 * EXPPOLY_32F_A0,
    1.6104903319492543081795206673574 * EXPPOLY_32F_A0,
    1.628027421857347766848218522014 * EXPPOLY_32F_A0,
    1.6457554781539648445187567247258 * EXPPOLY_32F_A0,
    1.6636765803267364350463364569764 * EXPPOLY_32F_A0,
    1.6817928305074290860622509524664 * EXPPOLY_32F_A0,
    1.7001063537185234695013625734975 * EXPPOLY_32F_A0,
    1.7186192981224779156293443764563 * EXPPOLY_32F_A0,
    1.7373338352737062489942020818722 * EXPPOLY_32F_A0,
    1.7562521603732994831121606193753 * EXPPOLY_32F_A0,
    1.7753764925265212525505592001993 * EXPPOLY_32F_A0,
    1.7947090750031071864277032421278 * EXPPOLY_32F_A0,
    1.8142521755003987562498346003623 * EXPPOLY_32F_A0,
    1.8340080864093424634870831895883 * EXPPOLY_32F_A0,
    1.8539791250833855683924530703377 * EXPPOLY_32F_A0,
    1.8741676341102999013299989499544 * EXPPOLY_32F_A0,
    1.8945759815869656413402186534269 * EXPPOLY_32F_A0,
    1.9152065613971472938726112702958 * EXPPOLY_32F_A0,
    1.9360617934922944505980559045667 * EXPPOLY_32F_A0,
    1.9571441241754002690183222516269 * EXPPOLY_32F_A0,
    1.9784560263879509682582499181312 * EXPPOLY_32F_A0,
};

static const double exp_prescale = 1.4426950408889634073599246810019 * (1 << EXPTAB_SCALE);
static const double exp_postscale = 1./(1 << EXPTAB_SCALE);
static const double exp_max_val = 3000.*(1 << EXPTAB_SCALE); // log10(DBL_MAX) < 3000

IPCVAPI_IMPL( CvStatus, icvExp_32f, ( const float *_x, float *y, int n ), (_x, y, n) )
{
    static const double
        EXPPOLY_32F_A4 = 1.000000000000002438532970795181890933776 / EXPPOLY_32F_A0,
        EXPPOLY_32F_A3 = .6931471805521448196800669615864773144641 / EXPPOLY_32F_A0,
        EXPPOLY_32F_A2 = .2402265109513301490103372422686535526573 / EXPPOLY_32F_A0,
        EXPPOLY_32F_A1 = .5550339366753125211915322047004666939128e-1 / EXPPOLY_32F_A0;

    #undef EXPPOLY
    #define EXPPOLY(x)  \
        (((((x) + EXPPOLY_32F_A1)*(x) + EXPPOLY_32F_A2)*(x) + EXPPOLY_32F_A3)*(x) + EXPPOLY_32F_A4)

    int i = 0;
    DBLINT buf[4];
    const Cv32suf* x = (const Cv32suf*)_x;

    if( !x || !y )
        return CV_NULLPTR_ERR;
    if( n <= 0 )
        return CV_BADSIZE_ERR;

    buf[0].i.lo = buf[1].i.lo = buf[2].i.lo = buf[3].i.lo = 0;

    for( ; i <= n - 4; i += 4 )
    {
        double x0 = x[i].f * exp_prescale;
        double x1 = x[i + 1].f * exp_prescale;
        double x2 = x[i + 2].f * exp_prescale;
        double x3 = x[i + 3].f * exp_prescale;
        int val0, val1, val2, val3, t;

        if( ((x[i].i >> 23) & 255) > 127 + 10 )
            x0 = x[i].i < 0 ? -exp_max_val : exp_max_val;

        if( ((x[i+1].i >> 23) & 255) > 127 + 10 )
            x1 = x[i+1].i < 0 ? -exp_max_val : exp_max_val;

        if( ((x[i+2].i >> 23) & 255) > 127 + 10 )
            x2 = x[i+2].i < 0 ? -exp_max_val : exp_max_val;

        if( ((x[i+3].i >> 23) & 255) > 127 + 10 )
            x3 = x[i+3].i < 0 ? -exp_max_val : exp_max_val;

        val0 = cvRound(x0);
        val1 = cvRound(x1);
        val2 = cvRound(x2);
        val3 = cvRound(x3);

        x0 = (x0 - val0)*exp_postscale;
        x1 = (x1 - val1)*exp_postscale;
        x2 = (x2 - val2)*exp_postscale;
        x3 = (x3 - val3)*exp_postscale;

        t = (val0 >> EXPTAB_SCALE) + 1023;
        t = (t | ((t < 2047) - 1)) & (((t < 0) - 1) & 2047);
        buf[0].i.hi = t << 20;

        t = (val1 >> EXPTAB_SCALE) + 1023;
        t = (t | ((t < 2047) - 1)) & (((t < 0) - 1) & 2047);
        buf[1].i.hi = t << 20;

        t = (val2 >> EXPTAB_SCALE) + 1023;
        t = (t | ((t < 2047) - 1)) & (((t < 0) - 1) & 2047);
        buf[2].i.hi = t << 20;

        t = (val3 >> EXPTAB_SCALE) + 1023;
        t = (t | ((t < 2047) - 1)) & (((t < 0) - 1) & 2047);
        buf[3].i.hi = t << 20;

        x0 = buf[0].d * icvExpTab[val0 & EXPTAB_MASK] * EXPPOLY( x0 );
        x1 = buf[1].d * icvExpTab[val1 & EXPTAB_MASK] * EXPPOLY( x1 );
        
        y[i] = (float)x0;
        y[i + 1] = (float)x1;
        
        x2 = buf[2].d * icvExpTab[val2 & EXPTAB_MASK] * EXPPOLY( x2 );
        x3 = buf[3].d * icvExpTab[val3 & EXPTAB_MASK] * EXPPOLY( x3 );
        
        y[i + 2] = (float)x2;
        y[i + 3] = (float)x3;
    }

    for( ; i < n; i++ )
    {
        double x0 = x[i].f * exp_prescale;
        int val0, t;

        if( ((x[i].i >> 23) & 255) > 127 + 10 )
            x0 = x[i].i < 0 ? -exp_max_val : exp_max_val;

        val0 = cvRound(x0);
        t = (val0 >> EXPTAB_SCALE) + 1023;
        t = (t | ((t < 2047) - 1)) & (((t < 0) - 1) & 2047);

        buf[0].i.hi = t << 20;
        x0 = (x0 - val0)*exp_postscale;

        y[i] = (float)(buf[0].d * icvExpTab[val0 & EXPTAB_MASK] * EXPPOLY(x0));
    }

    return CV_OK;
}


IPCVAPI_IMPL( CvStatus, icvExp_64f, ( const double *_x, double *y, int n ), (_x, y, n) )
{
    static const double
        A5 = .99999999999999999998285227504999 / EXPPOLY_32F_A0,
        A4 = .69314718055994546743029643825322 / EXPPOLY_32F_A0,
        A3 = .24022650695886477918181338054308 / EXPPOLY_32F_A0,
        A2 = .55504108793649567998466049042729e-1 / EXPPOLY_32F_A0,
        A1 = .96180973140732918010002372686186e-2 / EXPPOLY_32F_A0,
        A0 = .13369713757180123244806654839424e-2 / EXPPOLY_32F_A0;

    #undef EXPPOLY
    #define EXPPOLY(x)  (((((A0*(x) + A1)*(x) + A2)*(x) + A3)*(x) + A4)*(x) + A5)

    int i = 0;
    DBLINT buf[4];
    const Cv64suf* x = (const Cv64suf*)_x;

    if( !x || !y )
        return CV_NULLPTR_ERR;
    if( n <= 0 )
        return CV_BADSIZE_ERR;