📄 mgcdistvec2qdr2.cpp
字号:
// 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 "MgcEigen.h"
#include "MgcDistVec2Qdr2.h"
#include "MgcPolynomial.h"
#include "MgcRTLib.h"
//----------------------------------------------------------------------------
static void ComputeCoeff2 (MgcReal afA[2], MgcReal afB[2], MgcReal afD[2],
MgcReal fC, MgcPolynomial& rkPoly)
{
MgcReal afBPad[2] = { afB[0]+afA[0]*afD[0], afB[1]+afA[1]*afD[1] };
MgcReal afBSqr[2] = { afB[0]*afB[0], afB[1]*afB[1] };
MgcReal afDSqr[2] = { afD[0]*afD[0], afD[1]*afD[1] };
MgcReal fDPrd = afD[0]*afD[1];
MgcReal fDSum = afD[0]+afD[1];
rkPoly[0] = afA[0]*afBPad[0]+afA[1]*afBPad[1]+fC;
rkPoly[1] = -afBSqr[0]-afBSqr[1]+4.0*(afA[0]*afD[1]*afBPad[0]+
afA[1]*afD[0]*afBPad[1]+fC*fDSum);
rkPoly[2] = -afBSqr[0]*(afD[0]+4.0*afD[1])-afBSqr[1]*(afD[1]+4.0*afD[0])+
4.0*(afA[0]*afDSqr[1]*afBPad[0]+afA[1]*afDSqr[0]*afBPad[1]+
fC*(afDSqr[0]+afDSqr[1]+4.0*fDPrd));
MgcReal fTmp = -4.0*(afBSqr[0]*afD[1]+afBSqr[1]*afD[0]-4.0*fC*fDPrd);
rkPoly[3] = fDSum*fTmp;
rkPoly[4] = fDPrd*fTmp;
}
//----------------------------------------------------------------------------
MgcReal MgcSqrDistance (const MgcVector2& rkPoint, const MgcReal afQuad[6],
MgcVector2& rkClosest)
{
// eigendecomposition
MgcEigen kES(2);
kES.Matrix(0,0) = afQuad[3];
kES.Matrix(0,1) = 0.5*afQuad[5];
kES.Matrix(1,0) = kES.Matrix(0,1);
kES.Matrix(1,1) = afQuad[4];
kES.IncrSortEigenStuff2();
MgcReal** aafR = kES.GetEigenvector();
MgcReal afA[2], afB[2], afD[2], fC = afQuad[0];
int i;
for (i = 0; i < 2; i++)
{
afA[i] = aafR[0][i]*rkPoint.x + aafR[1][i]*rkPoint.y;
afB[i] = aafR[0][i]*afQuad[1] + aafR[1][i]*afQuad[2];
afD[i] = kES.GetEigenvalue(i);
}
MgcPolynomial kPoly(4);
ComputeCoeff2(afA,afB,afD,fC,kPoly);
int iCount;
MgcReal afRoot[4];
kPoly.GetAllRoots(iCount,afRoot);
if ( iCount > 0 )
{
MgcReal fMinDistSqr = MgcMath::INFINITY;
int iMinIndex = -1;
MgcReal afV[2], fDenom;
for (int iIndex = 0; iIndex < iCount; iIndex++)
{
// compute closest point for this root
for (i = 0; i < 2; i++)
{
fDenom = 1.0 + 2.0*afRoot[iIndex]*afD[i];
afV[i] = (afA[i]-afRoot[iIndex]*afB[i])/fDenom;
}
rkClosest.x = aafR[0][0]*afV[0] + aafR[1][0]*afV[1];
rkClosest.y = aafR[0][1]*afV[0] + aafR[1][1]*afV[1];
// compute squared distance from point to quadratic
MgcVector2 kDiff = rkClosest - rkPoint;
MgcReal fDistSqr = kDiff.SquaredLength();
if ( fDistSqr < fMinDistSqr )
{
fMinDistSqr = fDistSqr;
iMinIndex = iIndex;
}
}
for (i = 0; i < 2; i++)
{
fDenom = 1.0 + 2.0*afRoot[iMinIndex]*afD[i];
afV[i] = (afA[i]-afRoot[iMinIndex]*afB[i])/fDenom;
}
rkClosest.x = aafR[0][0]*afV[0] + aafR[1][0]*afV[1];
rkClosest.y = aafR[0][1]*afV[0] + aafR[1][1]*afV[1];
return fMinDistSqr;
}
else
{
// should not happen
assert( false );
return -1.0;
}
}
//----------------------------------------------------------------------------
MgcReal MgcDistance (const MgcVector2& rkPoint, const MgcReal afQuad[6],
MgcVector2& rkClosest)
{
return MgcMath::Sqrt(MgcSqrDistance(rkPoint,afQuad,rkClosest));
}
//----------------------------------------------------------------------------
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -