mgcquadraticnetwork.cpp

《3D游戏引擎设计》的源码

// Magic Software, Inc.
// http://www.magic-software.com
// Copyright (c) 2000, All Rights Reserved
//
// Source code from Magic Software is supplied under the terms of a license
// agreement and may not be copied or disclosed except in accordance with the
// terms of that agreement.  The various license agreements may be found at
// the Magic Software web site.  This file is subject to the license
//
// FREE SOURCE CODE
// http://www.magic-software.com/License/free.pdf

#include "MgcQuadraticNetwork.h"

//----------------------------------------------------------------------------
MgcQuadraticNetwork::MgcQuadraticNetwork (int iVertexQuantity,
    MgcVector2* akVertex, MgcReal* afF, MgcReal* afFx, MgcReal* afFy)
    :
    MgcTriangleNetwork(iVertexQuantity,akVertex)
{
    assert( afF && afFx && afFy );
    m_afF = afF;
    m_afFx = afFx;
    m_afFy = afFy;
    ProcessTriangles();
}
//----------------------------------------------------------------------------
MgcQuadraticNetwork::MgcQuadraticNetwork (int iVertexQuantity,
    MgcVector2* akVertex, MgcReal* afF)
    :
    MgcTriangleNetwork(iVertexQuantity,akVertex)
{
    assert( afF );
    m_afF = afF;
    EstimateDerivatives();
    ProcessTriangles();
}
//----------------------------------------------------------------------------
MgcQuadraticNetwork::MgcQuadraticNetwork (MgcTriangleNetwork& rkNet,
    MgcReal* afF, MgcReal* afFx, MgcReal* afFy)
    :
    MgcTriangleNetwork(rkNet)
{
    assert( afF && afFx && afFy );
    m_afF = afF;
    m_afFx = afFx;
    m_afFy = afFy;
    ProcessTriangles();
}
//----------------------------------------------------------------------------
MgcQuadraticNetwork::MgcQuadraticNetwork (MgcTriangleNetwork& rkNet,
    MgcReal* afF)
    :
    MgcTriangleNetwork(rkNet)
{
    assert( afF );
    m_afF = afF;
    EstimateDerivatives();
    ProcessTriangles();
}
//----------------------------------------------------------------------------
MgcQuadraticNetwork::~MgcQuadraticNetwork ()
{
    delete[] m_afF;
    delete[] m_afFx;
    delete[] m_afFy;
    delete[] m_akTData;
    delete[] m_akECenter;
}
//----------------------------------------------------------------------------
void MgcQuadraticNetwork::EstimateDerivatives ()
{
    m_afFx = new MgcReal[m_iVertexQuantity];
    m_afFy = new MgcReal[m_iVertexQuantity];
    MgcReal* afFz = new MgcReal[m_iVertexQuantity];
    memset(m_afFx,0,m_iVertexQuantity*sizeof(MgcReal));
    memset(m_afFy,0,m_iVertexQuantity*sizeof(MgcReal));
    memset(afFz,0,m_iVertexQuantity*sizeof(MgcReal));

    // accumulate normals at spatial locations (averaging process)
    int i;
    for (i = 0; i < m_iTriangleQuantity; i++)
    {
        Triangle& rkTri = m_akTriangle[i];

        // get three vertices of triangle
        int j0 = rkTri.m_aiVertex[0];
        int j1 = rkTri.m_aiVertex[1];
        int j2 = rkTri.m_aiVertex[2];

        // compute normal vector of triangle (with positive z-component)
        MgcReal fDx1 = m_akVertex[j1].x - m_akVertex[j0].x;
        MgcReal fDy1 = m_akVertex[j1].y - m_akVertex[j0].y;
        MgcReal fDz1 = m_afF[j1] - m_afF[j0];
        MgcReal fDx2 = m_akVertex[j2].x - m_akVertex[j0].x;
        MgcReal fDy2 = m_akVertex[j2].y - m_akVertex[j0].y;
        MgcReal fDz2 = m_afF[j2] - m_afF[j0];
        MgcReal fNx = fDy1*fDz2 - fDy2*fDz1;
        MgcReal fNy = fDz1*fDx2 - fDz2*fDx1;
        MgcReal fNz = fDx1*fDy2 - fDx2*fDy1;
        if ( fNz < 0.0 )
        {
            fNx = -fNx;
            fNy = -fNy;
            fNz = -fNz;
        }

        m_afFx[j0] += fNx;  m_afFy[j0] += fNy;  afFz[j0] += fNz;
        m_afFx[j1] += fNx;  m_afFy[j1] += fNy;  afFz[j1] += fNz;
        m_afFx[j2] += fNx;  m_afFy[j2] += fNy;  afFz[j2] += fNz;
    }

    // scale the normals to form (x,y,-1)
    for (i = 0; i < m_iVertexQuantity; i++)
    {
        if ( MgcMath::Abs(afFz[i]) > 1e-06 )
        {
            MgcReal fInv = -1.0/afFz[i];
            m_afFx[i] *= fInv;
            m_afFy[i] *= fInv;
        }
        else
        {
            m_afFx[i] = 0.0;
            m_afFy[i] = 0.0;
        }
    }

    delete[] afFz;
}
//----------------------------------------------------------------------------
void MgcQuadraticNetwork::ProcessTriangles ()
{
    // Add degenerate triangles to boundary triangles so that interpolation
    // at the boundary can be treated in the same way as interpolation in
    // the interior.

    // add quadratic data to triangle network
    m_akTData = new TriangleData[m_iTriangleQuantity];
    m_akECenter = new MgcVector2[m_iExtraTriangleQuantity];

    // compute centers of inscribed circles for triangles
    int iT;
    for (iT = 0; iT < m_iTriangleQuantity; iT++)
    {
        Triangle& rkTri = m_akTriangle[iT];
        MgcTriangleNetwork::ComputeInscribedCenter(
            m_akVertex[rkTri.m_aiVertex[0]],m_akVertex[rkTri.m_aiVertex[1]],
            m_akVertex[rkTri.m_aiVertex[2]],m_akTData[iT].m_kCenter);
    }

    // compute centers of edges on convex hull
    int iE = 0;
    for (iT = 0; iT < m_iTriangleQuantity; iT++)
    {
        Triangle& rkTri = m_akTriangle[iT];
        for (int j = 0; j < 3; j++)
        {
            if ( rkTri.m_aiAdjacent[j] >= m_iTriangleQuantity )
            {
                Triangle& rkETri = m_akExtraTriangle[iE];
                m_akECenter[iE] = 0.5*(
                    m_akVertex[rkTri.m_aiVertex[j]] +
                    m_akVertex[rkTri.m_aiVertex[(j+1) % 3]]);
                iE++;
            }
        }
    }

    // compute cross-edge intersections
    for (iT = 0; iT < m_iTriangleQuantity; iT++)
        ComputeCrossEdgeIntersections(iT);

    // compute Bezier coefficients
    for (iT = 0; iT < m_iTriangleQuantity; iT++)
        ComputeCoefficients(iT);
}
//----------------------------------------------------------------------------
void MgcQuadraticNetwork::ComputeCrossEdgeIntersections (int iT)
{
    Triangle& rkTri = m_akTriangle[iT];
    const MgcVector2& rkV0 = m_akVertex[rkTri.m_aiVertex[0]];
    const MgcVector2& rkV1 = m_akVertex[rkTri.m_aiVertex[1]];
    const MgcVector2& rkV2 = m_akVertex[rkTri.m_aiVertex[2]];

    MgcVector2 akU[3];
    for (int i = 0; i < 3; i++)
    {
        int iA = rkTri.m_aiAdjacent[i];
        if ( iA < m_iTriangleQuantity )
            akU[i] = m_akTData[iA].m_kCenter;
        else
            akU[i] = m_akECenter[iA - m_iTriangleQuantity];
    }

    MgcReal fM00, fM01, fM10, fM11, fR0, fR1, fInvDet;

    // intersection on edge <V0,V1>
    fM00 = rkV0.y - rkV1.y;
    fM01 = rkV1.x - rkV0.x;
    fM10 = m_akTData[iT].m_kCenter.y - akU[0].y;
    fM11 = akU[0].x - m_akTData[iT].m_kCenter.x;
    fR0  = fM00*rkV0.x + fM01*rkV0.y;
    fR1  = fM10*m_akTData[iT].m_kCenter.x + fM11*m_akTData[iT].m_kCenter.y;
    fInvDet = 1.0/(fM00*fM11 - fM01*fM10);
    m_akTData[iT].m_akIntersect[0].x = (fM11*fR0-fM01*fR1)*fInvDet;
    m_akTData[iT].m_akIntersect[0].y = (fM00*fR1-fM10*fR0)*fInvDet;

    // intersection on edge <V1,V2>
    fM00 = rkV1.y - rkV2.y;
    fM01 = rkV2.x - rkV1.x;
    fM10 = m_akTData[iT].m_kCenter.y - akU[1].y;
    fM11 = akU[1].x - m_akTData[iT].m_kCenter.x;
    fR0  = fM00*rkV1.x + fM01*rkV1.y;
    fR1  = fM10*m_akTData[iT].m_kCenter.x + fM11*m_akTData[iT].m_kCenter.y;
    fInvDet = 1.0/(fM00*fM11 - fM01*fM10);
    m_akTData[iT].m_akIntersect[1].x = (fM11*fR0-fM01*fR1)*fInvDet;
    m_akTData[iT].m_akIntersect[1].y = (fM00*fR1-fM10*fR0)*fInvDet;

    // intersection on edge <V0,V2>
    fM00 = rkV0.y - rkV2.y;
    fM01 = rkV2.x - rkV0.x;
    fM10 = m_akTData[iT].m_kCenter.y - akU[2].y;
    fM11 = akU[2].x - m_akTData[iT].m_kCenter.x;
    fR0  = fM00*rkV0.x + fM01*rkV0.y;