cvsnakes.cpp.svn-base

来自「非结构化路识别」· SVN-BASE 代码 · 共 537 行 · 第 1/2 页

            float minEcont = _CV_SNAKE_BIG;
            float minEcurv = _CV_SNAKE_BIG;
            float minEimg = _CV_SNAKE_BIG;
            float Emin = _CV_SNAKE_BIG;

            int offsetx = 0;
            int offsety = 0;
            float tmp;

            /* compute bounds */
            int left = MIN( pt[i].x, win.width >> 1 );
            int right = MIN( roi.width - 1 - pt[i].x, win.width >> 1 );
            int upper = MIN( pt[i].y, win.height >> 1 );
            int bottom = MIN( roi.height - 1 - pt[i].y, win.height >> 1 );

            maxEcont = 0;
            minEcont = _CV_SNAKE_BIG;
            for( j = -upper; j <= bottom; j++ )
            {
                for( k = -left; k <= right; k++ )
                {
                    int diffx, diffy;
                    float energy;

                    if( i == 0 )
                    {
                        diffx = pt[n - 1].x - (pt[i].x + k);
                        diffy = pt[n - 1].y - (pt[i].y + j);
                    }
                    else
                    {
                        diffx = pt[i - 1].x - (pt[i].x + k);
                        diffy = pt[i - 1].y - (pt[i].y + j);
                    }
                    Econt[(j + centery) * win.width + k + centerx] = energy =
                        (float) fabs( ave_d -
                                      cvSqrt( (float) (diffx * diffx + diffy * diffy) ));

                    maxEcont = MAX( maxEcont, energy );
                    minEcont = MIN( minEcont, energy );
                }
            }
            tmp = maxEcont - minEcont;
            tmp = (tmp == 0) ? 0 : (1 / tmp);
            for( k = 0; k < neighbors; k++ )
            {
                Econt[k] = (Econt[k] - minEcont) * tmp;
            }

            /*  Calculate Ecurv */
            maxEcurv = 0;
            minEcurv = _CV_SNAKE_BIG;
            for( j = -upper; j <= bottom; j++ )
            {
                for( k = -left; k <= right; k++ )
                {
                    int tx, ty;
                    float energy;

                    if( i == 0 )
                    {
                        tx = pt[n - 1].x - 2 * (pt[i].x + k) + pt[i + 1].x;
                        ty = pt[n - 1].y - 2 * (pt[i].y + j) + pt[i + 1].y;
                    }
                    else if( i == n - 1 )
                    {
                        tx = pt[i - 1].x - 2 * (pt[i].x + k) + pt[0].x;
                        ty = pt[i - 1].y - 2 * (pt[i].y + j) + pt[0].y;
                    }
                    else
                    {
                        tx = pt[i - 1].x - 2 * (pt[i].x + k) + pt[i + 1].x;
                        ty = pt[i - 1].y - 2 * (pt[i].y + j) + pt[i + 1].y;
                    }
                    Ecurv[(j + centery) * win.width + k + centerx] = energy =
                        (float) (tx * tx + ty * ty);
                    maxEcurv = MAX( maxEcurv, energy );
                    minEcurv = MIN( minEcurv, energy );
                }
            }
            tmp = maxEcurv - minEcurv;
            tmp = (tmp == 0) ? 0 : (1 / tmp);
            for( k = 0; k < neighbors; k++ )
            {
                Ecurv[k] = (Ecurv[k] - minEcurv) * tmp;
            }

            /* Calculate Eimg */
            for( j = -upper; j <= bottom; j++ )
            {
                for( k = -left; k <= right; k++ )
                {
                    float energy;

                    if( scheme == _CV_SNAKE_GRAD )
                    {
                        /* look at map and check status */
                        int x = (pt[i].x + k) >> 3;
                        int y = (pt[i].y + j) >> 3;

                        if( map[y * map_width + x] == 0 )
                        {
                            int l, m;							
                            CvSize g_roi;
                            uchar *source;

                            /* evaluate block location */
                            int upshift = y ? 1 : 0;
                            int leftshift = x ? 1 : 0;
                            int bottomshift = MIN( 1, roi.height - ((y + 1) << 3) );
                            int rightshift = MIN( 1, roi.width - ((x + 1) << 3) );

                            source =
                                src + ((y << 3) - upshift) * srcStep + (x << 3) - leftshift;
                            g_roi.height = upshift + 8 + bottomshift;
                            g_roi.width = leftshift + 8 + rightshift;

                            
							icvSobelInitAlloc( g_roi.width, cv8u, 3, CV_ORIGIN_TL, 1, 0,
                                                &pX );
							icvSobelInitAlloc( g_roi.width, cv8u, 3, CV_ORIGIN_TL, 0, 1,
                                                &pY );

                            icvSobel_8u16s_C1R( source, srcStep, dx, 20, &g_roi, pX, 0 );
                            icvSobel_8u16s_C1R( source, srcStep, dy, 20, &g_roi, pY, 0 );

                            icvFilterFree( &pX );
                            icvFilterFree( &pY );
                            for( l = 0; l < 8 + bottomshift; l++ )
                            {
                                for( m = 0; m < 8 + rightshift; m++ )
                                {
                                    gradient[((y << 3) + l) * roi.width + (x << 3) + m] =
                                        (float) (dx[(l + upshift) * 10 + m + leftshift] *
                                                 dx[(l + upshift) * 10 + m + leftshift] +
                                                 dy[(l + upshift) * 10 + m + leftshift] *
                                                 dy[(l + upshift) * 10 + m + leftshift]);
                                }
                            }
                            map[y * map_width + x] = 1;
                        }
                        Eimg[(j + centery) * win.width + k + centerx] = energy =
                            gradient[(pt[i].y + j) * roi.width + pt[i].x + k];
                    }
                    else
                    {
                        Eimg[(j + centery) * win.width + k + centerx] = energy =
                            src[(pt[i].y + j) * srcStep + pt[i].x + k];
                    }

                    maxEimg = MAX( maxEimg, energy );
                    minEimg = MIN( minEimg, energy );
                }
            }

            tmp = (maxEimg - minEimg);
            tmp = (tmp == 0) ? 0 : (1 / tmp);

            for( k = 0; k < neighbors; k++ )
            {
                Eimg[k] = (minEimg - Eimg[k]) * tmp;
            }

            /* locate coefficients */
            if( coeffUsage == CV_VALUE)
            {
                _alpha = *alpha;
                _beta = *beta;
                _gamma = *gamma;
            }
            else
            {                   
                _alpha = alpha[i];
                _beta = beta[i];
                _gamma = gamma[i];
            }

            /* Find Minimize point in the neighbors */
            for( k = 0; k < neighbors; k++ )
            {
                E[k] = _alpha * Econt[k] + _beta * Ecurv[k] + _gamma * Eimg[k];
            }
            Emin = _CV_SNAKE_BIG;
            for( j = -upper; j <= bottom; j++ )
            {
                for( k = -left; k <= right; k++ )
                {

                    if( E[(j + centery) * win.width + k + centerx] < Emin )
                    {
                        Emin = E[(j + centery) * win.width + k + centerx];
                        offsetx = k;
                        offsety = j;
                    }
                }
            }

            if( offsetx || offsety )
            {
                pt[i].x += offsetx;
                pt[i].y += offsety;
                moved++;
            }
        }
        converged = (moved == 0);
        if( (criteria.type & CV_TERMCRIT_ITER) && (iteration >= criteria.maxIter) )
            converged = 1;
        if( (criteria.type & CV_TERMCRIT_EPS) && (moved <= criteria.epsilon) )
            converged = 1;
    }

    icvFree( &Econt );
    icvFree( &Ecurv );
    icvFree( &Eimg );
    icvFree( &E );

    if( scheme == _CV_SNAKE_GRAD )
    {
        icvFree( &gradient );
        icvFree( &map );
        icvFree( &dx );
        icvFree( &dy );
    }
    return CV_OK;
}


CV_IMPL void
cvSnakeImage( const IplImage* src, CvPoint* points,
              int length, float *alpha,
              float *beta, float *gamma,
              CvCoeffType coeffUsage, CvSize win,
              CvTermCriteria criteria, int calcGradient )
{

    CV_FUNCNAME( "cvSnakeImage" );

    __BEGIN__;

    uchar *data;
    CvSize size;
    int step;

    if( src->nChannels != 1 )
        CV_ERROR( CV_BadNumChannels, "input image has more than one channel" );

    if( src->depth != IPL_DEPTH_8U )
        CV_ERROR( CV_BadDepth, icvUnsupportedFormat );


    cvGetImageRawData( src, &data, &step, &size );

    if( calcGradient )
    {
        IPPI_CALL( icvSnakeImageGrad8uC1R( data, step, size, points, length,
                              alpha, beta, gamma, coeffUsage, win, criteria ));
    }
    else
    {
        IPPI_CALL( icvSnakeImage8uC1R( data, step, size, points, length, alpha, beta,
                                       gamma, coeffUsage, win, criteria ));
    }
    __CLEANUP__;
    __END__;
}

/* end of file */