asmsearch.cpp

这是个人脸识别程序

if (sscanf(s, "%s %d %s", sLev, &nLev, sFirstShape) != 3)
	SysErr("Can't read level %d header from %s, line is %s", iLev, sAsmFile, s);

if (strcmp(sLev, "Lev") != 0 || strcmp(sFirstShape, "FirstShape") != 0)
	SysErr("Can't read level %d header from %s, line is %s", iLev, sAsmFile, s);

if (nLev != iLev)
	SysErr("Illegal level %d in level header from %s, expected %d, line is %s", nLev, sAsmFile, iLev, s);
}

//-----------------------------------------------------------------------------
// This gets called AFTER initializing the following for Asm[0]: nPoints EigVals EigVecs EigInverse nEigs

static void ReadAsmLev (tAsmLev &AsmLev,							// out
						  int iLev, 								// in
						  tAsmLev AsmLevs[], int nPoints,			// in
						  const char sAsmFile[], FILE *pAsmFile)	// in
{
CheckAsmLevHeader(sAsmFile, pAsmFile, iLev);

ReadProfSpecs(AsmLev.ProfSpecs, nPoints, sAsmFile, pAsmFile);

AsmLev.Covars.resize(nPoints);
AsmLev.SparseCovars.resize(nPoints);

for (int iPoint = 0; iPoint < nPoints; iPoint++)
	ReadMatVec(&AsmLev.Covars[iPoint], NULL, sAsmFile, pAsmFile, NULL, 0, 0, nSubProfsForProfSpec(AsmLev.ProfSpecs[iPoint]));

AsmLev.Profs.resize(nPoints);

for (iPoint = 0; iPoint < nPoints; iPoint++)
	ReadMatVec(&AsmLev.Profs[iPoint], NULL, sAsmFile, pAsmFile, NULL, 0, 0, nSubProfsForProfSpec(AsmLev.ProfSpecs[iPoint]));

// convert sparse covar matrices from Mat to SparseMat

for (iPoint = 0; iPoint < nPoints; iPoint++)
	{
	int nSubProfs = nSubProfsForProfSpec(AsmLev.ProfSpecs[iPoint]);
	for (int iSub = 0; iSub < nSubProfs; iSub++)
		if (fSparseMat(AsmLev.Covars[iPoint][iSub]))
			{
			if (AsmLev.SparseCovars[iPoint].size() < nSubProfs)
				AsmLev.SparseCovars[iPoint].resize(nSubProfs);
			CopyMatToSparseMat(AsmLev.SparseCovars[iPoint][iSub], AsmLev.Covars[iPoint][iSub]);
			AsmLev.Covars[iPoint][iSub].dim(1, 1);		// we no longer need Covars, just SparseCovars TODO tidy up, see GetProfileFit()
			}
	}
}

//-----------------------------------------------------------------------------
// return true if there are 2D profs in the given model level

template<typename ModelT>
static bool fTwoDProfs (const ModelT &Model, int nPoints)
{
for (int iPoint = 0; iPoint < nPoints; iPoint++)
	if (Model.ProfSpecs[iPoint] & PROF_2d)				// two dimensional profile?
		return TRUE;

return FALSE;
}

//-----------------------------------------------------------------------------
// Init active shape model by reading an ASM file and by using #defines in masmconf.hpp
// TODO rationalize use of #defines in CALLS to this function and use of #defines IN this function

void InitAsm (tAsmModel &Asm,				// out
			  const char sAsmFile[],		// in
			  bool fFancyShapeModel,		// in
			  int nShapeModelIters,			// in
			  eStartMethod StartMethod)		// in
{
static clock_t StartTime = clock();
lprintf("Init ASM %s ", sAsmFile);
FILE *pAsmFile = Fopen(sAsmFile, "r");

// header

int nLevs;
ReadAsmHeader(Asm.nPoints, nLevs, Asm.PyrRatio, Asm.PyrReduceMethod,
				Asm.fExplicitPrevNext, Asm.NormalizedProfLen, Asm.nScaledFaceWidth,
				Asm.fBilinearRescale, Asm.nTrimCovar, Asm.SigmoidScale,
				sAsmFile, pAsmFile);

Asm.nStartLev = nLevs - 1;

extern double SigmoidScaleGlobal;
SigmoidScaleGlobal = Asm.SigmoidScale;	//TODO ugly

// eigen values and vectors

Asm.AsmLevs[0].EigVals.sread(sAsmFile, pAsmFile);
Asm.EigVecs.sread(sAsmFile, pAsmFile);
Asm.EigInverse = Asm.EigVecs.t();		// inverse of EigVecs (EigVecs is orthogonal, so transpose=inverse)

// file mean shape

Asm.AsmLevs[0].MeanShape.sread(sAsmFile, pAsmFile);

ASSERT(Asm.nPoints == Asm.AsmLevs[0].MeanShape.nrows());
ASSERT(Asm.EigVecs.nrows() == 2 * Asm.nPoints);

Asm.StartMethod = StartMethod;
Asm.nShapeModelIters = nShapeModelIters;
Asm.fFancyShapeModel = fFancyShapeModel;

for (int iLev = nLevs-1; iLev >= 0; iLev--)		// for each level in the image pyramid
	ReadAsmLev(Asm.AsmLevs[iLev], iLev, Asm.AsmLevs, Asm.nPoints, sAsmFile, pAsmFile);

fclose(pAsmFile);

Asm.FileMeanShape = Asm.AsmLevs[0].MeanShape;

// initializations from global search #defines

bool fTwoDProfsLev0 = fTwoDProfs(Asm.AsmLevs[0], Asm.nPoints);	// true if there are ANY 2d profs on pyramid level 0

if (fTwoDProfsLev0)
	{
	Asm.nPixSearch = CONF_nPixSearch2d;
	Asm.nMaxSearchIters = CONF_nMaxSearchIters2d;
	Asm.nQualifyingDisplacments = CONF_nQualifyingDisplacments2d;
	}
else
	{
	Asm.nPixSearch = CONF_nPixSearch;
	Asm.nMaxSearchIters = CONF_nMaxSearchIters;
	Asm.nQualifyingDisplacments = CONF_nQualifyingDisplacments;
	}
for (iLev = 0; iLev < nLevs; iLev++)		// for each level in the image pyramid except level 0 (which was done above)
	{
	if (iLev == 0)
		{
		if (fTwoDProfsLev0)
			{
			if (CONF_nEigsLev0_2d)
				{
				Asm.AsmLevs[0].nEigs = CONF_nEigsLev0_2d;
				Asm.AsmLevs[0].BMax = CONF_BMaxLev0_2d;
				}
			else
				{
				Asm.AsmLevs[0].nEigs = CONF_nEigs;
				Asm.AsmLevs[0].BMax = CONF_BMax;
				}
			}
		else
			{
			if (CONF_nEigsLev0)
				{
				Asm.AsmLevs[0].nEigs = CONF_nEigsLev0;
				Asm.AsmLevs[0].BMax = CONF_BMaxLev0;
				}
			else
				{
				Asm.AsmLevs[0].nEigs = CONF_nEigs;
				Asm.AsmLevs[0].BMax = CONF_BMax;
				}
			}
		}
	else	// not iLev0
		{
		Asm.AsmLevs[iLev].nEigs = CONF_nEigs;
		Asm.AsmLevs[iLev].BMax = CONF_BMax;
		}
	if (CONF_nEigsFinal)
		{
		Asm.AsmLevs[iLev].nEigsFinal = CONF_nEigsFinal;
		Asm.AsmLevs[iLev].BMaxFinal = CONF_BMaxFinal;
		}
	else
		{
		Asm.AsmLevs[iLev].nEigsFinal = CONF_nEigs;
		Asm.AsmLevs[iLev].BMaxFinal = CONF_BMax;
		}
	Asm.AsmLevs[iLev].EigVals = Asm.AsmLevs[0].EigVals / pow(Asm.PyrRatio * Asm.PyrRatio, iLev);

	// In theory we shouldn't have to scale MeanShape here because ConformShapeToModel takes care of scaling
	// In practice we get better results when we do scale here
	// TODO look into this

	Asm.AsmLevs[iLev].MeanShape = Asm.AsmLevs[0].MeanShape / GetPyrScale(iLev, Asm.PyrRatio);
	}
double TimeTaken = (double)(clock() - StartTime) / CLOCKS_PER_SEC;
lprintf(" [Init ASM in %.2f secs]\n", TimeTaken);
}

//-----------------------------------------------------------------------------
// Generate a C file of an asm model
// You will have to hand tweak the C file later, this just gets you started

void GenerateCFileFromAsmModel(tAsmModel &m)
{
int nPoints = m.nPoints;
int nModel = m.iModel;
int iLev, i;
char sFile[SLEN];
char s[SLEN];
sprintf(sFile, "model%d.cpp", nModel+1);
lprintf("Generating %s\n", sFile);
FILE *pFile = Fopen(sFile, "w");
#define P fprintf(pFile,

P "// %s\n", sFile);
P "// This file was generated by GenerateCFileFromAsmModel.\n");
P "// The arrays are split into pieces so the C compiler doesn't choke.\n");
P "// Stephen Milborrow\n");
P "\n");
P "#include <all.hpp>\n");
P "#define NPOINTS %d\n", m.AsmLevs[0].ProfSpecs.size());
P "\n");
m.AsmLevs[0].EigVals.print("const double EigValsData[NPOINTS * 2] = {\n", "%g, ", sFile, pFile);
P "};\n");
m.EigVecs.print("static const double EigVecsData[NPOINTS * 2 * NPOINTS * 2] = {\n", "%g, ", sFile, pFile);
P "};\n");
m.FileMeanShape.print("const double FileMeanShapeData[NPOINTS * 2] = {\n", "%g, ", sFile, pFile);
P "};\n");
for (iLev = 0; iLev < 4; iLev++)
	{
	P "static unsigned ProfSpecs%dData[%d] = {", iLev, nPoints);
	for (i = 0; i < nPoints; i++)
		P "0x%x, ", m.AsmLevs[iLev].ProfSpecs[i]);
	P "};\n");
	}
for (iLev = 0; iLev < 4; iLev++)
	for (i = 0; i < nPoints; i++)
		{
		int nelems = m.AsmLevs[iLev].Covars[i][0].nelems();
		if (nelems > 1)
			{
			sprintf(s, "static double CovarLev%dPoint%dData[%d] = {\n", iLev, i, nelems);
			m.AsmLevs[iLev].Covars[i][0].print(s, "%g, ", sFile, pFile);
			}
		else	// sparse matrix
			{
			nelems = m.AsmLevs[iLev].SparseCovars[i][0].size();
			P "static tSparseElem CovarLev%dPoint%dData[%d] = {\n", iLev, i, nelems);
			for (int iElem = 0; iElem < nelems; iElem++)
				{
				tSparseElem *p = &m.AsmLevs[iLev].SparseCovars[i][0][iElem];
				P "%d, %d, %.1fF,\n", p->iRow, p->iCol, p->Val);
				}
			}
		P "};\n");
		sprintf(s, "static double ProfLev%dPoint%dData[%d] = {\n", iLev, i, m.AsmLevs[iLev].Profs[i][0].nelems());
		m.AsmLevs[iLev].Profs[i][0].print(s, "%g, ", sFile, pFile);
		P "};\n");
		}
P "static void InitProf(tAsmLev *pAsmLev, const int iPoint, const double ProfData[], const int nProfLen)\n");
P "{\n");
P "	Mat Prof(ProfData, 1, nProfLen);\n");
P "	pAsmLev->Profs[iPoint].resize(1);	//needed?\n");
P "	pAsmLev->Profs[iPoint][0] = Prof;   // iSub=0\n");	// Models with more than one subprofile aren't supported
P "}\n");
P "static void InitCovar(tAsmLev *pAsmLev, const int iPoint, const double CovarData[], const int nProfLen)\n");
P "{\n");
P "	Mat Covar(CovarData, nProfLen, nProfLen);\n");
P "	pAsmLev->Covars[iPoint].push_back(Covar);\n");
P "}\n");
P "void InitModel (tAsmModel &m)\n");
P "{\n");
P "int i;\n");
P "m.iModel = %d;\n", nModel);
P "strcpy(m.sFile, \"%s\");\n", m.sFile);
P "m.nPoints = %d;\n", m.nPoints);
P "m.nStartLev = %d;\n", m.nStartLev);
P "m.PyrRatio = %g;\n", m.PyrRatio);
P "m.PyrReduceMethod = %d;\n", m.PyrReduceMethod);
P "m.nMaxSearchIters = %d;\n", m.nMaxSearchIters);
P "m.nShapeModelIters = %d;\n", m.nShapeModelIters);
P "m.nQualifyingDisplacments = %d;\n", m.nQualifyingDisplacments);
P "m.nTrimCovar = %d;\n", m.nTrimCovar);
P "m.SigmoidScale = %g;\n", m.SigmoidScale);
P "m.nScaledFaceWidth = %d;\n", m.nScaledFaceWidth);
P "m.nPixSearch = %d;\n", m.nPixSearch);
P "m.fExplicitPrevNext = 0;  // revisit?\n");
P "m.fFancyShapeModel = 1;  // revisit?\n");
P "m.fBilinearRescale = 1;  // revisit?\n");
P "m.NormalizedProfLen = %g;\n", m.NormalizedProfLen);
P "m.StartMethod = SM_UseNnFaceAndEyesRealignEyes; // revisit?\n", m.StartMethod);
P "SHAPE FileMeanShape(FileMeanShapeData, NPOINTS, 2);\n");
P "m.FileMeanShape = FileMeanShape;\n");
P "SHAPE EigVecs(EigVecsData, NPOINTS * 2, NPOINTS * 2);\n");
P "m.EigVecs = EigVecs;\n");
P "m.EigInverse = m.EigVecs.t();\n");
P "Mat EigVals(EigValsData, NPOINTS * 2, 1);\n");
P "for (int iLev = 0; iLev < 4; iLev++)\n");
P "    {\n");
P "     m.AsmLevs[iLev].nEigs = %d;       // revisit?\n", CONF_nEigs);
P "     m.AsmLevs[iLev].nEigsFinal = %d;  // revisit?\n", CONF_nEigs);
P "     m.AsmLevs[iLev].BMax = %g;       // revisit?\n", CONF_BMax);
P "     m.AsmLevs[iLev].BMaxFinal = %g;  // revisit?\n", CONF_BMax);
P "     m.AsmLevs[iLev].MeanShape = m.FileMeanShape / GetPyrScale(iLev, m.PyrRatio);\n");
P "     m.AsmLevs[iLev].EigVals = EigVals / pow(m.PyrRatio * m.PyrRatio, iLev);\n");
P "     }\n");
P "tAsmLev *pAsmLev;\n");

for (iLev = 0; iLev < 4; iLev++)
	{
	P "\n");
	P "// level %d\n", iLev);
	P "\n");
	P "pAsmLev = &m.AsmLevs[%d];\n", iLev);
	P "pAsmLev->ProfSpecs.resize(NPOINTS);\n");
	P "for (i = 0; i < NPOINTS; i++) \n");
	P "    pAsmLev->ProfSpecs[i] = ProfSpecs%dData[i];\n", iLev);
	P "pAsmLev->Profs.resize(NPOINTS);\n");
	P "pAsmLev->Covars.resize(NPOINTS);\n");
	for (i = 0; i < nPoints; i++)
		{
		int nProfLen = 	m.AsmLevs[iLev].Profs[i][0].nelems();
		P "InitProf(pAsmLev,  %d, ProfLev%dPoint%dData, %d);\n", i, iLev, i, nProfLen);
		P "InitCovar(pAsmLev, %d, CovarLev%dPoint%dData, %d);\n", i, iLev, i, nProfLen);
		}
	}
P "}\n");
fclose(pFile);
#undef P
}

//-----------------------------------------------------------------------------
// If there are multiple sub-profiles within a profile this returns the sub-profile width

int nGetProfWidthFromModel (int iPoint, int iSub, const tAsmLev &AsmLev)
{
return nGetProfWidth(AsmLev.Profs[iPoint][iSub].ncols(), AsmLev.ProfSpecs[iPoint]);
}