wmlconvexhull3.cpp

3D Game Engine Design Source Code非常棒

// Magic Software, Inc.
// http://www.magic-software.com
// http://www.wild-magic.com
// Copyright (c) 2003.  All Rights Reserved
//
// The Wild Magic Library (WML) source code is supplied under the terms of
// the license agreement http://www.magic-software.com/License/WildMagic.pdf
// and may not be copied or disclosed except in accordance with the terms of
// that agreement.

#include "WmlConvexHull3.h"
using namespace Wml;

#include <algorithm>
using namespace std;

//----------------------------------------------------------------------------
template <class Real>
ConvexHull3<Real>::ConvexHull3 (int iVQuantity, const Vector3<Real>* akVertex)
    :
    m_kHullS(NULL,CreateTriangle)
{
    // uniformly scale and translate to [-1,1]^3 for numerical preconditioning
    Real fMin = akVertex[0].X(), fMax = fMin;
    int i;
    for (i = 0; i < iVQuantity; i++)
    {
        if ( akVertex[i].X() < fMin )
            fMin = akVertex[i].X();
        else if ( akVertex[i].X() > fMax )
            fMax = akVertex[i].X();

        if ( akVertex[i].Y() < fMin )
            fMin = akVertex[i].Y();
        else if ( akVertex[i].Y() > fMax )
            fMax = akVertex[i].Y();

        if ( akVertex[i].Z() < fMin )
            fMin = akVertex[i].Z();
        else if ( akVertex[i].Z() > fMax )
            fMax = akVertex[i].Z();
    }
    Real fHalfRange = ((Real)0.5)*(fMax - fMin);
    Real fInvHalfRange = ((Real)1.0)/fHalfRange;

    m_iVQuantity = iVQuantity;
    m_akVertex = new Vector3<Real>[m_iVQuantity];
    for (i = 0; i < m_iVQuantity; i++)
    {
        m_akVertex[i].X() = -(Real)1.0+fInvHalfRange*(akVertex[i].X()-fMin);
        m_akVertex[i].Y() = -(Real)1.0+fInvHalfRange*(akVertex[i].Y()-fMin);
        m_akVertex[i].Z() = -(Real)1.0+fInvHalfRange*(akVertex[i].Z()-fMin);
    }

    // Compute convex hull incrementally.  The first and second vertices in
    // the hull are managed separately until at least one triangle is formed.
    // At that point, a triangle mesh is used to store the hull.
    m_iHullType = HULL_POINT;
    m_kHullP.push_back(0);
    for (i = 1; i < m_iVQuantity; i++)
    {
        switch ( m_iHullType )
        {
        case HULL_POINT:
            MergePoint(i);
            break;
        case HULL_LINEAR:
            MergeLinear(i);
            break;
        case HULL_PLANAR:
            MergePlanar(i);
            break;
        case HULL_SPATIAL:
            MergeSpatial(i);
            break;
        }
    }

    // map convex hull indices back to original ordering
    if ( m_iHullType == HULL_SPATIAL )
    {
        const ETManifoldMesh::TMap& rkTMap = m_kHullS.GetTriangles();
        m_kHull.resize(3*rkTMap.size());
        i = 0;
        ETManifoldMesh::TMap::const_iterator pkTIter;
        for (pkTIter = rkTMap.begin(); pkTIter != rkTMap.end(); pkTIter++)
        {
            const Triangle* pkTri = (const Triangle*)pkTIter->second;
            for (int j = 0; j < 3; j++)
                m_kHull[i++] = pkTri->V[j];
        }
    }
    else
    {
        int iHQuantity = (int)m_kHullP.size();
        m_kHull.resize(iHQuantity);

        for (i = 0; i < iHQuantity; i++)
            m_kHull[i] = m_kHullP[i];
    }

    m_kHullP.clear();
    delete[] m_akVertex;
    m_akVertex = NULL;
}
//----------------------------------------------------------------------------
template <class Real>
ConvexHull3<Real>::~ConvexHull3 ()
{
    assert( m_akVertex == NULL );
}
//----------------------------------------------------------------------------
template <class Real>
int ConvexHull3<Real>::GetType () const
{
    return m_iHullType;
}
//----------------------------------------------------------------------------
template <class Real>
const vector<int>& ConvexHull3<Real>::GetConnectivity () const
{
    return m_kHull;
}
//----------------------------------------------------------------------------
template <class Real>
void ConvexHull3<Real>::MergePoint (int iP)
{
    // hull is {Q0}
    Vector3<Real> kDiff = m_akVertex[m_kHullP[0]] - m_akVertex[iP];
    if ( kDiff.Length() > Math<Real>::EPSILON )
    {
        m_iHullType = HULL_LINEAR;
        m_kHullP.push_back(iP);
    }
}
//----------------------------------------------------------------------------
template <class Real>
void ConvexHull3<Real>::MergeLinear (int iP)
{
    // hull is {Q0,Q1}
    const Vector3<Real>& rkQ0 = m_akVertex[m_kHullP[0]];
    const Vector3<Real>& rkQ1 = m_akVertex[m_kHullP[1]];
    const Vector3<Real>& rkP = m_akVertex[iP];

    // assert:  Q0 and Q1 are not collocated
    Vector3<Real> kEdge1 = rkQ1 - rkQ0;
    Vector3<Real> kEdge2 = rkP - rkQ0;
    m_kNormal = kEdge1.Cross(kEdge2);
    Real fLength = m_kNormal.Normalize();

    if ( fLength > Math<Real>::EPSILON )
    {
        // <Q0,Q1,P> is a triangle, updated hull is {Q0,Q1,P}
        m_iHullType = HULL_PLANAR;
        m_kHullP.push_back(iP);
        m_kOrigin = rkQ0;
        return;
    }

    // P is on line of <Q0,Q1>
    Real fE1dE2 = kEdge1.Dot(kEdge2);
    if ( fE1dE2 < 0.0f )
    {
        // order is <P,Q0,Q1>, updated hull is {P,Q1}
        m_kHullP[0] = iP;
        return;
    }

    Real fE1dE1 = kEdge1.SquaredLength();
    if ( fE1dE2 > fE1dE1 )
    {
        // order is <Q0,Q1,P>, updated hull is {Q0,P}
        m_kHullP[1] = iP;
        return;
    }

    // order is <Q0,P,Q1>, hull does not change
}
//----------------------------------------------------------------------------
template <class Real>
void ConvexHull3<Real>::MergePlanar (int iP)
{
    // hull is convex polygon in plane m_kNormal.Dot(X - m_kOrigin) = 0
    const Vector3<Real>& rkP = m_akVertex[iP];
    Real fOrder = m_kNormal.Dot(rkP - m_kOrigin);
    int i;
    if ( Math<Real>::FAbs(fOrder) <= Math<Real>::EPSILON )
    {
        // order essentially zero
        Vector3<Real>* pkQ0;
        Vector3<Real>* pkQ1;
        Vector3<Real> kEdge1, kEdge2;

        // find an edge visible to P
        int iSize = (int)m_kHullP.size(), iSizeM1 = iSize-1;
        int iL, iU;
        for (iL = iSizeM1, iU = 0; iU < iSize; iL = iU++)
        {
            pkQ0 = &m_akVertex[m_kHullP[iL]];
            pkQ1 = &m_akVertex[m_kHullP[iU]];
            kEdge1 = *pkQ1 - *pkQ0;
            kEdge2 = rkP - *pkQ0;
            fOrder = m_kNormal.Dot(kEdge1.UnitCross(kEdge2));
            if ( fOrder < (Real)0.0 )
                break;
        }
        if ( iU == iSize )
        {
            // P is inside the current hull
            return;
        }

        // Edge <L,U> is visible to P.  Perform a breadth-first search to
        // locate all remaining visible edges.

        // find the upper index for visible edges
        for (i = iU+1; i < iSize; iU = i++)
        {
            pkQ0 = &m_akVertex[m_kHullP[iU]];
            pkQ1 = &m_akVertex[m_kHullP[i]];
            kEdge1 = *pkQ1 - *pkQ0;
            kEdge2 = rkP - *pkQ0;
            fOrder = m_kNormal.Dot(kEdge1.UnitCross(kEdge2));
            if ( fOrder >= (Real)0.0 )
                break;
        }

        // If edge <size-1,0> was not visible, we had to search for the first
        // visible edge in which case the lower index L is already correct.
        // Thus, we only need to handle the case when L = size-1 and search
        // for smaller indices to find the lower index.
        if ( iL == iSizeM1 )
        {
            for (i = iL-1; i >= 0; iL = i--)
            {
                pkQ0 = &m_akVertex[m_kHullP[i]];
                pkQ1 = &m_akVertex[m_kHullP[iL]];
                kEdge1 = *pkQ1 - *pkQ0;
                kEdge2 = rkP - *pkQ0;
                fOrder = m_kNormal.Dot(kEdge1.UnitCross(kEdge2));
                if ( fOrder >= (Real)0.0 )
                    break;
            }