C++ HalfEdge::getNextEdge方法代码示例

T FacePerimeter:: getFacePerimeter ( const HalfEdge <T> *e_start )
{
        //Get perimeter of the Face given by start half edge
        T perimeter = 0;
        HalfEdge <T> *e = const_cast <HalfEdge <T> *> ( e_start );

        //There is a valid edge
        if ( e_start != NULL )
        {
                //Get point
                Node3DCartesian <T> *pi = e->getPoint();

                //Proces all edges of the Face
                do
                {
                        //Get point
                        Node3DCartesian <T> *pii = e->getNextEdge()->getPoint();

                        //Compute area
                        perimeter += EuclDistance::getEuclDistance2D ( pi->getX(), pi->getY(), pii->getX(), pii->getY() );

                        //Assign point
                        pi = pii;

                        //Assign edge
                        e = e->getNextEdge();

                }
                while ( e != e_start );
        }

        //Get area
        return perimeter;
}

unsigned int Face <T> ::countVertices() const
{
        //Count vertices of the Face
        unsigned int vertices = 0;

        //There is a valid edge
        if ( edge != NULL )
        {
                //Get next edge in cell
                HalfEdge <T> *e = edge;

                //Proces all edges of the Face
                do
                {
                        //Increment vertices
                        vertices ++;

                        //Increment edge
                        e = e->getNextEdge();

                }
                while ( e != edge );
        }

        return vertices;
}

void Face <T>::removeAdjacency()
{
        //Set twin edges of Faces adjacent to actual Face to NULL
        if ( edge != NULL )
        {
                //Get actual edge in cell
                HalfEdge <T> *e = edge;

                //Proces all edges of the Face
                do
                {
                        //Get twin edge
                        HalfEdge <T> *e_twin = e->getTwinEdge();

                        //Adjacent cell exists
                        if ( e_twin != NULL )
                        {
                                //Set pointer to NULL
                                e_twin -> setTwinEdge ( NULL );
                        }

                        //Increment edge
                        e = e->getNextEdge();

                }
                while ( e != edge );
        }
}

T DDT2D::globalCostFunction ( Container <HalfEdge <T> *> &half_edges )
{
        // Compute cost of the triangulation using selected criterion
        const unsigned int n = half_edges.size();
        T cost = 0;

        //Loop all edges
        for ( unsigned int i = 0; i < n; i++ )
        {

                //Take half edge
                HalfEdge <T> *e = half_edges [i];

                //Use simplex indentificator to eliminate T processing of the edge
                if ( !e->isSimplexEdge() )
                {
                        //Get the first triangle
                        const HalfEdge <T> *e12 = e->getNextEdge();
                        const HalfEdge <T> *e13 = e12->getNextEdge();

                        //Does a twin edge exist?
                        if ( e->getTwinEdge() )
                        {
                                //Get coordinates of the first triangle
                                const Point3DCartesian <T> *p11 = e->getPoint();
                                const Point3DCartesian <T> *p12 = e12->getPoint() ;
                                const Point3DCartesian <T> *p13 =  e13->getPoint() ;

                                //Get second triangle
                                HalfEdge <T> *e21 = e->getTwinEdge();
                                const HalfEdge <T> *e22 = e21->getNextEdge();
                                const HalfEdge <T> *e23 = e22->getNextEdge();

                                //Get coordinates
                                const Point3DCartesian <T> *p21 = p12;
                                const Point3DCartesian <T> *p22 = p11;
                                const Point3DCartesian <T> *p23 = e23->getPoint();

                                //Compute local criterion
                                T fi = ( *pcriterion ) ( p11, p12, p13, p21, p22, p23 );

                                //Global cost
                                cost += fi;

                                //Set both edges to be processed
                                e->setEdgeAsSimplex ( true );
                                e21->setEdgeAsSimplex ( true );
                        }
                }
        }

        //Reset attribute
        for ( unsigned int i = 0; i < n; i++ )
        {
                ( *half_edges ) [i]->setEdgeAsSimplex ( false );
        }

        //Return global cost of the triangulation
        return cost;
}

void TangentFunction::computeTangentFunctionFace ( const Face <T> *face, typename TTangentFunction <T>::Type & tf, T & rotation, const TTangentFunctionRotationMethod & rotation_method,
	const TTangentFunctionScaleMethod & scale_method )
{
        //Compute normalized tangent function for the Face  given by HalfEdge
        T turning_angle_sum = 0, normalize_length_sum = 0;
        const HalfEdge <T> *e_start = face->getHalfEdge();
        HalfEdge <T> *e = const_cast <HalfEdge <T> *> ( e_start ) ;

        //Get triplet of points
        Node3DCartesian <T> *pi, *piii = NULL;
        Node3DCartesian <T> *pii = e_start->getPoint();

        //Jump short segments
        for ( pi = e->getPreviousEdge()->getPoint(); EuclDistance::getEuclDistance2D ( pi, pii ) < MIN_POSITION_DIFF ; e = e->getPreviousEdge(), pi = e->getPreviousEdge()->getPoint() ) {}

        //Get third point for initial edge
        for ( piii = e->getNextEdge()->getPoint(); EuclDistance::getEuclDistance2D ( pii, piii ) < MIN_POSITION_DIFF ; e = e->getNextEdge(), piii = e->getNextEdge()->getPoint() ) {}

        //Compute perimeter of the Face
        const T face_perimeter = FacePerimeter::getFacePerimeter ( e_start );

        //Throw exception
        if ( face_perimeter == 0 )
        {
                throw ErrorMathZeroDevision <T> ( "ErrorMathZeroDevision: tangent function, can not normalize function (1 / perimeter), ", "perimeter = 0." );
        }

        //Compute initial angle: rotation invariant
        if ( rotation_method == RotationInvariant )
        {
                turning_angle_sum += 180 - Angle3Points::getAngle3Points ( pi, pii, piii );
        }

        //Compute initial angle: rotation dependent (do not compute the sum)
        else
        {
                turning_angle_sum = 360 - Angle3Points::getAngle3Points ( &Node3DCartesian <T> ( pii->getX() - 1, pii->getY() ), pii, piii );
        }

        //Add turning angle computed for normalized length of initial edge to map
        tf[normalize_length_sum] = turning_angle_sum;

        //Assign points from initial edge
        pi = pii;
        pii = piii;

        //Increment initial edge
        e = e_start->getNextEdge();

        //Process all edges of the Face (2, n-1)
        do
        {
                //Get third point
                piii = e->getNextEdge()->getPoint();

                //Segment is too short
                if ( EuclDistance::getEuclDistance2D ( pii, piii ) > MIN_POSITION_DIFF )
                {
                        //Compute angle: rotation invariant
                        if ( rotation_method == RotationInvariant )
                        {
                                turning_angle_sum += 180 - Angle3Points::getAngle3Points ( pi, pii, piii );
                        }

                        //Compute initial angle: rotation dependent (do not compute the sum)
                        else
                        {
                                turning_angle_sum = 360 - Angle3Points::getAngle3Points ( &Node3DCartesian <T> ( pii->getX() - 1, pii->getY() ), pii, piii );
                        }

                        //Compute normalized length sum
                        normalize_length_sum += ( scale_method == ScaleDependent ? EuclDistance::getEuclDistance2D ( pi, pii ) / face_perimeter : EuclDistance::getEuclDistance2D ( pi, pii ) );

                        //Add turning angle computed for normalized length to map
                        tf[normalize_length_sum] = turning_angle_sum;

                        //Assign points
                        pi = pii;
                        pii = piii;
                }

                //Increment edge
                e = e->getNextEdge();

        } while ( e != e_start );

        //Add last value: we visit first point again
	const T last = ( scale_method == ScaleDependent ? 1.0 : face_perimeter );
        if ( rotation_method == RotationInvariant )
        {
		tf[ last] = turning_angle_sum + tf.begin()->second;
        }

        else
        {
                tf[last] = tf.begin()->second;
        }
}

void DDT2D::DDTLOP ( Container <Node3DCartesian <T> *> &nl, Container <HalfEdge <T> *> &half_edges, unsigned short swap_criterion_selected, const bool print_message, const bool print_exception, std::ostream * output )
{
        // Data Depending triangulation using selected local criterion
        bool swap_exist = true;

        //Set iterations
        unsigned iterations = 0;

        //Create Delaunay triangulation
        DT2D::DT ( nl, half_edges, print_message );

        //Get number of HalfEdges
        const unsigned int n = half_edges->size();

        //Set local swap criterion
        setSwapCriterion ( swap_criterion_selected );

        //Global cost before swapping
        const T global_cost_old = globalCostFunction ( half_edges );
        T global_cost = global_cost_old;

        //Print info
        if ( print_message )
        {
                *output << "> Starting DDTLOP... " ;
        }

        //Initialize counters
        unsigned int counter = MAX_INT, counter_old = MAX_INT;

        try
        {
                //Run until swap exists or decrease number of swaps between two loops
                do
                {
                        //We suppose ordered set of triangles, no swap will be required
                        swap_exist = false;

                        //Remember old counter
                        counter_old = counter;

                        //Assign new counter value
                        counter = 0;

                        //Loop all edges
                        for ( unsigned int i = 0; i < n; i++ )
                        {

                                //Take half edge
                                HalfEdge <T> *e = ( *half_edges ) [i];

                                //Use simplex flag to eliminate T processing of the edge
                                if ( !e->isSimplexEdge() )
                                {
                                        //Does twin edge exist?
                                        if ( e->getTwinEdge() )
                                        {
                                                // Test of convexity for quadrilateral
                                                HalfEdge <T> *e12 = e->getNextEdge();
                                                HalfEdge <T> *e13 = e12->getNextEdge();
                                                HalfEdge <T> *e21 = e->getTwinEdge();
                                                HalfEdge <T> *e22 = e21->getNextEdge();
                                                HalfEdge <T> *e23 = e22->getNextEdge();

                                                // Get nodes, counterclockwise set of nodes
                                                const Node3DCartesian <T> *p1 = e->getPoint();
                                                const Node3DCartesian <T> *p2 = e23->getPoint();
                                                const Node3DCartesian <T> *p3 = e12->getPoint();
                                                const Node3DCartesian <T> *p4 = e13->getPoint();

                                                //Is convex (non convex can not be swapped)
                                                if ( ConvexQuadrilateral::isStrictlyConvex ( p1, p2, p3, p4 ) == 1 )
                                                {

                                                        //Set twin edge to be processed
                                                        e->getTwinEdge()->setEdgeAsSimplex ( true );

                                                        //Get first triangle
                                                        e12 = e->getNextEdge();
                                                        e13 = e12->getNextEdge();

                                                        //Get coordinates (cast to parent using static_cast)
                                                        const Point3DCartesian <T> *p11 = e->getPoint();
                                                        const Point3DCartesian <T> *p12 = e12->getPoint();
                                                        const Point3DCartesian <T> *p13 = e13->getPoint();

                                                        //Get second triangle
                                                        e21 = e->getTwinEdge();
                                                        e22 = e21->getNextEdge();
                                                        e23 = e22->getNextEdge();

                                                        //Get coordinates
                                                        const Point3DCartesian <T> *p21 = p12;
                                                        const Point3DCartesian <T> *p22 = p11;
                                                        const Point3DCartesian <T> *p23 = e23->getPoint();

                                                        //Compute local criterion
                                                        const T c1 = ( *pcriterion ) ( p11, p12, p13, p21, p22, p23 );

                                                        //Calculation local criterion from swapped diagonal
//.........这里部分代码省略.........

float DDT2D::getInitialTemperature ( Container <HalfEdge <T> *> &half_edges )
{
        // Compute cost of the triangulation using selected criterion
        const unsigned int n = half_edges.size();
        T max_difference = 0; ;

        //Loop all edges
        for ( unsigned int i = 0; i < n; i++ )
        {

                //Take half edge
                HalfEdge <T> *e = half_edges [i];

                if ( !e->isSimplexEdge() )
                {
                        //Get next edges
                        const HalfEdge <T> *e12 = e->getNextEdge();
                        const HalfEdge <T> *e13 = e12->getNextEdge();

                        //Does a twin edge exist?
                        if ( e->getTwinEdge() )
                        {
                                //Get coordinates of the first triangle (cast to parent using static_cast)
                                const Point3DCartesian <T> *p11 = e->getPoint();
                                const Point3DCartesian <T> *p12 = e12->getPoint();
                                const Point3DCartesian <T> *p13 = e13->getPoint();

                                //Get second triangle
                                HalfEdge <T> *e21 = e->getTwinEdge();
                                HalfEdge <T> *e22 = e21->getNextEdge();
                                HalfEdge <T> *e23 = e22->getNextEdge();

                                //Get coordinates
                                const Point3DCartesian <T> *p21 = p12;
                                const Point3DCartesian <T> *p22 = p11;
                                const Point3DCartesian <T> *p23 = e23->getPoint();

                                //Only strictly convex quadrilaterals
                                if ( ConvexQuadrilateral::isStrictlyConvex ( p11,  p23, p12,  p13 ) == 1 )
                                {
                                        //Compute local criterion  for adjacent triangles
                                        const T fi1 = ( *pcriterion ) ( p11, p12, p13, p21, p22, p23 );

                                        //Compute local criterion  for swapped triangles
                                        const T fi2 = ( *pcriterion ) ( p13, p23, p12, p23, p13, p11 );

                                        //Find max difference
                                        if ( fabs ( fi2 - fi1 ) > max_difference )
                                        {
                                                max_difference = fabs ( fi2 - fi1 );
                                        }
                                }

                                //Set both edges to be processed
                                e->setEdgeAsSimplex ( true );
                                e21->setEdgeAsSimplex ( true );
                        }
                }
        }

        //Reset attribute
        for ( unsigned int i = 0; i < n; i++ )
        {
                ( *half_edges ) [i]->setEdgeAsSimplex ( false );
        }

        //Return result
        return ( float ) ( 2 * max_difference );
}

void DDT2D::setNewState ( float tk, unsigned int & good_swap, unsigned int n, T & global_cost, Container <HalfEdge <T> *> &half_edges )
{
        //Set new state or recover old state
        unsigned int i = int ( n * rand() / ( RAND_MAX + 1.0 ) );

        //i can not be grater than n
        if ( i >= n )
        {
                i = n - 1;
        }

        //Get random edge
        HalfEdge <T> *e = ( *half_edges ) [i];

        //Cost difference
        T cost_difference = 0.0001 * MAX_FLOAT;

        //Compute cost difference (convex quadrilateral), else set difference to MAX_FLOAT
        getCostDifference ( e, cost_difference, 0, 1 );

        //Set new state, perform swap of the edge
        if ( cost_difference < 0 )
        {

                //First triangle
                HalfEdge <T> *e12 = e->getNextEdge();
                HalfEdge <T> *e13 = e12->getNextEdge();

                //Get the second triangle T2
                HalfEdge <T> *e21 = e->getTwinEdge();
                HalfEdge <T> *e22 = e21->getNextEdge();
                HalfEdge <T> *e23 = e22->getNextEdge();

                //Swap diagonal
                DT2D::swapDiagonal ( e, e12, e13, e21, e22, e23 );

                //Increment good swap
                good_swap ++;

                //Assign old global cost
                global_cost += cost_difference;
        }

        //Decide if set a new state or old state
        else
        {

                //Generate random number (0, 1)
                const float theta = ( T ) rand() / ( T ) RAND_MAX;

                //Compare with Boltzmann function and decide about new state
                if ( theta < exp ( - cost_difference / tk ) )
                {

                        //Set new state (acceptable increasing of the cost), perform the swap

                        //First triangle
                        HalfEdge <T> *e12 = e->getNextEdge();
                        HalfEdge <T> *e13 = e12->getNextEdge();

                        //Get the second triangle T2
                        HalfEdge <T> *e21 = e->getTwinEdge();
                        HalfEdge <T> *e22 = e21->getNextEdge();
                        HalfEdge <T> *e23 = e22->getNextEdge();

                        //Swap diagonal
                        DT2D::swapDiagonal ( e, e12, e13, e21, e22, e23 );

                        //Assign old global cost
                        global_cost += cost_difference;
                }
        }
}

unsigned short PointFacePosition:: getPointFacePosition ( const Point3DCartesian <T> * q, const Face <T> *face )
{
        /*Returns position of the point to the Face
        0: q is strictly interior
        1: q is strictly exterior
        2: q is on the edge but not an endpoint
        3: q is a vertex
        */

        //There is a valid Face
        if ( face != NULL )
        {
                unsigned short l_inters = 0, r_inters = 0;

                HalfEdge <T> *e = face->getHalfEdge();
                const HalfEdge <T> *e_start = e;

                //Get point P[i-1]
                const Node3DCartesian <T> *pi1 = e->getPreviousEdge()->getPoint();

                //Process all edges
                do
                {
                        //Get point P[i]
                        const Node3DCartesian <T> *pi = e->getPoint();

                        //Test point is identical with end points: d(q, pi) < POSITION_ROUND_ERROR
                        if ( ( q->getX() - pi->getX() ) * ( q->getX() - pi->getX() ) + ( q->getY() - pi->getY() ) * ( q->getY() - pi->getY() ) < POSITION_ROUND_ERROR * POSITION_ROUND_ERROR )
                        {
                                return 3;
                        }

                        //Perform  tests: different position of segment end points to ray
                        bool t1 = ( pi->getY() > q->getY() + POSITION_ROUND_ERROR ) != ( pi1->getY() > q->getY() + POSITION_ROUND_ERROR );
                        bool t2 = ( pi->getY() < q->getY() - POSITION_ROUND_ERROR ) != ( pi1->getY() < q->getY() - POSITION_ROUND_ERROR );

                        //Segment intersects left or ray
                        if ( t1 || t2 )
                        {
                                //Compute intersection
                                T x = ( ( pi->getX() - q->getX() ) * ( pi1->getY() - q->getY() ) - ( pi1->getX() - q->getX() ) * ( pi->getY() - q->getY() ) )
                                      / double ( ( pi1->getY() - pi->getY() ) );

                                //Aditional test to avoid roundness error: point too close to edge
                                if ( PointLineDistance::getPointLineSegmentDistance2D ( q, pi, pi1 ) < POSITION_ROUND_ERROR )
                                {
                                        return 2;
                                }

                                //Increment number of right intersections
                                if ( t1 && x > 0 )
                                {
                                        r_inters ++;
                                }

                                //Increment  number of left intersections
                                if ( t2 && x < 0 )
                                {
                                        l_inters ++;
                                }
                        }

                        //Assign point
                        pi1 = pi;

                        //Increment edge
                        e = e->getNextEdge();

                }
                while ( e != e_start );

                //Point lies on the edge, but not at endpoints
                if ( ( l_inters % 2 ) != ( r_inters % 2 ) )
                {
                        return 2;
                }

                //Points is strictly interior
                if ( l_inters % 2 == 1 )
                {
                        return 0;
                }

                //Point is strictly exterior
                return 1;
        }

        //Throw exception
        else
        {
                throw ErrorBadData ( "ErrorBadData: no face incident to edge, ", "can not test point vs. face position" );
        }
}

void TurningFunction::computeTurningFunctionFace ( const Face <T> *face, typename TTurningFunction <T>::Type & tf, T & rotation, const TTurningFunctionRotationMethod & rotation_method,
                const TTurningFunctionScaleMethod & scale_method )
{
        //Compute normalized turning function for the Face  given by HalfEdge
        T turning_angle_sum = 0, normalize_length_sum = 0;
        const HalfEdge <T> *e_start = face->getHalfEdge();
        HalfEdge <T> *e = const_cast <HalfEdge <T> *> ( e_start ) ;

        //Get triplet of points
        Node3DCartesian <T> *pi, *piii = NULL;
        Node3DCartesian <T> *pii = e_start->getPoint();

        //Jump short segments
        for ( pi = e->getPreviousEdge()->getPoint(); EuclDistance::getEuclDistance2D ( pi->getX(), pi->getY(), pii->getX(), pii->getY() ) < MIN_POSITION_DIFF ; e = e->getPreviousEdge(), pi = e->getPreviousEdge()->getPoint() );

        //Get third point for initial edge
        for ( piii = e->getNextEdge()->getPoint(); EuclDistance::getEuclDistance2D ( pii->getX(), pii->getY(), piii->getX(), piii->getY() ) < MIN_POSITION_DIFF ; e = e->getNextEdge(), piii = e->getNextEdge()->getPoint() );

        //Compute perimeter of the Face
        const T face_perimeter = FacePerimeter::getFacePerimeter ( e_start );

        //Throw exception
        if ( fabs ( face_perimeter ) < MIN_FLOAT )
        {
                throw MathZeroDevisionException <T> ( "MathZeroDevisionException: turning function, can not normalize function (1 / perimeter), ", "perimeter = 0.", face_perimeter );
        }

        //Compute initial angle: rotation invariant
        if ( rotation_method == RotationInvariant )
        {
                turning_angle_sum += 180 - Angle3Points::getAngle3Points ( pi, pii, piii );
        }

        //Compute initial angle: rotation dependent (do not compute the sum)
        else
        {
                //Create temporary point P (x_temp, y_temp) on x axis
                const T dx = std::max ( fabs ( piii->getX() - pii->getX() ), fabs ( piii->getY() - pii->getY() ) );
                Node3DCartesian <T> n_temp ( pii->getX() - 2 * dx , pii->getY() );

                turning_angle_sum += Bearing::getBearing ( pii, &n_temp );
        }

        //Add turning angle computed for normalized length of initial edge to map
        tf[normalize_length_sum] = turning_angle_sum;

        //Assign points from initial edge
        pi = pii;
        pii = piii;

        //Increment initial edge
        e = e_start->getNextEdge();

        //Process all edges of the Face (2, n-1)
        do
        {
                //Get third point
                piii = e->getNextEdge()->getPoint();

                //Segment is too short
                if ( EuclDistance::getEuclDistance2D ( pii->getX(), pii->getY(), piii->getX(), piii->getY() ) > MIN_POSITION_DIFF )
                {
                        //Compute angle: rotation invariant
                        turning_angle_sum += 180.0 - Angle3Points::getAngle3Points ( pi, pii, piii );

                        //Compute normalized length sum
                        normalize_length_sum += ( scale_method == ScaleInvariant ? EuclDistance::getEuclDistance2D ( pi->getX(), pi->getY(), pii->getX(), pii->getY() ) / face_perimeter :
                                                  EuclDistance::getEuclDistance2D ( pi->getX(), pi->getY(), pii->getX(), pii->getY() ) );

                        //Add turning angle computed for normalized length to map
                        tf[normalize_length_sum] = turning_angle_sum;

                        //Assign points
                        pi = pii;
                        pii = piii;
                }

                //Increment edge
                e = e->getNextEdge();

        }
        while ( e != e_start );

        //Add last value: we visit first point again
        const T last_key = ( scale_method == ScaleInvariant ? 1.0 : face_perimeter );
        tf[last_key] = ( rotation_method == RotationInvariant ? turning_angle_sum + tf.begin()->second : tf.begin()->second );
}