// $masm\genfeats.cpp 1.5 milbo$ print search features for pos analysis with R// Warning: this is raw research code -- expect it to be quite messy.// milbo petaluma oct 06//-----------------------------------------------------------------------------// This program is free software; you can redistribute it and/or modify// it under the terms of the GNU General Public License as published by// the Free Software Foundation; either version 2 of the License, or// (at your option) any later version.//// This program is distributed in the hope that it will be useful,// but WITHOUT ANY WARRANTY; without even the implied warranty of// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the// GNU General Public License for more details.//// A copy of the GNU General Public License is available at// http://www.r-project.org/Licenses///-----------------------------------------------------------------------------#include "all.hpp"static FILE *pgFeatFile;		// -G flagstatic char sgFeatFile[SLEN];#define MAX_OFFSET 10#define MAX_PROF   20typedef struct tFeatData	{	bool fUsed;	double Fit;	int iMax;				// posn of biggest feature	double Max;				// value of biggset feature	Mat Prof[1];			// copy of profile	}tFeatData;static tFeatData gFeats[1 + 2 * MAX_OFFSET];const int igFeatPoint = 30;#define TWOD 1//-----------------------------------------------------------------------------bool fGenFeats (void){return pgFeatFile != 0;}//-----------------------------------------------------------------------------void InitGenFeat (bool fGenFeatsFlag, tAsmModel &Model, const char sRefShapeFile[], int nModels){if (!fGenFeatsFlag)	return;sprintf(sgFeatFile, "out/feat%2.2d.tab", igFeatPoint);lprintf("Opening %s\n", sgFeatFile);pgFeatFile = Fopen(sgFeatFile, "w");Model.StartMethod = SM_UseRefShapeAndAlignedMeanShape;ASSERT(sRefShapeFile[0]);ASSERT(nModels == 1);ASSERT(CONF_xRefShapeOffset == 0 && CONF_yRefShapeOffset == 0);for (int i = 0; i < 1 + 2 * MAX_OFFSET; i++)	gFeats[i].fUsed = false;}//-----------------------------------------------------------------------------void AccumulateFeats (const int ix, const int iy, int iSub, double Fit, const Vec &Prof, unsigned ProfSpec){if (!pgFeatFile)	return;if (iy != 0)	return;if (ix != 0)	return;}//-----------------------------------------------------------------------------static void GenFeats1 (bool fTwoD,				SHAPE &Shape,				tSearchImages &SearchImgs,				const tAsmLev *pAsmLev,				int iLev, const tLand Lands[],				const int iPoint,				int nPixSearch,				int ix, int iy){SHAPE AlignedMeanShape(Shape);	// TODO probably unneededconst unsigned ProfSpec = (fTwoD? pAsmLev->ProfSpecs[iPoint]: (PROF_Grad|PROF_Flat));ASSERT(nSubProfsForProfSpec(ProfSpec) == 1);const int iSub = 0;const int nProfWidth = nGetProfWidthFromModel(iPoint, iSub, *pAsmLev);int nelems = nProfWidth * (fTwoD? nProfWidth: 1);Vec Prof(nelems, ROWVEC);if (fTwoD)	Get2dProf(Prof, ProfSpec, iSub, SearchImgs.Img, SearchImgs.Grads, Shape, iPoint, ix, iy, nProfWidth, false);else	{	PrepareProf1D(SearchImgs.Img, Shape, ProfSpec, Lands, AlignedMeanShape, iPoint, nProfWidth + 2 * nPixSearch, 0, 0);	Get1dProf(Prof, ProfSpec, SearchImgs.Img, iPoint, ix, 0, false);	}Vec Diag(nProfWidth);int iStep = (fTwoD? nProfWidth: 1);for (int i = 0; i < nProfWidth; i++)	Diag(i) = Prof(i * iStep);// print profile statsint iMaxAbs, iMaxNeg, iMaxPos;double MaxAbs = -FLT_MAX, MaxPos = -FLT_MAX, MaxNeg = FLT_MAX;for (i = 0; i < nProfWidth; i++)	{	double x = Diag(i);	if (fabs(x) > MaxAbs)		{		MaxAbs = fabs(x);		iMaxAbs = i;		}	if (x > MaxPos)		{		MaxPos = x;		iMaxPos = i;		}	if (x < MaxNeg)		{		MaxNeg = x;		iMaxNeg = i;		}	}Fprintf(pgFeatFile, "% 7d %7.2f % 7d %7.2f % 7d %7.2f ", iMaxAbs, MaxAbs, iMaxNeg, MaxNeg, iMaxPos, MaxPos);#if TWODProf.viewAsRow().print(NULL, "%7.2f", sgFeatFile, pgFeatFile);#elseDiag.t().print(NULL, "%7.2f", sgFeatFile, pgFeatFile);#endif}//-----------------------------------------------------------------------------// Generate features for igFeatPoint for the given Shape and image// If -G flag not used, then this just returns// If iLev!=0 then this just returnsvoid GenFeats (SHAPE &Shape,										// io				tSearchImages &SearchImgs,							// in				const tAsmModel &Model, 							// in				int iLev, const tLand Lands[],				const char sImageBase[], int nStartLev)				// in{if (!pgFeatFile)	return;if (iLev != 0)	return;SHAPE AlignedMeanShape(Shape);	// TODO probably unneededint nUsedPoints = nGetNbrUsedPoints(Shape, Shape.nrows());if (fgExplicitPrevNext && (nUsedPoints != Model.nPoints))		// AlignedMeanShape is only needed if don't have all landmarks	SysErr("GenFeats");static bool fWriteHeader = true;if (fWriteHeader)	{	fWriteHeader = false;	// xp is predicted offset, y is actual offset	Fprintf(pgFeatFile, " y      xp    ");	Fprintf(pgFeatFile, "Fit iMaxAbs  MaxAbs iMaxNeg  MaxNeg iMaxPos  MaxPos    ");#if TWOD	for (int i = -169/2; i < 0; i++)		Fprintf(pgFeatFile, "P_%-3d  ", -i);	for (i = 0; i < 169/2+1; i++)		Fprintf(pgFeatFile, "P%-3d   ", i);	Fprintf(pgFeatFile, "\n");#else	Fprintf(pgFeatFile, "P_6    P_5    P_4    P_3    P_2    P_1     P0     P1     P2     P3     P4     P5     P6\n");#endif	}double BestFit = DBL_MAX;int ixBest = 0;int iy = 0;int iPoint = igFeatPoint;const tAsmLev *pAsmLev = &Model.AsmLevs[iLev];const unsigned ProfSpec = pAsmLev->ProfSpecs[iPoint];const bool fTwoD = (pAsmLev->ProfSpecs[iPoint] & PROF_2d) != 0;ASSERT(fTwoD);const int nPixSearch = Model.nPixSearch;for (int ix = -nPixSearch; ix <= nPixSearch; ix++)		// 1D: offset along whisker 2D: offset horizontally	{	// call model-specific fitness routine	double Fit = GetProfileFit(SearchImgs, iPoint, ix, iy, *pAsmLev, Shape, Lands, ProfSpec, Model.nTrimCovar);	// Test for a new best fit	// We test against ix<=0 and iy<=0 so if there are exact matches then ixBest=0 (no change)	if ((ix <= 0) && (iy <= 0)? Fit <= BestFit:  Fit < BestFit)		{		ixBest = ix;		BestFit = Fit;		}	}// Cycle through estimated offsets.  Do it this way because for// a linear model, we want no noise on xp (predictor), although the y (response) can be noisystatic int iOffset;if (++iOffset > nPixSearch)	iOffset = -nPixSearch;Fprintf(pgFeatFile, "% 2d      % 2d %6.2f ", iOffset-ixBest, iOffset, BestFit);	// xp xe FitGenFeats1(TWOD, Shape, SearchImgs, pAsmLev, iLev, Lands, iPoint, nPixSearch, iOffset, iy);fflush(pgFeatFile);	//TODO}