masm.cpp

这是个人脸识别程序

		SumShapes(NewAvShape, Shapes[iShape], AltShape);
		}
	// constrain new average shape by aligning it to the ref shape

	AlignShape(NewAvShape, RefShape);
	double Distance = AvShape.distSquared(NewAvShape);
	if (Distance < CONF_MaxShapeAlignDist)
		break;	// converged
	AvShape = NewAvShape;
	}
if (iPass > CONF_nMaxAlignPasses)
	Err("Didn't align in %d passes", CONF_nMaxAlignPasses);

// align all the shapes to the average shape

for (iShape = 0; iShape < nShapes; iShape++)
	AlignShape(Shapes[iShape], AvShape);		//TODO was was CentralizeShape(Shapes[iShape])

#if CONF_fShowDebugImages
ShowCentralizedShapes(Shapes, AvShape);
#endif
}

//-----------------------------------------------------------------------------
// Print the angle wrt horizontal of the line joining the eye pupils in Shape

static void PrintEyeAngle (SHAPE &Shape, char sFormat[])
{
const tLand *pLand = pSelectLand(Shape.nrows());		// landmark information

if (iTranslatedPoint(MLEye) < Shape.nrows()
	&& iTranslatedPoint(MREye) < Shape.nrows()
	&& fPointUsed(Shape, iTranslatedPoint(MLEye))
	&& fPointUsed(Shape, iTranslatedPoint(MREye)))
	{
	lprintf(sFormat,
		Degrees(atan2(Shape(iTranslatedPoint(MREye), VY) - Shape(iTranslatedPoint(MLEye), VY),
					  Shape(iTranslatedPoint(MREye), VX) - Shape(iTranslatedPoint(MLEye), VX))));
	}
}

//-----------------------------------------------------------------------------
static void AlignTrainingShapes (tShapes &Shapes,	// io: Shapes field and AvShape updated
								 int iRefShape)		// in
{
char *sShapeName = &Shapes.TagStrings[iRefShape][FNAME_OFFSET];
lprintf("Aligning shapes to %s",  sShapeName);

DASSERT(iRefShape < Shapes.nShapes);

int nPoints = Shapes.Shapes[iRefShape].nrows();
Shapes.AvShape.dimKeep(nPoints, 2);

SHAPE RefShape(Shapes.Shapes[iRefShape]);
#if CONF_fStraightenBeforeAligning
StraightenAndCentralizeShape(RefShape);
#else
CentralizeShape(RefShape);
#endif
if (CONF_GenSubModel == SM_All || CONF_GenSubModel == SM_20)
	PrintEyeAngle(Shapes.AvShape, " (eye angle %.2f degrees)");

AlignShapes(Shapes.Shapes, Shapes.AvShape, RefShape);

if (CONF_GenSubModel == SM_All || CONF_GenSubModel == SM_20)
	PrintEyeAngle(Shapes.AvShape, ", aligned eye angle %.2f degrees");

lprintf(" %s\n", sSecsElapsed());

// sanity check: not necessary but prevents slip-ups during classifier evaluation

#if CONF_nReducedNbrPoints == 0
if (Shapes.nUsedLandmarks[iRefShape] != CONF_nPointsXm2vts &&
	Shapes.nUsedLandmarks[iRefShape] != CONF_nPoints84)
	{
	Warn("The reference shape (%s shape number %d) has an unusual number %d of landmarks",
			sShapeName, iRefShape, Shapes.nUsedLandmarks[iRefShape]);
	}
#else
if (Shapes.nUsedLandmarks[iRefShape] != CONF_nReducedNbrPoints)
	Err("The reference shape (%s shape number %d) has %d landmarks, expected same number of landmarks as XMTVTS or M4",
		sShapeName, iRefShape, Shapes.nUsedLandmarks[iRefShape]);
#endif
}

//-----------------------------------------------------------------------------
// if CONF_fMasmGaussianNoise is true, return a gaussian rand value with a std dev of Range
//									(i.e. 68% iof the values will be between -Range and +Range)
//
// else return a rand value between -Range and +Range, uniform distribution

static double Rand1 (double Range)
{
#if CONF_fMasmGaussianNoise
return RandGauss(Range);
#else
return Range * (RandDouble(2) - 1);
#endif
}

//-----------------------------------------------------------------------------
// This calculates the eigen vecs and vals for Shapes.
// There is a twist: we only the shapes that have the same number of used landmarks as the first shape.
// TODO One day it would be nice to include all shapes

static void CalcShapeEigs (Mat &EigVals, Mat &EigVecs, 						   									// out
						 const ShapeVec &Shapes, SHAPE &AvShape, size_t nShapes, const IntVec &nUsedLandmarks)	// in
{
int	nPoints = AvShape.nrows(); int nShapesWithAllLandmarks = 0;

for (int iShape = 0; iShape < nShapes; iShape++)
	if (nUsedLandmarks[iShape] == nUsedLandmarks[0])
		nShapesWithAllLandmarks++;

if (nShapesWithAllLandmarks != nShapes)
	Warn("Only %d of %d shapes used in CalcShapeEigs because of missing landmarks", nShapesWithAllLandmarks, nShapes);

Mat ShapeDelta(nShapesWithAllLandmarks, 2 * nPoints);

VecView AvShapeAsRow(AvShape.viewAsRow()); 	// view all coords for shape as a single row vector, first x coords then y coords

int iShape1 = 0;
for (iShape = 0; iShape < nShapes; iShape++)
	if (nUsedLandmarks[iShape] == nUsedLandmarks[0])
		{
		SHAPE Shape(Shapes[iShape]);
#if CONF_fRecentralizeShapes			// TODO why does centralizing here give different results (A+3.6% B-4.4% m+1.3%)
#if CONF_fStraightenBeforeAligning
		StraightenAndCentralizeShape(Shape);
#else
		CentralizeShape(Shape);
#endif
#endif
		if (CONF_GenSubModel == SM_All &&
				CONF_xStretchShape != 0.0 || CONF_xRhsStretchShape != 0.0 || CONF_ShapeNoise != 0.0)
			{
			// add noise to the shape

			double xStretch = Rand1(CONF_xStretchShape);		// rand val, 68% of values between +- CONF_xStretchShape
			double xDist    = Rand1(CONF_xRhsStretchShape);
			for (int iRow = 0; iRow < Shape.nrows(); iRow++)
				{
				if (Shape(iRow, VX) > 0)						// stretch rhs of shape horizontally
					Shape(iRow, VX) *= (1 + xDist);
				Shape(iRow, VX) *= (1 + xStretch);				// stretch entire shape horizontally
				Shape(iRow, VX) += Rand1(CONF_ShapeNoise);		// add noise to both x and y
				Shape(iRow, VY) += Rand1(CONF_ShapeNoise);
				}
			}
		ShapeDelta.row(iShape1) = Shape.viewAsRow() - AvShapeAsRow;
		iShape1++;
		}

DASSERT(iShape1 == nShapesWithAllLandmarks);

Mat Covar((ShapeDelta.t() * ShapeDelta) / nShapesWithAllLandmarks);	// dimensions are 2*nPoints x 2*nPoints

EigVecs = GetEigsForSymMat(Covar, EigVals);							// returns eig vecs sorted on eig vals

// For improved consistency between _DEBUG and release versions, set
// very small eig vals and their corresponding eig vecs to 0.
// These small values are the unpredictable result of numerical errors but it's not
// clear why the GSL generates slightly different vals for _DEBUG and release versions.

for (int iEig = 0; iEig < EigVals.nrows(); iEig++)
	if (EigVals(iEig) < EigVals(0) / 1E6)
		{
		EigVals(iEig) = 0;
		EigVecs.col(iEig).fill(0);
		}
}

//-----------------------------------------------------------------------------
// Print first n eig vals, up to and including the first small value. But no more than 10 of them.

static void ShowShapeEigs (const Mat &EigVals)
{
char s[PLEN];
int iPrint = 0;

double MinEig = EigVals(0) / 1000;
while (iPrint < EigVals.nelems() && EigVals(iPrint) > MinEig && iPrint < 10-1)
	iPrint++;

if (iPrint < EigVals.nelems() + 1)	// show first small eig value, for context
	iPrint++;

sprintf(s, "%s%d eig vals: ", (iPrint < EigVals.nelems()? "First ": ""), iPrint);

EigVals.t().print(s, "%.1f ", NULL, NULL, iPrint);

if (fabs(EigVals(0)) < 0.1)			// number is somewhat arbitrary
	Warn("Eigen data invalid (not enough variation in shape file)");
}

//-----------------------------------------------------------------------------
// Reserve space in the shape Stats

static void AllocStats (tStats &Stats, 			// io
						char sProfSpecs[],		// in
						const tShapes &Shapes)	// in
{
DASSERT(fOdd(CONF_nProfWidth1d));
DASSERT(fOdd(CONF_nProfWidth2d));

int nPoints = Shapes.nPoints;
int nMaxProfsPerPoint = Shapes.nShapes;	// max nbr of profs per landmark, not all used if some landmarks rejected by Shapes.TagStrings[iShape]

Stats.nProfs.resize(CONF_nLevs);
Stats.Covars.resize(CONF_nLevs);
Stats.Profs.resize(CONF_nLevs);
Stats.AvProfs.resize(CONF_nLevs);

static unsigned LastProfSpec = 0;
sProfSpecs[0] = 0;
for (int iLev = 0; iLev < CONF_nLevs; iLev++)
	{
	Stats.nProfs[iLev].resize(nPoints, 0);
	Stats.Covars[iLev].resize(nPoints, 0);
	Stats.Profs[iLev].resize(nPoints, 0);
	Stats.AvProfs[iLev].resize(nPoints, 0);

	for (int iPoint = 0; iPoint < nPoints; iPoint++)
		if (fPointUsed(Shapes.AvShape, iPoint))	// not all points are used if we are generating a submodel
			{
			// Here we make the assumption that all sub profs need no more space than SubProf0.
			// This can waste space but is otherwise harmless.

			unsigned ProfSpec = GetGenProfSpec(iLev, iPoint);
			int nSubProfs = nSubProfsForProfSpec(ProfSpec);
			int nSubProfPoints = nGenelemsPerSubProf(GetSubProfSpec(ProfSpec, 0));

			AllocMatVec(Stats.Covars[iLev][iPoint],  nSubProfs, nSubProfPoints,    nSubProfPoints);
			AllocMatVec(Stats.Profs[iLev][iPoint],   nSubProfs, nMaxProfsPerPoint, nSubProfPoints);
			AllocMatVec(Stats.AvProfs[iLev][iPoint], nSubProfs, 1, 		  		   nSubProfPoints);

			if (LastProfSpec != ProfSpec)
				{
				char s[SLEN]; sprintf(s, "%d,%d:%8.8x ", iLev, iPoint, ProfSpec);
				strcat(sProfSpecs, s);
				LastProfSpec = ProfSpec;
				}
			}
	}
}

//-----------------------------------------------------------------------------
// Compact Profs by keeping elements that are in used, discarding the rest.
// fUsed tells us which elements are used.

#if CONF_fMultiEdit || CONF_fCondense || CONF_fReduce
static void CompactProfs (int &nProfs, Mat &Profs, Vec *pProfs, CharVec *pLabs,	// io
						  const BoolVec &fUse, char sTitle[], bool fVerbose=false)		// in
{
int ncols = Profs.ncols();
int iProfNew = 0;

for (int iProf = 0; iProf < nProfs; iProf++)
	{
	if (fUse[iProf])	// profile used?
		{
		for (int i = 0; i < ncols; i++)
			Profs(iProfNew, i) = Profs(iProf, i);

		iProfNew++;
		}
	}
Profs.dimKeep(iProfNew, ncols);

DASSERT(iProfNew > 0);
nProfs = iProfNew;
}
#endif CONF_fMultiEdit || CONF_fCondense || CONF_fReduce

//-----------------------------------------------------------------------------
// This routine is for debugging.  It shows the image with the shape and
// whiskers superimposed.  We call this just before sampling profiles so we can
// see that what we are sampling is correct i.e the whiskers are in the right places etc.

#if CONF_fShowSampledImageWithShapes
static void ShowSampledImageWithShapes (int iLev, int iShape, Image &ScaledImg, 
				SHAPE &ScaledShape, SHAPE &ScaledAlignedAvShape, const tShapes &Shapes)
{
if (CONF_fShowSampledImageAllLevs || iLev == 0)
	{
	double Scale = pow(2,iLev);
	#define MAG 1.0
	RgbImage RgbImg(MAG * ScaledImg.width, MAG * ScaledImg.height);
	Scale *= MAG;

	DrawShape1(RgbImg, ScaledImg, ScaledShape, NULL /* &ScaledAlignedAvShape */,
		C_DRED, Scale, IM_NO_CONNECT_DOTS,
		(MAG > 1? IM_ANNOTATE: IM_NO_ANNOTATE), IM_NO_WHISKERS, 0, IM_TRANSPARENT, IM_NO_CROP_TO_FACE,
		gLandTabAll, &ScaledAlignedAvShape, IM_NO_DRAW_CIRCLES);

#if CONF_fCropImageToShape
	CropImageToFace(RgbImg, ScaledAlignedAvShape, Scale);
#endif
	char sDrive[_MAX_DRIVE], sDir[_MAX_DIR], sFname[_MAX_FNAME], sExt[_MAX_EXT], sPath[SLEN];
	_splitpath(&Shapes.TagStrings[iShape][FNAME_OFFSET], sDrive, sDir, sFname, sExt);
	sprintf(sPath, "%s/%s_%d", CONF_sOutDir, sFname, iLev);
	WriteBmp(RgbImg, sPath, VERBOSE);
	}
}
#endif

//-----------------------------------------------------------------------------
static void CollectProfForOneLandmark (tStats &Stats,								// io
					int iShape, int iLev, int iPoint, const tLand *pLand,			// in
					const Image &ScaledImg, 										// in
					const SHAPE &ScaledShape, const SHAPE &ScaledAlignedAvShape,	// in
					const MatVec &Grads)											// in
{
unsigned ProfSpec = GetGenProfSpec(iLev, iPoint);

if ((ProfSpec & PROF_2d) == 0)		// 1D prof? if so save time for multiple iOffsets by pre-preparing profile
	PrepareProf1D(ScaledImg, ScaledShape, ProfSpec, pLand, ScaledAlignedAvShape,
					iPoint, CONF_nProfWidth1d, 0, 0);

int iProf = Stats.nProfs[iLev][iPoint];		// start where we left off last time
const int iOrthOffset = 0;					// constant for now (will need multiple iOffsets for nearest neighbors when I implement it TODO)
const int iOffset = 0;						// ditto
bool fSuccess = true;

if (ProfSpec & PROF_2d)						// two dimensional prof?
	{
#if CONF_fDebug2dProfs
	lprintf("\n2dProf at iPoint %d: ", iPoint);
	// draw a square around the landmark and show the image
	static int iBmp = 0; char s[SLEN]; sprintf(s, "out/Img%2.2d.bmp", iBmp++);
	SHAPE SubShape(5, 2);
	int x = ScaledShape(iPoint, VX), y = ScaledShape(iPoint, VY), Half = (CONF_nProfWidth2d-1)/2 + 1;
	SubShape(0, VX) = x - Half; SubShape(0, VY) = y - Half;
	SubShape(1, VX) = x - Half; SubShape(1, VY) = y + Half;
	SubShape(2, VX) = x + Half; SubShape(2, VY) = y + Half;
	SubShape(3, VX) = x + Half; SubShape(3, VY) = y - Half;