mgcminimize1d.cpp

3D Game Engine Design Source Code非常棒

// 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 "MgcMinimize1D.h"
#include "MgcRTLib.h"

//----------------------------------------------------------------------------
MgcMinimize1D::MgcMinimize1D (Function oF, unsigned int uiMaxLevel,
    unsigned int uiMaxBracket, void* pvUserData)
{
    assert( m_oF );
    m_oF = oF;
    m_uiMaxLevel = uiMaxLevel;
    m_uiMaxBracket = uiMaxBracket;
    m_pvUserData = pvUserData;
}
//----------------------------------------------------------------------------
void MgcMinimize1D::GetMinimum (MgcReal fT0, MgcReal fT1, MgcReal fTInitial,
    MgcReal& rfTMin, MgcReal& rfFMin)
{
    assert( fT0 <= fTInitial && fTInitial <= fT1 );

    m_fTMin = MgcMath::INFINITY;
    m_fFMin = MgcMath::INFINITY;

    MgcReal fF0 = m_oF(fT0,m_pvUserData);
    MgcReal fFInitial = m_oF(fTInitial,m_pvUserData);
    MgcReal fF1 = m_oF(fT1,m_pvUserData);

    GetMinimum(fT0,fF0,fTInitial,fFInitial,fT1,fF1,m_uiMaxLevel);

    rfTMin = m_fTMin;
    rfFMin = m_fFMin;
}
//----------------------------------------------------------------------------
void MgcMinimize1D::GetMinimum (MgcReal fT0, MgcReal fF0, MgcReal fTm,
    MgcReal fFm, MgcReal fT1, MgcReal fF1, unsigned int uiLevel)
{
    if ( fF0 < m_fFMin )
    {
        m_fTMin = fT0;
        m_fFMin = fF0;
    }

    if ( fFm < m_fFMin )
    {
        m_fTMin = fTm;
        m_fFMin = fFm;
    }

    if ( fF1 < m_fFMin )
    {
        m_fTMin = fT1;
        m_fFMin = fF1;
    }

    if ( uiLevel-- == 0 )
        return;

    if ( (fT1 - fTm)*(fF0 - fFm) > (fTm - fT0)*(fFm - fF1) )
    {
        // quadratic fit has positive second derivative at midpoint

        if ( fF1 > fF0 )
        {
            if ( fFm >= fF0 )
            {
                // increasing, repeat on [t0,tm]
                GetMinimum(fT0,fF0,fTm,fFm,uiLevel);
            }
            else
            {
                // not monotonic, have a bracket
                GetBracketedMinimum(fT0,fF0,fTm,fFm,fT1,fF1,uiLevel);
            }
        }
        else if ( fF1 < fF0 )
        {
            if ( fFm >= fF1 )
            {
                // decreasing, repeat on [tm,t1]
                GetMinimum(fTm,fFm,fT1,fF1,uiLevel);
            }
            else
            {
                // not monotonic, have a bracket
                GetBracketedMinimum(fT0,fF0,fTm,fFm,fT1,fF1,uiLevel);
            }
        }
        else
        {
            // constant, repeat on [t0,tm] and [tm,t1]
            GetMinimum(fT0,fF0,fTm,fFm,uiLevel);
            GetMinimum(fTm,fFm,fT1,fF1,uiLevel);
        }
    }
    else
    {
        // quadratic fit has nonpositive second derivative at midpoint

        if ( fF1 > fF0 )
        {
            // repeat on [t0,tm]
            GetMinimum(fT0,fF0,fTm,fFm,uiLevel);
        }
        else if ( fF1 < fF0 )
        {
            // repeat on [tm,t1]
            GetMinimum(fTm,fFm,fT1,fF1,uiLevel);
        }
        else
        {
            // repeat on [t0,tm] and [tm,t1]
            GetMinimum(fT0,fF0,fTm,fFm,uiLevel);
            GetMinimum(fTm,fFm,fT1,fF1,uiLevel);
        }
    }
}
//----------------------------------------------------------------------------
void MgcMinimize1D::GetMinimum (MgcReal fT0, MgcReal fF0, MgcReal fT1,
    MgcReal fF1, unsigned int uiLevel)
{
    if ( fF0 < m_fFMin )
    {
        m_fTMin = fT0;
        m_fFMin = fF0;
    }

    if ( fF1 < m_fFMin )
    {
        m_fTMin = fT1;
        m_fFMin = fF1;
    }

    if ( uiLevel-- == 0 )
        return;

    MgcReal fTm = 0.5*(fT0+fT1);
    MgcReal fFm = m_oF(fTm,m_pvUserData);

    if ( fF0 - 2.0*fFm + fF1 > 0.0 )
    {
        // quadratic fit has positive second derivative at midpoint

        if ( fF1 > fF0 )
        {
            if ( fFm >= fF0 )
            {
                // increasing, repeat on [t0,tm]
                GetMinimum(fT0,fF0,fTm,fFm,uiLevel);
            }
            else
            {
                // not monotonic, have a bracket
                GetBracketedMinimum(fT0,fF0,fTm,fFm,fT1,fF1,uiLevel);
            }
        }
        else if ( fF1 < fF0 )
        {
            if ( fFm >= fF1 )
            {
                // decreasing, repeat on [tm,t1]
                GetMinimum(fTm,fFm,fT1,fF1,uiLevel);
            }
            else
            {
                // not monotonic, have a bracket
                GetBracketedMinimum(fT0,fF0,fTm,fFm,fT1,fF1,uiLevel);
            }
        }
        else
        {
            // constant, repeat on [t0,tm] and [tm,t1]
            GetMinimum(fT0,fF0,fTm,fFm,uiLevel);
            GetMinimum(fTm,fFm,fT1,fF1,uiLevel);
        }
    }
    else
    {
        // quadratic fit has nonpositive second derivative at midpoint

        if ( fF1 > fF0 )
        {
            // repeat on [t0,tm]
            GetMinimum(fT0,fF0,fTm,fFm,uiLevel);
        }
        else if ( fF1 < fF0 )
        {
            // repeat on [tm,t1]
            GetMinimum(fTm,fFm,fT1,fF1,uiLevel);
        }
        else
        {
            // repeat on [t0,tm] and [tm,t1]
            GetMinimum(fT0,fF0,fTm,fFm,uiLevel);
            GetMinimum(fTm,fFm,fT1,fF1,uiLevel);
        }
    }
}
//----------------------------------------------------------------------------
void MgcMinimize1D::GetBracketedMinimum (MgcReal fT0, MgcReal fF0,
    MgcReal fTm, MgcReal fFm, MgcReal fT1, MgcReal fF1, unsigned int uiLevel)
{
    for (unsigned int uiI = 0; uiI < m_uiMaxBracket; uiI++)
    {
        // update minimum value
        if ( fFm < m_fFMin )
        {
            m_fTMin = fTm;
            m_fFMin = fFm;
        }

        // test for convergence
        const MgcReal fEps = 1e-08, fTol = 1e-04;
        if ( MgcMath::Abs(fT1-fT0) <= 2.0*fTol*MgcMath::Abs(fTm) + fEps )
            break;

        // compute vertex of interpolating parabola
        MgcReal fDT0 = fT0 - fTm, fDT1 = fT1 - fTm;
        MgcReal fDF0 = fF0 - fFm, fDF1 = fF1 - fFm;
        MgcReal fTmp0 = fDT0*fDF1, fTmp1 = fDT1*fDF0;
        MgcReal fDenom = fTmp1 - fTmp0;
        if ( MgcMath::Abs(fDenom) < fEps )
           return;

        MgcReal fTv = fTm + 0.5*(fDT1*fTmp1-fDT0*fTmp0)/fDenom;
        assert( fT0 <= fTv && fTv <= fT1 );
        MgcReal fFv = m_oF(fTv,m_pvUserData);

        if ( fTv < fTm )
        {
            if ( fFv < fFm )
            {
                fT1 = fTm;
                fF1 = fFm;
                fTm = fTv;
                fFm = fFv;
            }
            else
            {
                fT0 = fTv;
                fF0 = fFv;
            }
        }
        else if ( fTv > fTm )
        {
            if ( fFv < fFm )
            {
                fT0 = fTm;
                fF0 = fFm;
                fTm = fTv;
                fFm = fFv;
            }
            else
            {
                fT1 = fTv;
                fF1 = fFv;
            }
        }
        else
        {
            // vertex of parabola is already at middle sample point
            GetMinimum(fT0,fF0,fTm,fFm,uiLevel);
            GetMinimum(fTm,fFm,fT1,fF1,uiLevel);
        }
    }
}
//----------------------------------------------------------------------------