masm.cpp

这是个人脸识别程序

//
// Warning messages are printed if the inversion process is not 100% standard
// i.e. if any of the above processing is done.

static void GenInvertedCovarMat (Mat &Covar,							// out
						 		 Mat &Profs, const Mat &AvProfs,		// in: except that one prof is deleted if linear dependence
								 int nProfs, unsigned ProfSpec)			// in
{
char cStatus = 0;

if (CONF_fDropProfToPreventSingCovar && fProfileNormalizationCausesLinearDependence(ProfSpec))
	{
	static bool fWarned = false;
	if (!fWarned)
		lprintf("\nDropping prof to prevent linear dependence ");
	fWarned = true;
	ASSERT(nProfs == Profs.nrows());
	nProfs--;
	Profs.dimKeep(nProfs, Profs.ncols());		// drop final row of Profs
	cStatus = 'd';
	}
DASSERT(nProfs >= 1);

Covar = (Profs.t() * Profs - AvProfs.t() * AvProfs * nProfs) / nProfs;

int nAllowedZeroEigs = 0;
if (fCheckCovarMat(Covar, nAllowedZeroEigs))
	Covar.invertMe();
else if (CONF_fDropEigToPreventSingCovar)
	{
	int nColinear = 1;
	lprintf("\nDropping eigs ");
	InvertMatRemovingSmallEigs(Covar, nColinear); // remove small eigenvalues in Covar and then invert using principal components
	if (nColinear)
		lprintf("dropped %d eig(s) ", nColinear);
	nAllowedZeroEigs = 1;
	cStatus = 's';
	}
else if (CONF_RidgeLambdaForSingCovars != 0)
	{
	InvertMatWithRidgeRegresssion(Covar);	// do ridge regression on Covar and then invert it
	cStatus = 'r';
	}
if (CONF_MaxEigRatio != 0)
	RemoveSmallEigs(Covar, nAllowedZeroEigs, CONF_MaxEigRatio);

if (!fCheckCovarMat(Covar, nAllowedZeroEigs))
	{															// inverted matrix still not good (even after above processing)
	cStatus = 'I';
	Covar.identityMe();											// turn covariance matrix into an identity matrix
	}
if (cStatus != 0)
	lprintf("%c", cStatus);
}

//-----------------------------------------------------------------------------
// This generates the (inverted) Covar and AvProfs for one landmark on one level

static void GenAsmStats (Mat &Covar, Mat &AvProf, 			// out
						 Mat &Profs, 						// in: except that one prof is deleted if linear dependence
						 int nProfs, unsigned ProfSpec) 	// in
{
GenAvProf(AvProf, Profs, nProfs);
GenInvertedCovarMat(Covar, Profs, AvProf, nProfs, ProfSpec);
}

//-----------------------------------------------------------------------------
// Dump feature iFeature for all level 0 profiles for all shapes
// iFeature is an index in the profile vector, feature at landmark is iFeature=CONF_nProfWidth1d/2+1
// I used this to generate density histograms to check for normality with R.

#if CONF_fDumpLandmark
static void DumpLandmark (int iLandmark, const char sShapeFile[], const tStats &Stats, const tShapes &Shapes)
{
char sDrive[_MAX_DRIVE], sDir[_MAX_DIR], sFileBase[_MAX_FNAME], sExt[_MAX_EXT];
_splitpath(sShapeFile, sDrive, sDir, sFileBase, sExt);	// generate name from .shape file name
char s[SLEN]; sprintf(s, "out/%c%c%clandmark%2.2d.tab",
				sFileBase[0], (sFileBase[1]? sFileBase[1]: '_'),
				(sFileBase[2]? sFileBase[2]: '_'), iLandmark);
FILE *pFile = Fopen(s, "w");
lprintf("Generating %s ", s);
const int iLev = 0;
const int iSub = 0;
int ncols = Stats.Profs[iLev][iLandmark][iSub].ncols();
Fprintf(pFile, "# %s\n", s);
Fprintf(pFile, "     ");
for (int iCol = 0; iCol < ncols; iCol++)	// print header
	{
	sprintf(s, "x%d", iCol+1);
	Fprintf(pFile, "%8s ", s);
	}
Fprintf(pFile, "\n");
for (int iProf = 0; iProf < Stats.nProfs[iLev][iLandmark]; iProf++)
	{
	Fprintf(pFile, "%4d ", iProf+1);
	for (iCol = 0; iCol < ncols; iCol++)
		Fprintf(pFile, "%8.4g ", Stats.Profs[iLev][iLandmark][iSub](iProf, iCol));
	Fprintf(pFile, "\n");
	}
fclose(pFile);
lprintf("\n");
}
#endif

//-----------------------------------------------------------------------------
// Generate profile stats to matrices Stats.AvProfs, Stats.AvProfs, Stats.Covars

static void GenProfStats (tStats &Stats, 		// io
							tShapes &Shapes,	// in: but will be modified if CONF_nGenPrescaleFaceWidth != 0
							const char sShapeFile[])
{
int iPoint;

// Phase1: collect profiles i.e. fill in Stats.Profs

if (CONF_fAlignImageToAlignedShape)
	lprintf("Translate images to aligned shape frame before taking profiles\n");

lprintf("Reading images ");
logprintf("\n");

CollectProfs(Stats, Shapes);

// Phase2: release some memory by resizing the profile data to the size we are actually using
// We over-allocated (to nMaxProfsPerPoint) in AllocStats() earlier.

for (int iLev = 0; iLev < CONF_nLevs; iLev++)
	for (iPoint = 0; iPoint <  Shapes.nPoints; iPoint++)
		if (fPointUsed(Shapes.AvShape, iPoint))				// not all points are used if we are generating a submodel
			{
			int nSubProfs = nSubProfsForProfSpec(GetGenProfSpec(iLev, iPoint));
			for (int iSub = 0; iSub < nSubProfs; iSub++)
				{
				const int ncols = Stats.Profs[iLev][iPoint][iSub].ncols();
				Stats.Profs[iLev][iPoint][iSub].dimKeep(Stats.nProfs[iLev][iPoint], ncols);
				}
			}

// Phase3: calculate ASM summary stats i.e. fill in Stats.AvProfs and Stats.Covars from Stats.Profs

#if CONF_fDumpLandmark
DumpLandmark(MLEyeInner, sShapeFile, Stats, Shapes);
Err("Early exit after dumping landmark");
#endif

int nTrimmedMats = 0, nBadTrimmedMats = 0;

lprintf("Generating stats ");
InitPacifyUser(Shapes.nPoints);
for (iPoint = 0; iPoint < Shapes.nPoints; iPoint++)
	{
	PacifyUser(iPoint);
	if (fPointUsed(Shapes.AvShape, iPoint))				// not all points are used if we are generating a submodel
		for (int iLev = 0; iLev < CONF_nLevs; iLev++)
			{
			unsigned ProfSpec = GetGenProfSpec(iLev, iPoint);
			int nSubProfs = nSubProfsForProfSpec(ProfSpec);
			for (int iSub = 0; iSub < nSubProfs; iSub++)
				{
				Mat *pCovar = &Stats.Covars[iLev][iPoint][iSub];

				GenAsmStats(*pCovar, Stats.AvProfs[iLev][iPoint][iSub],
							Stats.Profs[iLev][iPoint][iSub], Stats.nProfs[iLev][iPoint], ProfSpec);

				if (ProfSpec & PROF_2d)
					{
					// do matrix trimming as specified by CONF_nGenTrimCovar and CONF_GenCovarMinRatio

					double AvTrace = pCovar->trace() / pCovar->ncols();	// average value of diagonal elements
					double Min = CONF_GenCovarMinRatio * AvTrace;		// min allowed val is a fraction of average diag element

					pCovar->clearSmallElems(CONF_GenCovarMinRatio * pCovar->trace() / pCovar->ncols());

					if (CONF_nGenTrimCovar)
						{
						// trim covars that are for pixels that are more than CONF_nGenTrimCovar apart

						nTrimmedMats++;
						if (!fIterativeTrimMat(*pCovar, CONF_nGenTrimCovar, ProfSpec))
							nBadTrimmedMats++;

						// The covar matrix now has lots of 0 elements and processing will be faster
						// if we store it as a sparse matrix.
						// This also deletes any very small elems that are an artifact of fIterativeTrimMat.

						ConvertSymmetricMatToSparseFormat(*pCovar, AvTrace / 1e6);
						}
					}
				}
			}
	}
lprintf("0\n");

if (nBadTrimmedMats)
	lprintf("%d of %d trimmed matrices are not pos def\n", nBadTrimmedMats, nTrimmedMats);
}

//-----------------------------------------------------------------------------
// Write covariance matrices for all points for one level iLev

static void WriteCovarMats (const tStats &Stats, int iLev, int nPoints, const char sFile[], FILE *pAsmFile)	// in
{
for (int iPoint = 0; iPoint < nPoints; iPoint++)
	{
	int nSubProfs = nSubProfsForProfSpec(GetGenProfSpec(iLev, iPoint));
	for (int iSub = 0; iSub < nSubProfs; iSub++)
		if (Stats.nProfs[iLev][iPoint] == 0)													// not all points are used if we are generating a submodel
			Fprintf(pAsmFile, "# Covar Lev %d Point %d Sub %d\n{ 0 0 }\n", iLev, iPoint, iSub);	// generate a dummy matrix
		else
			{
			const Mat *pCovar = &Stats.Covars[iLev][iPoint][iSub];
			char sComment[SLEN]; sprintf(sComment, "Covar Lev %d Point %d Sub %d", iLev, iPoint, iSub);
			pCovar->write(sFile, pAsmFile, QUIET, sComment, "%g");
			}
	}
}

//-----------------------------------------------------------------------------
// Write profiles for all points for one level iLev

static void WriteProfiles (const tStats &Stats, int iLev, int nPoints, const char sFile[], FILE *pFile)
{
for (int iPoint = 0; iPoint < nPoints; iPoint++)
	{
	if (Stats.nProfs[iLev][iPoint] == 0)													// not all points are used if we are generating a submodel
		{
		int nSubProfs = nSubProfsForProfSpec(GetGenProfSpec(iLev, iPoint));
		for (int iSub = 0; iSub < nSubProfs; iSub++)
			Fprintf(pFile, "# Lev %d Prof %d Sub %d\n{ 0 0 }\n", iLev, iPoint, iSub);		// generate a dummy matrix
		}
	else
		{
		char sHeader[SLEN]; sprintf(sHeader, "Lev %d Prof %d Sub", iLev, iPoint);
		WriteMatVec(Stats.AvProfs[iLev][iPoint], sFile, pFile, NULL, sHeader, NULL, QUIET, "%g");
		}
	}
}

//-----------------------------------------------------------------------------
// Create a new ASM file i.e. write results to disk

static void WriteAsmFile (const char sFile[], const char sShapeFile[],  Vec &EigVals, Mat &EigVecs,		// in
						  const tStats &Stats,															// in
						  const tShapes &Shapes, char *sTagRegExp, unsigned Mask1, unsigned Mask2,		// in
						  const char sProfSpecs[], int iRefShape)
{
FILE *pFile = Fopen(sFile, "w");			// will give an msg and exit if can't open

int nPoints = Shapes.nPoints;

char drive[_MAX_DRIVE], dir[_MAX_DIR], fname[_MAX_FNAME], ext[_MAX_EXT];
_splitpath(sFile, drive, dir, fname, ext);

Fprintf(pFile, "ASM%d [%s]\n", CONF_nAsmFileVersion, fname);

// parameters that the search routines need to know about

Fprintf(pFile, "Points %d\n", nPoints);
Fprintf(pFile, "Levs %d\n", CONF_nLevs);
Fprintf(pFile, "PyramidRatio %g\n", double(CONF_PyrRatio));
Fprintf(pFile, "PyramidReduceMethod %d\n", CONF_PyrReduceMethod);
Fprintf(pFile, "ExplicitPrevNext %d\n", fgExplicitPrevNext);
Fprintf(pFile, "NormalizedProfLen %g\n", double(gNormalizedProfLen));
Fprintf(pFile, "ScaledFaceWidth %d\n", CONF_nGenPrescaleFaceWidth);
Fprintf(pFile, "BilinearRescale %d\n", CONF_fGenPrescaleBilinear);
Fprintf(pFile, "TrimCovar %d\n", CONF_nGenTrimCovar);
extern double SigmoidScaleGlobal;
Fprintf(pFile, "SigmoidScale %g\n", SigmoidScaleGlobal);
Fprintf(pFile, "\n");

time_t ltime; time(&ltime);
char *sTime = ctime(&ltime);
sTime[strlen(sTime)-1] = 0;	// remove \n
Fprintf(pFile, "# Generated by masm [%s]\n", sTime);	// put time in [] so it is ignored by mdiff.exe

// parameters from the masm command line

Fprintf(pFile, "# nShapes %d\n", Shapes.nShapes);
Fprintf(pFile, "# ShapeFile %s\n", sShapeFile);
Fprintf(pFile, "# TagRegExp \"%s\"\n", sTagRegExp);
Fprintf(pFile, "# AttrMask1 %4.4x [%s]\n", Mask1, sGetAtFaceString(Mask1, false));
Fprintf(pFile, "# AttrMask2 %4.4x [%s]\n", Mask2, sGetAtFaceString(Mask2, false));
Fprintf(pFile, "# ProfSpecs %s\n", sProfSpecs);

// parameters from the shape file

Fprintf(pFile, "# ImageDirs %s\n", Shapes.sImageDirs);

// parameters generated by masm

Fprintf(pFile, "# iRefShape %d\n", iRefShape);		// shape found by iGetRefShapeIndex

// parameters from masmconf.hpp, in alphabetical order

Fprintf(pFile, "# fAlignImageToAlignedShape %d\n", CONF_fAlignImageToAlignedShape);
Fprintf(pFile, "# fDropEigToPreventSingCovar %d\n", CONF_fDropEigToPreventSingCovar);
Fprintf(pFile, "# fDropProfToPreventSingCovar %d\n", CONF_fDropProfToPreventSingCovar);
Fprintf(pFile, "# fEqualizeImages %d\n", CONF_fEqualizeImages);
Fprintf(pFile, "# fMasmGaussianNoise %d\n", CONF_fMasmGaussianNoise);
Fprintf(pFile, "# nRandomSeed %d\n", CONF_nRandomSeed);
Fprintf(pFile, "# fRecentralizeShapes %d\n", CONF_fRecentralizeShapes);
Fprintf(pFile, "# fShuffleLandmarks %d\n", CONF_fShuffleLandmarks);
Fprintf(pFile, "# fStraightenBeforeAligning %d\n", CONF_fStraightenBeforeAligning);
Fprintf(pFile, "# fUseOnlyUnobscuredFeats %d\n", CONF_fUseOnlyUnobscuredFeats);
Fprintf(pFile, "# GenSubModel %d\n", CONF_GenSubModel);
Fprintf(pFile, "# MaxEigRatio %g\n", (double)CONF_MaxEigRatio);
Fprintf(pFile, "# MaxShapeAlignDist %g\n", (double)CONF_MaxShapeAlignDist);
Fprintf(pFile, "# nGenLandmarksUsedInternally %d\n", CONF_nGenLandmarksUsedInternally);
Fprintf(pFile, "# nProfWidth1d %d\n", CONF_nProfWidth1d);
Fprintf(pFile, "# nProfWidth2d %d\n", CONF_nProfWidth2d);
Fprintf(pFile, "# nRandSeedForShuffle %d\n", CONF_nRandSeedForShuffle);
Fprintf(pFile, "# nReducedNbrPoints %d\n", CONF_nReducedNbrPoints);
Fprintf(pFile, "# nSynthesizedEyePoints %d\n", CONF_nSynthesizedEyePoints);
Fprintf(pFile, "# RidgeLambdaForSingCovars %g\n", (double)CONF_RidgeLambdaForSingCovars);
Fprintf(pFile, "# ShapeNoise %g\n", (double)CONF_ShapeNoise);
Fprintf(pFile, "# xRhsStretchShape %g\n", (double)CONF_xRhsStretchShape);
Fprintf(pFile, "# xStretchShape %g\n", (double)CONF_xStretchShape);
Fprintf(pFile, "\n");

EigVals.write(sFile, pFile, QUIET, "EigVals");
EigVecs.write(sFile, pFile, QUIET, "EigVecs");
Shapes.AvShape.write(sFile, pFile, QUIET, "AvShape");

lprintf("Checking stats and writing %s ", sFile);
for (int iLev = CONF_nLevs-1; iLev >= 0; iLev--)
	{
	Fprintf(pFile, "Lev %d FirstShape xExtent %d yExtent %d\n",
		iLev,
		int(xShapeExtent(Shapes.OrgShapes[0])/(iLev+1)),
		int(yShapeExtent(Shapes.OrgShapes[0])/(iLev+1)));

	Fprintf(pFile, "# Lev %d ProfSpecs\n", iLev, nPoints);
	WriteProfSpecs(iLev, nPoints, pFile);

	WriteCovarMats(Stats, iLev, nPoints, sFile, pFile);
	WriteProfiles(Stats, iLev, nPoints, sFile, pFile);