wmlextractsurfacetetra.cpp

来自「3D Game Engine Design Source Code非常棒」· C++ 代码 · 共 1,297 行 · 第 1/3 页

    for (int iV = 0, iNextVertex = 0; iV < iVQuantity; iV++)
    {
        // keep only unique vertices
        pair<VMapIterator,bool> kResult = kVMap.insert(
            make_pair(rkVA[iV],iNextVertex));
        if ( kResult.second == true )
            iNextVertex++;
    }

    // use a hash table to generate unique storage
    TMap kTMap;
    TMapIterator pkTIter;
    for (int iT = 0, iNextTriangle = 0; iT < iTQuantity; iT++)
    {
        // replace old vertex indices by new ones
        TriangleKey& rkTri = rkTA[iT];
        pkVIter = kVMap.find(rkVA[rkTri.V[0]]);
        assert( pkVIter != kVMap.end() );
        rkTri.V[0] = pkVIter->second;
        pkVIter = kVMap.find(rkVA[rkTri.V[1]]);
        assert( pkVIter != kVMap.end() );
        rkTri.V[1] = pkVIter->second;
        pkVIter = kVMap.find(rkVA[rkTri.V[2]]);
        assert( pkVIter != kVMap.end() );
        rkTri.V[2] = pkVIter->second;

        // keep only unique triangles
        pair<TMapIterator,bool> kResult = kTMap.insert(
            make_pair(rkTri,iNextTriangle));
        if ( kResult.second == true )
            iNextTriangle++;
    }

    // pack the vertices
    rkVA.resize(kVMap.size());
    for (pkVIter = kVMap.begin(); pkVIter != kVMap.end(); pkVIter++)
        rkVA[pkVIter->second] = pkVIter->first;

    // pack the triangles
    rkTA.resize(kTMap.size());
    for (pkTIter = kTMap.begin(); pkTIter != kTMap.end(); pkTIter++)
        rkTA[pkTIter->second] = pkTIter->first;
}
//----------------------------------------------------------------------------
void ExtractSurfaceTetra::OrientTriangles (vector<Vector3f>& rkVA,
    vector<TriangleKey>& rkTA, bool bSameDir)
{
    for (int i = 0; i < (int)rkTA.size(); i++)
    {
        TriangleKey& rkTri = rkTA[i];

        // get triangle vertices
        Vector3f kV0 = rkVA[rkTri.V[0]];
        Vector3f kV1 = rkVA[rkTri.V[1]];
        Vector3f kV2 = rkVA[rkTri.V[2]];
        
        // construct triangle normal based on current orientation
        Vector3f kEdge1 = kV1 - kV0;
        Vector3f kEdge2 = kV2 - kV0;
        Vector3f kNormal = kEdge1.Cross(kEdge2);

        // get the image gradient at the vertices
        Vector3f kGrad0 = GetGradient(kV0);
        Vector3f kGrad1 = GetGradient(kV1);
        Vector3f kGrad2 = GetGradient(kV2);

        // compute the average gradient
        Vector3f kGradAvr = (kGrad0 + kGrad1 + kGrad2)/3.0f;
        
        // compute the dot product of normal and average gradient
        float fDot = kGradAvr.Dot(kNormal);

        // choose triangle orientation based on gradient direction
        int iSave;
        if ( bSameDir )
        {
            if ( fDot < 0.0f )
            {
                // wrong orientation, reorder it
                iSave = rkTri.V[1];
                rkTri.V[1] = rkTri.V[2];
                rkTri.V[2] = iSave;
            }
        }
        else
        {
            if ( fDot > 0.0f )
            {
                // wrong orientation, reorder it
                iSave = rkTri.V[1];
                rkTri.V[1] = rkTri.V[2];
                rkTri.V[2] = iSave;
            }
        }
    }
}
//----------------------------------------------------------------------------
void ExtractSurfaceTetra::ComputeNormals (const vector<Vector3f>& rkVA,
    const vector<TriangleKey>& rkTA, vector<Vector3f>& rkNA)
{
    // maintain a running sum of triangle normals at each vertex
    int iVQuantity = (int)rkVA.size();
    int iTQuantity = (int)rkTA.size();
    rkNA.resize(iVQuantity);
    int i, j;
    for (i = 0; i < iVQuantity; i++)
        rkNA[i] = Vector3f::ZERO;

    for (i = 0, j = 0; i < iTQuantity; i++)
    {
        const TriangleKey& rkT = rkTA[i];
        Vector3f kV0 = rkVA[rkT.V[0]];
        Vector3f kV1 = rkVA[rkT.V[1]];
        Vector3f kV2 = rkVA[rkT.V[2]];

        // construct triangle normal
        Vector3f kEdge1 = kV1 - kV0;
        Vector3f kEdge2 = kV2 - kV0;
        Vector3f kNormal = kEdge1.Cross(kEdge2);

        // maintain the sum of normals at each vertex
        rkNA[rkT.V[0]] += kNormal;
        rkNA[rkT.V[1]] += kNormal;
        rkNA[rkT.V[2]] += kNormal;
    }

    // The normal vector storage was used to accumulate the sum of
    // triangle normals.  Now these vectors must be rescaled to be
    // unit length.
    for (i = 0; i < iVQuantity; i++)
        rkNA[i].Normalize();
}
//----------------------------------------------------------------------------
int ExtractSurfaceTetra::AddVertex (VtxMap& rkVMap, int iXNumer, int iXDenom,
    int iYNumer, int iYDenom, int iZNumer, int iZDenom)
{
    Vertex kVertex(iXNumer,iXDenom,iYNumer,iYDenom,iZNumer,iZDenom);
    VtxMapIterator pkVIter = rkVMap.find(kVertex);
    if ( pkVIter != rkVMap.end() )
    {
        // Vertex already in map, just return its unique index.
        return pkVIter->second;
    }
    else
    {
        // Vertex not in map, insert it and assign it a unique index.
        int i = m_iNextIndex++;
        rkVMap.insert(make_pair(kVertex,i));
        return i;
    }
}
//----------------------------------------------------------------------------
void ExtractSurfaceTetra::AddEdge (VtxMap& rkVMap, ESet& rkESet, int iXNumer0,
    int iXDenom0, int iYNumer0, int iYDenom0, int iZNumer0, int iZDenom0,
    int iXNumer1, int iXDenom1, int iYNumer1, int iYDenom1, int iZNumer1,
    int iZDenom1)
{
    int iV0 = AddVertex(rkVMap,iXNumer0,iXDenom0,iYNumer0,iYDenom0,iZNumer0,
        iZDenom0);

    int iV1 = AddVertex(rkVMap,iXNumer1,iXDenom1,iYNumer1,iYDenom1,iZNumer1,
        iZDenom1);

    rkESet.insert(EdgeKey(iV0,iV1));
}
//----------------------------------------------------------------------------
void ExtractSurfaceTetra::AddTriangle (VtxMap& rkVMap, ESet& rkESet,
    TSet& rkTSet, int iXNumer0, int iXDenom0, int iYNumer0, int iYDenom0,
    int iZNumer0, int iZDenom0, int iXNumer1, int iXDenom1, int iYNumer1,
    int iYDenom1, int iZNumer1, int iZDenom1, int iXNumer2, int iXDenom2,
    int iYNumer2, int iYDenom2, int iZNumer2, int iZDenom2)
{
    int iV0 = AddVertex(rkVMap,iXNumer0,iXDenom0,iYNumer0,iYDenom0,iZNumer0,
        iZDenom0);

    int iV1 = AddVertex(rkVMap,iXNumer1,iXDenom1,iYNumer1,iYDenom1,iZNumer1,
        iZDenom1);

    int iV2 = AddVertex(rkVMap,iXNumer2,iXDenom2,iYNumer2,iYDenom2,iZNumer2,
        iZDenom2);

    // nothing to do if triangle already exists
    TriangleKey kT(iV0,iV1,iV2);
    if ( rkTSet.find(kT) != rkTSet.end() )
        return;

    // prevent double-sided triangles
    int iSave = kT.V[1];
    kT.V[1] = kT.V[2];
    kT.V[2] = iSave;
    if ( rkTSet.find(kT) != rkTSet.end() )
        return;

    rkESet.insert(EdgeKey(iV0,iV1));
    rkESet.insert(EdgeKey(iV1,iV2));
    rkESet.insert(EdgeKey(iV2,iV0));

    // compute triangle normal assuming counterclockwise ordering
    Vector3f kV0(
        iXNumer0/(float)iXDenom0,
        iYNumer0/(float)iYDenom0,
        iZNumer0/(float)iZDenom0);

    Vector3f kV1(
        iXNumer1/(float)iXDenom1,
        iYNumer1/(float)iYDenom1,
        iZNumer1/(float)iZDenom1);

    Vector3f kV2(
        iXNumer2/(float)iXDenom2,
        iYNumer2/(float)iYDenom2,
        iZNumer2/(float)iZDenom2);

    Vector3f kE0 = kV1 - kV0;
    Vector3f kE1 = kV2 - kV0;
    Vector3f kN = kE0.Cross(kE1);

    // choose triangle orientation based on gradient direction
    Vector3f kCentroid = (kV0+kV1+kV2)/3.0f;
    Vector3f kGrad = GetGradient(kCentroid);
    if ( kGrad.Dot(kN) <= 0.0f )
        rkTSet.insert(TriangleKey(iV0,iV1,iV2));
    else
        rkTSet.insert(TriangleKey(iV0,iV2,iV1));
}
//----------------------------------------------------------------------------
void ExtractSurfaceTetra::ProcessTetrahedron (VtxMap& rkVM,
    ESet& rkES, TSet& rkTS, int iX0, int iY0, int iZ0, int iF0,
    int iX1, int iY1, int iZ1, int iF1, int iX2, int iY2, int iZ2, int iF2,
    int iX3, int iY3, int iZ3, int iF3)
{
    int iXN0, iYN0, iZN0, iD0;
    int iXN1, iYN1, iZN1, iD1;
    int iXN2, iYN2, iZN2, iD2;
    int iXN3, iYN3, iZN3, iD3;

    if ( iF0 != 0 )
    {
        // convert to case +***
        if ( iF0 < 0 )
        {
            iF0 = -iF0;
            iF1 = -iF1;
            iF2 = -iF2;
            iF3 = -iF3;
        }

        if ( iF1 > 0 )
        {
            if ( iF2 > 0 )
            {
                if ( iF3 > 0 )
                {
                    // ++++
                    return;
                }
                else if ( iF3 < 0 )
                {
                    // +++-
                    iD0 = iF0 - iF3;
                    iXN0 = iF0*iX3 - iF3*iX0;
                    iYN0 = iF0*iY3 - iF3*iY0;
                    iZN0 = iF0*iZ3 - iF3*iZ0;
                    iD1 = iF1 - iF3;
                    iXN1 = iF1*iX3 - iF3*iX1;
                    iYN1 = iF1*iY3 - iF3*iY1;
                    iZN1 = iF1*iZ3 - iF3*iZ1;
                    iD2 = iF2 - iF3;
                    iXN2 = iF2*iX3 - iF3*iX2;
                    iYN2 = iF2*iY3 - iF3*iY2;
                    iZN2 = iF2*iZ3 - iF3*iZ2;
                    AddTriangle(rkVM,rkES,rkTS,
                        iXN0,iD0,iYN0,iD0,iZN0,iD0,
                        iXN1,iD1,iYN1,iD1,iZN1,iD1,
                        iXN2,iD2,iYN2,iD2,iZN2,iD2);
                }
                else
                {
                    // +++0
                    AddVertex(rkVM,iX3,1,iY3,1,iZ3,1);
                }
            }
            else if ( iF2 < 0 )
            {
                iD0 = iF0 - iF2;
                iXN0 = iF0*iX2 - iF2*iX0;
                iYN0 = iF0*iY2 - iF2*iY0;
                iZN0 = iF0*iZ2 - iF2*iZ0;
                iD1 = iF1 - iF2;
                iXN1 = iF1*iX2 - iF2*iX1;
                iYN1 = iF1*iY2 - iF2*iY1;
                iZN1 = iF1*iZ2 - iF2*iZ1;

                if ( iF3 > 0 )
                {
                    // ++-+
                    iD2 = iF3 - iF2;
                    iXN2 = iF3*iX2 - iF2*iX3;
                    iYN2 = iF3*iY2 - iF2*iY3;
                    iZN2 = iF3*iZ2 - iF2*iZ3;
                    AddTriangle(rkVM,rkES,rkTS,
                        iXN0,iD0,iYN0,iD0,iZN0,iD0,
                        iXN1,iD1,iYN1,iD1,iZN1,iD1,
                        iXN2,iD2,iYN2,iD2,iZN2,iD2);
                }
                else if ( iF3 < 0 )
                {
                    // ++--
                    iD2 = iF0 - iF3;
                    iXN2 = iF0*iX3 - iF3*iX0;
                    iYN2 = iF0*iY3 - iF3*iY0;
                    iZN2 = iF0*iZ3 - iF3*iZ0;
                    iD3 = iF1 - iF3;
                    iXN3 = iF1*iX3 - iF3*iX1;
                    iYN3 = iF1*iY3 - iF3*iY1;
                    iZN3 = iF1*iZ3 - iF3*iZ1;
                    AddTriangle(rkVM,rkES,rkTS,
                        iXN0,iD0,iYN0,iD0,iZN0,iD0,
                        iXN1,iD1,iYN1,iD1,iZN1,iD1,
                        iXN2,iD2,iYN2,iD2,iZN2,iD2);
                    AddTriangle(rkVM,rkES,rkTS,
                        iXN1,iD1,iYN1,iD1,iZN1,iD1,
                        iXN3,iD3,iYN3,iD3,iZN3,iD3,
                        iXN2,iD2,iYN2,iD2,iZN2,iD2);
                }
                else
                {
                    // ++-0
                    AddTriangle(rkVM,rkES,rkTS,
                        iXN0,iD0,iYN0,iD0,iZN0,iD0,
                        iXN1,iD1,iYN1,iD1,iZN1,iD1,
                        iX3,1,iY3,1,iZ3,1);
                }
            }
            else
            {
                if ( iF3 > 0 )
                {
                    // ++0+
                    AddVertex(rkVM,iX2,1,iY2,1,iZ2,1);
                }
                else if ( iF3 < 0 )
                {
                    // ++0-
                    iD0 = iF0 - iF3;
                    iXN0 = iF0*iX3 - iF3*iX0;
                    iYN0 = iF0*iY3 - iF3*iY0;
                    iZN0 = iF0*iZ3 - iF3*iZ0;
                    iD1 = iF1 - iF3;
                    iXN1 = iF1*iX3 - iF3*iX1;
                    iYN1 = iF1*iY3 - iF3*iY1;
                    iZN1 = iF1*iZ3 - iF3*iZ1;
                    AddTriangle(rkVM,rkES,rkTS,
                        iXN0,iD0,iYN0,iD0,iZN0,iD0,
                        iXN1,iD1,iYN1,iD1,iZN1,iD1,
                        iX2,1,iY2,1,iZ2,1);
                }
                else
                {
                    // ++00
                    AddEdge(rkVM,rkES,iX2,1,iY2,1,iZ2,1,iX3,1,iY3,1,iZ3,1);
                }
            }
        }
        else if ( iF1 < 0 )
        {
            if ( iF2 > 0 )
            {
                iD0 = iF0 - iF1;
                iXN0 = iF0*iX1 - iF1*iX0;
                iYN0 = iF0*iY1 - iF1*iY0;
                iZN0 = iF0*iZ1 - iF1*iZ0;
                iD1 = iF2 - iF1;
                iXN1 = iF2*iX1 - iF1*iX2;
                iYN1 = iF2*iY1 - iF1*iY2;
                iZN1 = iF2*iZ1 - iF1*iZ2;

                if ( iF3 > 0 )
                {
                    // +-++
                    iD2 = iF3 - iF1;
                    iXN2 = iF3*iX1 - iF1*iX3;
                    iYN2 = iF3*iY1 - iF1*iY3;
                    iZN2 = iF3*iZ1 - iF1*iZ3;
                    AddTriangle(rkVM,rkES,rkTS,
                        iXN0,iD0,iYN0,iD0,iZN0,iD0,
                        iXN1,iD1,iYN1,iD1,iZN1,iD1,
                        iXN2,iD2,iYN2,iD2,iZN2,iD2);
                }
                else if ( iF3 < 0 )
                {
                    // +-+-
                    iD2 = iF0 - iF3;
                    iXN2 = iF0*iX3 - iF3*iX0;
                    iYN2 = iF0*iY3 - iF3*iY0;
                    iZN2 = iF0*iZ3 - iF3*iZ0;
                    iD3 = iF2 - iF3;
                    iXN3 = iF2*iX3 - iF3*iX2;
                    iYN3 = iF2*iY3 - iF3*iY2;
                    iZN3 = iF2*iZ3 - iF3*iZ2;
                    AddTriangle(rkVM,rkES,rkTS,
                        iXN0,iD0,iYN0,iD0,iZN0,iD0,
                        iXN1,iD1,iYN1,iD1,iZN1,iD1,
                        iXN2,iD2,iYN2,iD2,iZN2,iD2);
                    AddTriangle(rkVM,rkES,rkTS,
                        iXN1,iD1,iYN1,iD1,iZN1,iD1,
                        iXN3,iD3,iYN3,iD3,iZN3,iD3,
                        iXN2,iD2,iYN2,iD2,iZN2,iD2);
                }
                else
                {
                    // +-+0
                    AddTriangle(rkVM,rkES,rkTS,
                        iXN0,iD0,iYN0,iD0,iZN0,iD0,
                        iXN1,iD1,iYN1,iD1,iZN1,iD1,
                        iX3,1,iY3,1,iZ3,1);
                }
            }
            else if ( iF2 < 0 )
            {
                iD0 = iF1 - iF0;
                iXN0 = iF1*iX0 - iF0*iX1;
                iYN0 = iF1*iY0 - iF0*iY1;
                iZN0 = iF1*iZ0 - iF0*iZ1;
                iD1 = iF2 - iF0;
                iXN1 = iF2*iX0 - iF0*iX2;
                iYN1 = iF2*iY0 - iF0*iY2;
                iZN1 = iF2*iZ0 - iF0*iZ2;

                if ( iF3 > 0 )
                {
                    // +--+
                    iD2 = iF1 - iF3;