#include "Graph.h"

void Graph::loadMeshes(char* fileFormat, bool dynamicProgramming, int version, bool spectralInitToMatchMore)
{
	//initializes all meshes to be used in collection

	if (! spectralInitToMatchMore && (N < 5 || N > 9))
	{
		cout << "\n\n\n\nWARNING: # extremes must be in [5, 9] for efficiency and/or accuracy reasons (our original version at coarse level)\n\n\n\n";
//		exit(0);
	}
	if (spectralInitToMatchMore && N > 500)
	{
		cout << "\n\n\n\nWARNING: # samples must less than 500 for accuracy reasons (our dense matching extension)\n\n\n\n";
//		exit(0);
	}
	cout << "................LOADING MESHES OF THE COLLECTION................\n";
	nMeshes = stopIdx - startIdx + 1; //# of meshes in use, start and stop inclusive
	if (nMeshes < 2)
	{
		cout << "# of meshes in collection must be greater than or equal to 2 but now " << nMeshes << endl;
		exit(0);
	}
	meshes = new Mesh*[nMeshes];
	outFiles = new FILE*[nMeshes-1];
	char fName[250], algoName[250], dName[250], orderingName[250] = "", vName[50], rName[50];
	if (dynamicProgramming && version == 0)
	{
		sprintf(vName, "EXPONENTIAL_ALGO\\");
//		topK = (nMeshes > 6 ? 8 : topK); //to make the exhaustive search tractable keep topK tiny no matter what constructor says if 7 or 8 shapes will be matched; not in this public release
		if (nMeshes > 8)
		{
			cout << "# of meshes should be less than 9 to make the exhaustive search tractable\n";
			exit(0);
		}
	}
	else
		//sprintf(vName, "version%d\\", version);
		sprintf(vName, "DYNAMIC_PROGRAMMING_ALGO\\"); //no version stuff in this public release
	sprintf(algoName, (dynamicProgramming ? vName : "SHORTEST_PATH_ALGO\\"));

	grdTruth = true;
	for (int m = 0; m < nMeshes; m++)
	{
		meshes[m] = new Mesh(startIdx+m, fileFormat, m);
//note that dName and orderingName are unused in this public release
		if (meshes[m]->theInput == SCAPE)
			sprintf(dName, "SCAPE\\");
		else if (meshes[m]->theInput == TOSCA_CAT)
			sprintf(dName, "TOSCA\\cat\\");
		else if (meshes[m]->theInput == TOSCA_CENTAUR)
			sprintf(dName, "TOSCA\\centaur\\");
		else if (meshes[m]->theInput == TOSCA_DAVID)
			sprintf(dName, "TOSCA\\david\\");
		else if (meshes[m]->theInput == TOSCA_DOG)
			sprintf(dName, "TOSCA\\dog\\");
		else if (meshes[m]->theInput == TOSCA_GORILLA)
			sprintf(dName, "TOSCA\\gorilla\\");
		else if (meshes[m]->theInput == TOSCA_HORSE)
			sprintf(dName, "TOSCA\\horse\\");
		else if (meshes[m]->theInput == TOSCA_VICTORIA)
			sprintf(dName, "TOSCA\\victoria\\");
		else if (meshes[m]->theInput == TOSCA_WOLF)
			sprintf(dName, "TOSCA\\wolf\\");
		else if (meshes[m]->theInput == TOSCA_DAVIDVICTORIA)
		{
			sprintf(dName, "TOSCA\\david-victoria\\");
			if (startIdx != 0 || stopIdx != 13)
			{
				cout << "for ground-truth evaluations, david-victoria pair must be loaded from startIdx=0 to stopIdx=13\n";
				exit(0);
			}
		}
		if (meshes[m]->ordering > 0)
			sprintf(orderingName, "another order%d\\", meshes[m]->ordering);
		sprintf(rName, spectralInitToMatchMore ? "dense\\" : "coarse\\"); //resolution name
		if (m > 0) //2 consecutive meshes achieved, so set the output file for their map
		{
			//sprintf(fName, "output\\%s%sN=%d\\%s%s%d-%d.dat", algoName, dName, N, orderingName, rName, startIdx+m-1, startIdx+m);
			sprintf(fName, "output\\%s\\%s%d-%d.dat", algoName, rName, startIdx+m-1, startIdx+m); //public release format
			outFiles[m-1] = fopen(fName, "w");
		}
	}
}
void Graph::loadMeshTriplets(char* fileFormat, bool dynamicProgramming, int version, bool spectralInitToMatchMore)
{
	//initializes 3 meshes to be used in collection

	cout << "loadMeshTriplets() not made ready for public release; play with file naming stuff here to integrate this which is easy and similar to loadMeshTriplets() but requires recompilation\n";
	exit(0);

	if (N < 5 || N > 9)
	{
		cout << "# extremes must be in [5, 9] for efficiency and/or accuracy reasons\n";
		exit(0);
	}
	cout << "................LOADING TRIPLET OF MESHES................\n";
	nMeshes = 3;
	meshes = new Mesh*[nMeshes];
	outFiles = new FILE*[nMeshes-1];
	char fName[250], algoName[250], dName[250], orderingName[250] = "", vName[50], rName[50];
	sprintf(vName, "version%d\\", version);
	sprintf(algoName, (dynamicProgramming ? vName : "SHORTEST_PATH_ALGO\\"));

	grdTruth = true;
	for (int m = 0; m < nMeshes; m++)
	{
		int idx;
		if (m == 0)
			idx = startIdx;
		else if (m == 1)
			idx = middleIdx;
		else
			idx = stopIdx;
		meshes[m] = new Mesh(idx, fileFormat, m);
		if (meshes[m]->theInput == SCAPE)
			sprintf(dName, "SCAPE\\");
		else if (meshes[m]->theInput == TOSCA_CAT)
			sprintf(dName, "TOSCA\\cat\\");
		else if (meshes[m]->theInput == TOSCA_CENTAUR)
			sprintf(dName, "TOSCA\\centaur\\");
		else if (meshes[m]->theInput == TOSCA_DAVID)
			sprintf(dName, "TOSCA\\david\\");
		else if (meshes[m]->theInput == TOSCA_DOG)
			sprintf(dName, "TOSCA\\dog\\");
		else if (meshes[m]->theInput == TOSCA_GORILLA)
			sprintf(dName, "TOSCA\\gorilla\\");
		else if (meshes[m]->theInput == TOSCA_HORSE)
			sprintf(dName, "TOSCA\\horse\\");
		else if (meshes[m]->theInput == TOSCA_VICTORIA)
			sprintf(dName, "TOSCA\\victoria\\");
		else if (meshes[m]->theInput == TOSCA_WOLF)
			sprintf(dName, "TOSCA\\wolf\\");
		else if (meshes[m]->theInput == TOSCA_DAVIDVICTORIA)
		{
			sprintf(dName, "TOSCA\\david-victoria\\");
			if (startIdx != 0 || stopIdx != 13)
			{
				cout << "for ground-truth evaluations, david-victoria pair must be loaded from startIdx=0 to stopIdx=13\n";
				exit(0);
			}
		}
		if (meshes[m]->ordering > 0)
			sprintf(orderingName, "another order%d\\", meshes[m]->ordering);
		sprintf(rName, spectralInitToMatchMore ? "dense\\" : ""); //resolution name
		if (m > 0) //2 consecutive meshes achieved, so set the output file for their map
		{
			sprintf(fName, "output\\%s%sN=%d\\%striplets\\%s%d-%d-%d %d-%d.dat", algoName, dName, N, orderingName, rName, startIdx, middleIdx, stopIdx, (m == 1 ? startIdx : middleIdx), idx);
			sprintf(fName, "output\\%s%sN=%d\\%s%s%d-%d.dat", algoName, dName, N, orderingName, rName, startIdx+m-1, startIdx+m);
			outFiles[m-1] = fopen(fName, "w");
		}
	}
}
/*void Graph::deleteMeshTriplets()
{
	//for parallel runs, accumulated mesh info (12K vertices etc) fills up memory

	for (int m = 0; m < nMeshes; m++)
	{
		meshes[m]->verts.clear();
		meshes[m]->tris.clear();
		meshes[m]->edges.clear();
		delete meshes[m];		
	}
}*/

void Graph::addGVertices(int m1, int m2, bool spectralInitToMatchMore)
{
	//add all graph vertices concerning extremity maps b/w meshes[m1] & meshes[m2]

	//maps and their distortions sorted w.r.t. distortions
	OutputPkt allMaps;
	if (spectralInitToMatchMore)
	{
		//N evenly-spaced samples to be matched (N < 500)
//		N = 80; //80 used in our PAMI spectral initialization
		int D = 6; //in D-dimensional spectral domain
		meshes[m1]->computeEvenlySpacedSamples(N, D);
		meshes[m2]->computeEvenlySpacedSamples(N, D);
		//min-weight perfect matching (blossom) b/w samples gives each graph vertex as a mapping along w/ its distortion, so gvert is a struct keeping mapping and its distortion
		allMaps = meshes[m1]->spectralMatchWith(meshes[m2], N, D, topK);
	}
	else
	{
		//N extreme vertices to be matched (N < 10)
		meshes[m1]->computeExtremes(N);
		meshes[m2]->computeExtremes(N);
		//permutations b/w extremes gives each graph vertex as a mapping along w/ its distortion, so gvert is a struct keeping mapping and its distortion
		allMaps = meshes[m1]->matchWith(meshes[m2], N, topK);
	}

	for (int i = 0; i < topK; i++)
	{
		gverts.push_back( new GVertex(allMaps.candidateMapsSorted[i], allMaps.distortionsSorted[i], m1, m2) );
//cout << "gverts[" << (int) gverts.size()-1 << "] added: " << allMaps.distortionsSorted[i] << endl;
	}

	//on 12.5K scape w/ N=9, 1 call to computeExtremes() a sec, to matchWith() ~2a secs where a = ~1; overall build() (over 72 shapes) 134 secs, computePaths2() for dynamic programming 15 secs, 149 in total
	//on 50K tosca_david w/ N=9, 1 call to computeExtremes() ~3a sec, to matchWith() a secs where a = ~1; overall build() (over 7 shapes) 31 secs, computePaths2() for dynamic programming 1 secs, 32 in total
	//																				 70a secs if i don't do the tooBigDistortion trick, i.e. sort all N! instad of just potentially good maps which sorts a constant amount of shapes <1000
}
void Graph::build(bool dynamicProgramming, bool spectralInitToMatchMore)
{
	//builds the directed graph on which multiple correspondences will be obtained via shortest paths algorithm; just need gverts for DYNAMIC_PROGRAMMING_ALGO*

	cout << "\n\n................BUILDING THE GRAPH OVER THE COLLECTION................\n";

	//the maps will always be NxN one-to-one b/w N=7 extremes; this makes N!^^ * (nMeshes-1) nodes in this graph each representing a map b/w 2 consecutive meshes
	//^^ since some maps are simply too bad, i.e. the ones w/ distortions above average of disotortions or the ones w/ distortions above the topK least distortions, do not use all N! maps in vain; use only the topK ones which replaces N! w/ a constant topK = typically 20

	cout << "graph vertices..\n";
	for (int m = 0; m < nMeshes-1; m++) //maps b/w 0-1, 1-2, 2-3, 3-4, .., (nMeshes-2)-(nMeshes-1)
		addGVertices(m, m+1, spectralInitToMatchMore);

	if (dynamicProgramming)
		return; //just need gverts in this case as the connections b/w them'll be implicitly achieved later by filling the dynamic programming table

	if (nMeshes > 2) //multiple shape correspondence problem that needs this graph
	{
		cout << "graph edges..\n";
		float weight;
		//connect gverts layer by layer via graph edges; every topK (not N! 'cos for efficiency topK mappings b/w N extremes used only) gverts make 1 layer and layers go from left to right; i'll connect all gverts of layer_i w/ all verts of layer_i+1, i.e. complete bipartite graph
		for (int layer = 1; layer < nMeshes-1; layer++) //nMeshes-1 vertical layers for nMeshes meshes as each gvert in a layer is related to a map b/w 2 consecutive meshes, e.g. for nMeshes=3, there'll be 2 layers, and just 1 edge set that connects layer1 and 2
			for (int gv = 0; gv < topK; gv++) //topK gverts at each layer
			{
				int gv1 = (layer-1)*topK + gv;
				for (int ggvv = 0; ggvv < topK; ggvv++) //connect gv1 from layer_i to all gverts in layer_i+1 (layer-1 just-above becomes layer just-below for i+1 effect)
				{
					int gv2 = layer*topK + ggvv;
#ifdef SIMPLE_WEIGHTING
					weight = (gverts[gv1]->distortion + gverts[gv2]->distortion) / 2;
#else
					weight = (gverts[gv1]->distortion + gverts[gv2]->distortion) / 2; //something more complicated than this!!!!!!!!!!!!!
#endif
//float weight = (gv1 == 60 && gv2 == 80 ? 0.99f : (gverts[gv1]->distortion + gverts[gv2]->distortion) / 2);
					gedges.push_back( new GEdge(gv1, gv2, weight) ); //weight of this edge is the avg distortion of the connected maps (represented by gverts); averaging is redundant 'cos div by 2 on all edges does not change shortest paths result
//if (gv1==40)cout << "ggggggggggggedges[" << (int) gedges.size()-1 << "] added: " << gv1 << "->" << gv2 << "\t" << (gverts[gv1]->distortion + gverts[gv2]->distortion) / 2 << endl;

//use max() instead of average in 3rd argument to GEdge??????????????????????????????????????????????????????????

					
					//this newcomer directed edge is from gv1 to gv2, so incident to gv1 only; by simply setting gv2.gedgeList.push_back((int) gedges.size() - 1) the whole directed graph becomes undirected, but i need directed graph here
					gverts[gv1]->gedgeList.push_back((int) gedges.size() - 1);
				}
			}
		//there's also a virtual start graph vertex sgv from which solve()'ll start dijkstra; sgv is connected to the first layer (layer1) w/ the same weights of say 0
		vector< int > dummy;
		int sgvIdx = (int) gverts.size();
		gverts.push_back( new GVertex(dummy, -1.0f, -1, -1) );
		for (int gv = 0; gv < topK; gv++) //sgvIdx's connections to the first layer, i.e. the first topK gverts
		{
			gedges.push_back( new GEdge(sgvIdx, gv, 0.0f) );
			gverts[sgvIdx]->gedgeList.push_back((int) gedges.size() - 1); //edge added just-above ((int) gedges.size() - 1) is incident to sgvIdx as it is directed from sgvIdx to gv
		}
		//finally, there's also a virtual target graph vertex tgv to which shortest path from sgv will go; the last layer (layer_nMeshes-1) is connected to tgv w/ the same weights of say 0
		int tgvIdx = (int) gverts.size();
		gverts.push_back( new GVertex(dummy, -1.0f, -1, -1) );
		for (int gv = 0; gv < topK; gv++) //last layer's connections to tgvIdx
		{
			gedges.push_back( new GEdge((nMeshes-2)*topK + gv, tgvIdx, 0.0f) );
			gverts[ (nMeshes-2)*topK + gv ]->gedgeList.push_back((int) gedges.size() - 1); //edge added just-above ((int) gedges.size() - 1) is incident to this last layer graph vert ((nMeshes-2)*topK + gv) as it is directed from (nMeshes-2)*topK + gv to tgvIdx
		}
	}
	else //conventional pairwise shape correspondence problem that does not even need this graph; simply select the min-distortion map b/w available 2 meshes
	{
		//select the min-distortion map
		float minDisto = INF;
		int minDistoMapIdx;
		for (int map = 0; map < (int) gverts.size(); map++)
			if (gverts[map]->distortion < minDisto)
			{
				minDisto = gverts[map]->distortion;
				minDistoMapIdx = map;
			}
		//result to file and cmd prompt
		cout << "Optimal min-distortion pairwise map b/w 2 meshes written to \\output folder and here:\n";
		printMap(0, minDistoMapIdx, true);		
	}
}

bool Graph::solve()
{
	//w/ shortest paths solution, finds the optimal path from the first mesh to the last, which implies intermediate paths and hence maps as well, where optimality is defined as minimization of sum of pairwise distortions on the final path b/w the first and the last mesh (so effect of simply selecting consecutive min-distortion maps independently)

	if (nMeshes == 2)
	{
		cout << "\n\nNo graph solving 'cos this is not a multiple correspondence instance\n\n"; //build() created and fprinted that optimal map
		return false;
	}
	cout << "\n\n................SOLVING THE GRAPH OVER THE COLLECTION................\n";

	//original Dijkstra to compute shortest paths from sourceVert to all others; i'm interested in the shortest path from sgv (2nd from last in gvert) to tgv (last in gverts)	
	int sourceVert = (int) gverts.size() - 2, targetVert = (int) gverts.size() - 1; //the 2nd from last is the virtual sgv is the source of our shortest paths computation; recall last is the virtual tgv
	gverts[sourceVert]->dSpanningUnn = new float[ (int) gverts.size() ]; //shortest distances from sourceVert to all other gverts
	FibHeap* h = new FibHeap(); //whole dijkstra w/ Fibonacci heap takes O(VlgV + e); w/ minHeap takes O(VlgV + elgV), w/ simple array takes O(V^2)	

	//initialization
	for (int v = 0; v < (int) gverts.size(); v++)
	{
		gverts[v]->dSpanned = (v == sourceVert ? 0.0f : INF);
		FibHeapNode* hn = new FibHeapNode(gverts[v]->dSpanned, v); //v is a handle from headNode to verts,
		gverts[v]->heapNode = hn;								  //and hn is an handle from verts to heapNode
		h->Insert(hn);
		gverts[v]->predecessor = -1;
		//no delete hn; here 'cos it'll be extracted, used, and deleted below
	}

	while (true)
	{
		//extract vert w/ min dSpanned value
		FibHeapNode* minNode = h->ExtractMin(); //fast O(lgV) priority queue algo
		if (minNode == NULL)
			break; //no more verts w/ unknown shortest dist; so, we're done
		int minDidx = minNode->element;
		float minD = minNode->key;
		gverts[sourceVert]->dSpanningUnn[minDidx] = minD; //distance of the shortest path from sourceVert to minDidx
		delete minNode; //new hn above is deleted here to save memo [info: delete minNode; automatically calls minNode->~FibHeapNode(), i.e. its destructor]

		//relax each edge incident to minDidx
		for (int ve = 0; ve < (int) gverts[minDidx]->gedgeList.size(); ve++)
		{
			int e = gverts[minDidx]->gedgeList[ve];
			int va = gedges[e]->gv2;
			if (gverts[minDidx]->dSpanned + gedges[e]->length < gverts[va]->dSpanned)
			{
				gverts[va]->dSpanned = gverts[minDidx]->dSpanned + gedges[e]->length; //relaxation
				h->DecreaseKey(gverts[va]->heapNode, gverts[va]->dSpanned); //worst O(lgV), amortized O(1)
				gverts[va]->predecessor = minDidx;
			}
		}
	} //end of while (true)
	delete h;

	//gverts[sourceVert]->dSpanningUnn[targetVert] and totalDistortion not necessarily the same 'cos dSpanningUnn considers edge weights which are averages of the endpnt maps, whereas totalDistortion considers distortions of the selected mappings hence to be compared w/ DYNAMIC_PROGRAMMING_ALGO* results
	cout << "distance of the computed shortest path (recall edges weighted by average distortions of endpoints): " << gverts[sourceVert]->dSpanningUnn[targetVert] << endl
		 << "path itself (given as id's of the graph vertices in the path):\n";
	//the chain of predecessors originating at a vertex v (= tgvIdx in our case) runs backwards along a shortest path from source vertex (sgvIdx) to v; use recursion to extract this path
	extractMapsOnTheShortestPath(sourceVert, targetVert, targetVert);
	//fprint the maps
	float totalDistortion = 0.0f, shoPathGrdTruthDistortion = 0.0f;
	int optPathSymFlips = 0;
	for (int i = 0; i < (int) mapsOnTheOptimalPath.size(); i++)
	{		
		//ground-truth stuff for the next 2 lines
		gverts[ mapsOnTheOptimalPath[i] ]->grdTruthDistortion = (grdTruth ? getGrdTruthDistortion(gverts[ mapsOnTheOptimalPath[i] ]->sidx, gverts[ mapsOnTheOptimalPath[i] ]->tidx, gverts[ mapsOnTheOptimalPath[i] ]->map) : 0.0f);
		shoPathGrdTruthDistortion += gverts[ mapsOnTheOptimalPath[i] ]->grdTruthDistortion;
		optPathSymFlips += (symFlipped(gverts[ mapsOnTheOptimalPath[i] ]->grdTruthDistortion) ? 1 : 0);
		printMap(i, mapsOnTheOptimalPath[i]);
		totalDistortion += gverts[ mapsOnTheOptimalPath[i] ]->distortion;
	}
	cout << "\ntotal & average isometric distortion on the shortest path: " << totalDistortion << " & " << totalDistortion / (int) mapsOnTheOptimalPath.size() << endl
		 << "\ntotal & average ground-truth distortion on this shortest path: " << shoPathGrdTruthDistortion << " & " << shoPathGrdTruthDistortion / (int) mapsOnTheOptimalPath.size() << endl;

	////////////////////////// composition of all pairwise maps for distortionSum computation //////////////////////////
	//heavily adapted from findTheOptimalPath() to generate total or average distortion over all (shortest path above + compositions) pairwise maps inplied by the mapsOnTheOptimalPath = shortest path; then i can compare this distortionSumMin vs. findTheOptimalPath().distortionSumMin
	//mapsOnTheOptimalPath is a vector/collection of maps
if ((int) mapsOnTheOptimalPath.size() != nMeshes-1) {cout << "# maps in the shortest path wrong\n"; exit(0);}
	//sum of isometric distortions imposed by the path of maps: mapsOnTheOptimalPath[0], mapsOnTheOptimalPath[1], .., mapsOnTheOptimalPath[nMeshes-2]
	float distortionSum = 0.0f, overallGrdTruthDistortion = 0.0f, tmp;
	int nAdds = 0;
	overallSymFlips = 0;
	for (int target = nMeshes-1; target > 0; target--)
		/////////////// 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1] ////////////////////////
		for (int src = 0; src < target; src++) //isometric distortions to be added will be b/w meshes 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1]
		{
			vector< int > composedMap;
			//compose the map b/w src and target
			for (int match = 0; match < N; match++) //for each of N matches in gverts[ path->maps[src] ]->map
			{
				int i = gverts[ mapsOnTheOptimalPath[src] ]->map[2*match], ci = gverts[ mapsOnTheOptimalPath[src] ]->map[2*match + 1]; //match i-ci where i from meshes[src], ci from meshes[src+1]
				for (int a = src+1; a < target; a++) //propagate ci info to the map at right
				{
					//path->maps[a] defined for all a values here unlike computePaths2() which requires a replacement w/ tryMap at the last iteration
					//look for the match ci-cii in the map at right, namely gverts[ path->maps[a] ]->map
					int cii = -1;
					for (int match2 = 0; match2 < N; match2++) //for each of N matches in gverts[ path->maps[a] ]->map
						if (gverts[ mapsOnTheOptimalPath[a] ]->map[2*match2] == ci)
						{
							cii = gverts[ mapsOnTheOptimalPath[a] ]->map[2*match2 + 1];
							break;
						}
//if (cii == -1) {cout << "cii cannot remain -1 during map composition\n"; exit(0);}
					ci = cii; //propagate new ci info to the map at right
				}
				//i-ci as the (match+1)'th match of composedMap b/w meshes[src]-meshes[target]
				composedMap.push_back(i);
				composedMap.push_back(ci);
			} //end of match
			//isometric distortion of the composed map b/w meshes[src]-meshes[target]
			distortionSum += meshes[src]->getCorrespListDistortion(meshes[target], composedMap);
			nAdds++;
			tmp = (grdTruth ? getGrdTruthDistortion(src, target, composedMap) : 0.0f); //meshes[src]->getCorrespListGrdTruthDistortion(composedMap) call makes unresolved external symbol error for no reason, very weird
			overallGrdTruthDistortion += tmp;
			overallSymFlips += (symFlipped(tmp) ? 1 : 0);
		} //end of src
	////////////////////////// composition of all pairwise maps for distortionSum computation ends //////////////////////////

	//distortionSumMin is actually minimized and can be compared w/ solve().distortionSumMin and we expect distortionSumMin here to be better/smaller
	cout << "\ntotal & average isometric distortion over ALL (shortest path + compositions) pairwise distortions: " << distortionSum << " & " << distortionSum / nAdds << " (for C(" << nMeshes << ", 2) = " << nAdds << " pairwise maps)\n\n"
		 << "\ntotal & average ground-truth distortion over ALL (shortest path + compositions) pairwise distortions: " << overallGrdTruthDistortion << " & " << overallGrdTruthDistortion / nAdds << " (for C(" << nMeshes << ", 2) = " << nAdds << " pairwise maps)\n\n"
		 << "\nmaps written to \\output\\SHORTEST_PATH_ALGO folder\n\n"
		 << ">>>>>> SUMMARY (optimal path isometric distortion) & (over ALL isometric distortions): " << totalDistortion / (int) mapsOnTheOptimalPath.size() << " & " << distortionSum / nAdds << "\n"
		 << ">>>>>> SUMMARY (optimal path ground-truth distortion) & (over ALL ground-truth distortions): " << shoPathGrdTruthDistortion / (int) mapsOnTheOptimalPath.size() << " & " << overallGrdTruthDistortion / nAdds << "\n"
		 << ">>>>>> SUMMARY (optimal path # symmetric flips) & (over ALL # symmetric flips): " << optPathSymFlips << " / " << (int) mapsOnTheOptimalPath.size() << " & " << overallSymFlips << " / " << nAdds << "\n\n";

//print specific map, 2nd argument (called argument in fnctn call, parameter in fnctn declaration) on some file pointed by new outFiles[0]
//outFiles[0] = fopen("output\\testMap.dat", "w"); printMap(0, 64);

	//to be used for triple runs only, i.e. if MESH_TRIPLETS_ONLY is defined
	overallIsoDistOfTriplet = distortionSum / nAdds;
	overallGrdTruthDistOfTriplet = overallGrdTruthDistortion / nAdds;

	return true;
}
void Graph::extractMapsOnTheShortestPath(int virtualSource, int v, int virtualTarget)
{
	if (v == virtualSource)
		;//cout << v << "\t"; //redundant to keep virtual source graph vertex's dummy map in mapsOnTheOptimalPath or print it
	else if (gverts[v]->predecessor != -1)
	{
		extractMapsOnTheShortestPath(virtualSource, gverts[v]->predecessor, virtualTarget);
		if (v != virtualTarget) //redundant to keep virtual target graph vertex's dummy map in mapsOnTheOptimalPath
		{
			cout << v << "\t";
			mapsOnTheOptimalPath.push_back(v);
		}
	}
}

void Graph::computePaths(int target)
{
	//computes a set of N! paths from meshes[0] (first mesh of the collection) to mesh meshes[target]; when i say N!, i actually mean just topK 'cos i use topK mappings at all times
	
	if (target > 1) //general case which appends N! pairwise maps b/w target-1 and target to N! existing cool paths stored via dynamic programming at pathsTo[target-1] list
		for (int p = 0; p < (int) pathsTo[target-1].size(); p++) //for each of N! existing path b/w meshes[0] to meshes[target-1]
		{
			//path is a vector/collection of maps where maps[i] represents gverts[ maps[i] ] as a path vertex/map
			GPath* path = pathsTo[target-1][p]; //(p+1)'th path from meshes[0] to target-1, so of size target-1 for sure: when target = 4, 1 map b/w 0-1, 1 b/w 1-2, 1 b/w 2-3 hence target-1 = 3 maps b/w 0-3
if ((int) path->maps.size() != target-1) {cout << "# maps in path wrong\n"; exit(0);}
			//among N! pairwise maps to be tried below b/w target-1 and target, i'll select the one that fits best to path
			int bestFitIdx; //keeping idx is simpler than keeping the pointer GVertex* bestFit
			//best fitting tests are based on distortion summation
			float distortionSumMin = INF; //best path to be added into pathesTo[target] minimizes the sum of distortions of the composed maps on it
			//try N! pairwise maps b/w target-1 and target
			for (int m = 0; m < topK; m++)
			{
				//map b/w target-1 and target to be appended to the path
				vector< int > tryMap = gverts[(target-1)*topK + m]->map;
				//sum of isometric distortions imposed by the path of maps: gverts[ path->maps[0] ]->map, gverts[ path->maps[1] ]->map, .., tryMap								
				float distortionSum = 0.0f; //based on the current tryMap
				for (int src = 0; src < target-1; src++) //isometric distortions to be added will be b/w meshes 0-target, 1-target, 2-target, .., (target-1)-target, where for all cases except the last i need to compose maps
				{
					vector< int > composedMap;
					//compose the map b/w src and target
					for (int match = 0; match < N; match++) //for each of N matches in gverts[ path->maps[src] ]->map
					{
						int i = gverts[ path->maps[src] ]->map[2*match], ci = gverts[ path->maps[src] ]->map[2*match + 1]; //match i-ci where i from meshes[src], ci from meshes[src+1]
						for (int a = src+1; a < target; a++) //propagate ci info to the map at right
						{
							vector< int > map = (a < target-1 ? gverts[ path->maps[a] ]->map : tryMap); //in the last iteration, path->maps[a] undefined and therefore needs to be replaced w/ tryMap
							//look for the match ci-cii in the map at right, namely gverts[ path->maps[a] ]->map
							int cii = -1;
							for (int match2 = 0; match2 < N; match2++) //for each of N matches in gverts[ path->maps[a] ]->map
								if (map[2*match2] == ci)
								{
									cii = map[2*match2 + 1];
									break;
								}
if (cii == -1) {cout << "cii cannot remain -1 during map composition\t" << target << endl; exit(0);} //INF constant may be too small to pass the first < test below whose result may make cii = -1; just increase INF or do some /= target normalization
							ci = cii; //propagate new ci info to the map at right
						}
						//i-ci as the (match+1)'th match of composedMap b/w meshes[src]-meshes[target]
						composedMap.push_back(i);
						composedMap.push_back(ci);
					} //end of match
					//isometric distortion of the composed map b/w meshes[src]-meshes[target]
					distortionSum += meshes[src]->getCorrespListDistortion(meshes[target], composedMap);
//if (p==0 && ((target == 2 && src == 0) || (target == 3 && src == 0))) cout << meshes[target]->verts[0]->coords[0] << " & " << meshes[src]->verts[0]->coords[0] << "\t" << meshes[src]->getCorrespListDistortion(meshes[target], composedMap) << " for " << target << src << endl;//order4
//if (p==0 && ((target == 3 && src == 0) || (target == 3 && src == 2))) cout << meshes[target]->verts[0]->coords[0] << " & " << meshes[src]->verts[0]->coords[0] << "\t" << meshes[src]->getCorrespListDistortion(meshes[target], composedMap) << " for " << target << src << endl;//order5
				} //end of src
				//isometric distortion of the map b/w meshes[target-1]-meshes[target] which is not composed 'cos it's simply ready
				distortionSum += meshes[target-1]->getCorrespListDistortion(meshes[target], tryMap);
				//distortionSum /= target; //target += done above so scale distortionSum in [0,1] which is required for large sequences as my INF initialization above may fail after accumulations; not needed 'cos i currently use m == 0 trick in that test


//finally there is the overall distortion of this current path (to which tryMap is appended)
distortionSum += path->overallDistortion;


				//current distortionSum based on the current tryMap
				if (m == 0 || distortionSum < distortionSumMin)
				{
					distortionSumMin = distortionSum;
					bestFitIdx = (target-1)*topK + m; //idx used in the retrieval of this best fit map tryMap
				}
			} //end of m
			//(p+1)'th path of N! pathes from 0 to target: path + bestFitIdx
			GPath* bestPath = new GPath();
			for (int j = 0; j < (int) path->maps.size(); j++)
				bestPath->maps.push_back(path->maps[j]); //copy the generator path
			bestPath->maps.push_back(bestFitIdx); //append the best fit map


//finally there is the overall distortion of this newly created bestPath to be used in the next iteration (on p) during distortionSum += above; the new overall distortion is all C(bestPath->maps.size+1, 2) pairwise distortions = 
//all pairwise distortions including the new mesh target (0-target, 1-target, 2-target, .., (target-1)-target) which is stored in distortionSumMin plus all pairwise distortions on the generating path which is stored in path->overallDistortion again via dynamic programming
bestPath->overallDistortion = //(distortionSumMin + path->overallDistortion);// / 2; //distortionSumMin in [0,1] plus path->overallDistortion in [0,1] so div by 2 to ensure bestPath->overallDistortion in [0,1] (redundant but anyway); no normalization anymore 'cos INF test is stable currently thanks to m==0 check
								distortionSumMin; //path->overallDistortion already added above


			pathsTo[target].push_back(bestPath); //N! additions'll be done here for sure 'cos p runs until (int) pathsTo[target-1].size() = N!
//cout << p << ": generatorPathDist " << path->overallDistortion << "\t+ newDist " << distortionSumMin << " = " << bestPath->overallDistortion << " [" << bestFitIdx << "]\n";
		}
	else //initial call b/w meshes[0] and [1] cannot use pathsTo[] so write here a modified version of the block just-above
		//path is a vector/collection of maps where maps[i] represents gverts[ maps[i] ] as a path vertex/map; for the initial 0-1 case here, each of N! paths will consist of 1-size maps, i.e. path'll have only 1 vertex/map on it
		for (int p = 0; p < topK; p++)
		{
			//(p+1)'th path of N! pathes from 0 to target=1: path (none existing) + m (bestFitIdx)
			GPath* bestPath = new GPath();
			//no generator path to copy
			bestPath->maps.push_back(p); //append the best fit map, no need to call it best fit actually 'cos there's no existing path that p is tried to be fit
			bestPath->overallDistortion = gverts[p]->distortion; //bestPath consists of 1 map and hence the overall distortion is the corresponding gverts[p].distortion
			pathsTo[target].push_back(bestPath); //N! additions'll be done here for sure 'cos p runs until (int) pathsTo[target-1].size() = N!
//cout << p << ": " << bestPath->overallDistortion << endl;
		}
}
void Graph::computePaths2(int target)
{
	//same as computePaths() except pathsTo[] filled slightly differently; here we try N! possible paths for each map to be appended

	///////////// slight difference explanation ////////////////////
	//w/ computePaths(), dynamic programming keeps a pathsTo[i][1..N!] table that stores N! optimal paths from meshes[0] to meshes[i] and appends the best of N! candidate maps to each of them; may be missing some maps w/ this tactic, so-called the missing maps problem
	//w/ computePaths2(), pathsTo[i][1..N!] table is created slightly differently (yet still have the missing maps problem); above, i was trying N! candidate maps from the next layer L as the last map/gvert on the current path pathsTo[i][1..N!]; this way some maps at L may be left unused; now, proceed as follows:
	//for each map m \in L
	//	try m w/ each path p \in pathsTo[i][1..N!]
	//		append m to the best-fit p
	//this way, all maps in L are assigned to a path for sure, none left unused; however, this time maps at L-1 may be left unmatched; so, we again cannot get the optimal solution that minimizes all pairwise distortions in the end; as an approx solution, use DYNAMIC_PROGRAMMING_ALGO_V3 which maintains topP best-fit p instead of just 1 here
	///////////// slight difference explanation ends ////////////////////
	
	if (target > 1) //general case which appends N! pairwise maps b/w target-1 and target to N! existing cool paths stored via dynamic programming at pathsTo[target-1] list
		for (int m = 0; m < topK; m++) //for each N! pairwise maps b/w target-1 and target from the next layer (next w.r.t. the current end point of the pathsTo[target-1]; so, this next'll be added to one of the pathsTo[target-1])
		{
			//map b/w target-1 and target to be appended to the best-fit path after tryings below
			vector< int > tryMap = gverts[(target-1)*topK + m]->map;
			//decide the best-fit path for the current tryMap
			int bestFitIdx; //keeping idx is simpler than keeping the pointer GPath* bestFit
			//best fitting tests are based on distortion summation
			float distortionSumMin = INF; //best path to be added into pathesTo[target] minimizes the sum of distortions of the composed maps on it
			for (int p = 0; p < (int) pathsTo[target-1].size(); p++) //for each of N! existing path b/w meshes[0] to meshes[target-1]
			{
				//path is a vector/collection of maps where maps[i] represents gverts[ maps[i] ] as a path vertex/map
				GPath* path = pathsTo[target-1][p]; //(p+1)'th path from meshes[0] to target-1, so of size target-1 for sure: when target = 4, 1 map b/w 0-1, 1 b/w 1-2, 1 b/w 2-3 hence target-1 = 3 maps b/w 0-3
//if ((int) path->maps.size() != target-1) {cout << "# maps in path wrong\n"; exit(0);}				
				//sum of isometric distortions imposed by the path of maps: gverts[ path->maps[0] ]->map, gverts[ path->maps[1] ]->map, .., tryMap								
				float distortionSum = 0.0f; //based on the current path
				for (int src = 0; src < target-1; src++) //isometric distortions to be added will be b/w meshes 0-target, 1-target, 2-target, .., (target-1)-target, where for all cases except the last i need to compose maps
				{
					vector< int > composedMap;
					//compose the map b/w src and target
					for (int match = 0; match < N; match++) //for each of N matches in gverts[ path->maps[src] ]->map
					{
						int i = gverts[ path->maps[src] ]->map[2*match], ci = gverts[ path->maps[src] ]->map[2*match + 1]; //match i-ci where i from meshes[src], ci from meshes[src+1]
						for (int a = src+1; a < target; a++) //propagate ci info to the map at right
						{
							vector< int > map = (a < target-1 ? gverts[ path->maps[a] ]->map : tryMap); //in the last iteration, path->maps[a] undefined and therefore needs to be replaced w/ tryMap
							//look for the match ci-cii in the map at right, namely gverts[ path->maps[a] ]->map
							int cii = -1;
							for (int match2 = 0; match2 < N; match2++) //for each of N matches in gverts[ path->maps[a] ]->map
								if (map[2*match2] == ci)
								{
									cii = map[2*match2 + 1];
									break;
								}
//if (cii == -1) {cout << "cii cannot remain -1 during map composition\t" << target << endl; exit(0);} //INF constant may be too small to pass the first < test below whose result may make cii = -1; just increase INF or do some /= target normalization
							ci = cii; //propagate new ci info to the map at right
						}
						//i-ci as the (match+1)'th match of composedMap b/w meshes[src]-meshes[target]
						composedMap.push_back(i);
						composedMap.push_back(ci);
					} //end of match
					//isometric distortion of the composed map b/w meshes[src]-meshes[target]
					distortionSum += meshes[src]->getCorrespListDistortion(meshes[target], composedMap);
				} //end of src
				//isometric distortion of the map b/w meshes[target-1]-meshes[target] which is not composed 'cos it's simply ready
				distortionSum += meshes[target-1]->getCorrespListDistortion(meshes[target], tryMap);
				//distortionSum /= target; //target += done above so scale distortionSum in [0,1] which is required for large sequences as my INF initialization above may fail after accumulations; not needed 'cos i currently use m == 0 trick in that test


//finally there is the overall distortion of this current path (to which tryMap is appended)
distortionSum += path->overallDistortion;


				//current distortionSum based on the current tryMap
				if (p == 0 || distortionSum < distortionSumMin)
				{
					distortionSumMin = distortionSum;
					bestFitIdx = p; //idx used in the retrieval of this best fit path path
				}
			} //end of p
			//pathsTo[target][bestFitIdx] = pathsTo[target-1][bestFitIdx] + tryMap as bestFitIdx and tryMap fit very good
			GPath* bestPath = new GPath();
			for (int j = 0; j < (int) pathsTo[target-1][bestFitIdx]->maps.size(); j++)
				bestPath->maps.push_back(pathsTo[target-1][bestFitIdx]->maps[j]); //copy the generator path
			bestPath->maps.push_back((target-1)*topK + m); //append the tryMap which is retrieved w/ (target-1)*topK + m above


//finally there is the overall distortion of this newly created bestPath to be used in the next iteration (on p) during distortionSum += above; the new overall distortion is all C(bestPath->maps.size+1, 2) pairwise distortions = 
//all pairwise distortions including the new mesh target (0-target, 1-target, 2-target, .., (target-1)-target) which is stored in distortionSumMin plus all pairwise distortions on the generating path which is stored in path->overallDistortion again via dynamic programming
bestPath->overallDistortion = //(distortionSumMin + pathsTo[target-1][bestFitIdx]->overallDistortion);// / 2; //distortionSumMin in [0,1] plus path->overallDistortion in [0,1] so div by 2 to ensure bestPath->overallDistortion in [0,1] (redundant but anyway); no normalization anymore 'cos INF test is stable currently thanks to m==0 check
							  distortionSumMin; //path->overallDistortion already added above


			pathsTo[target].push_back(bestPath); //N! additions'll be done here for sure 'cos p runs until (int) pathsTo[target-1].size() = N!
//if(target==3)cout << (target-1)*topK + m << " appendedTo: expandingPathDist " << pathsTo[target-1][bestFitIdx]->overallDistortion << " + newDist " << distortionSumMin << " = " << bestPath->overallDistortion << " [" << bestFitIdx << "]\n";
		}
	else //initial call b/w meshes[0] and [1] cannot use pathsTo[] so write here a modified version of the block just-above
		//path is a vector/collection of maps where maps[i] represents gverts[ maps[i] ] as a path vertex/map; for the initial 0-1 case here, each of N! paths will consist of 1-size maps, i.e. path'll have only 1 vertex/map on it
		for (int p = 0; p < topK; p++)
		{
			//(p+1)'th path of N! pathes from 0 to target=1: path (none existing) + m (bestFitIdx)
			GPath* bestPath = new GPath();
			//no generator path to copy
			bestPath->maps.push_back(p); //append the best fit map, no need to call it best fit actually 'cos there's no existing path that p is tried to be fit
			bestPath->overallDistortion = gverts[p]->distortion; //bestPath consists of 1 map and hence the overall distortion is the corresponding gverts[p].distortion
			pathsTo[target].push_back(bestPath); //N! additions'll be done here for sure 'cos p runs until (int) pathsTo[target-1].size() = N!
//cout << p << ": " << bestPath->overallDistortion << endl;
		}
}
/*void Graph::computePaths3(int target)
{
	//same as computePaths2() except pathsTo[] filled slightly differently; here we try N!^^ possible paths for each map to be appended and maintains P best-fit paths instead of just 1 best-fit map as in computePaths2(), hence alleviating the missing maps problem; P=1 makes computePaths3() and computePaths2() exactly the same
	//^^when i say N!, i actually mean just topK 'cos i use topK of all N! mappings at all times

	if (P < 1 || P > 100 || P > topK)
	{
		cout << P << " = # best-fit paths for a map to be appended is not in [1, 100] or greater than: " << topK << endl;
		exit(0);
	}
	int nPathsToAMesh = P*topK;
	if (target > 1) //general case which appends N! pairwise maps b/w target-1 and target to N! existing cool paths stored via dynamic programming at pathsTo[target-1] list
		for (int m = 0; m < topK; m++) //for each N! pairwise maps b/w target-1 and target from the next layer (next w.r.t. the current end point of the pathsTo[target-1]; so, this next'll be added to one of the pathsTo[target-1])
		{
			//map b/w target-1 and target to be appended to the best-fit path after tryings below
			vector< int > tryMap = gverts[(target-1)*topK + m]->map;			
if ((int) pathsTo[target-1].size() != nPathsToAMesh) {cout << "# best-fit paths entering meshes[" << target-1 << "] is wrong\n"; exit(0);}
			float* allDistortionSums = new float[nPathsToAMesh]; //top P of this many distortions (and their corresponding paths) will be chosen after sorting performed on this array
			int* correspondingPaths = new int[nPathsToAMesh]; //and this corresponding paths array
			for (int p = 0; p < (int) pathsTo[target-1].size(); p++) //for each of P*N! existing path b/w meshes[0] to meshes[target-1]
			{
				//path is a vector/collection of maps where maps[i] represents gverts[ maps[i] ] as a path vertex/map
				GPath* path = pathsTo[target-1][p]; //(p+1)'th path from meshes[0] to target-1, so of size target-1 for sure: when target = 4, 1 map b/w 0-1, 1 b/w 1-2, 1 b/w 2-3 hence target-1 = 3 maps b/w 0-3
if ((int) path->maps.size() != target-1) {cout << "# maps in path wrong\n"; exit(0);}
				//sum of isometric distortions imposed by the path of maps: gverts[ path->maps[0] ]->map, gverts[ path->maps[1] ]->map, .., tryMap								
				float distortionSum = 0.0f; //based on the current path
				for (int src = 0; src < target-1; src++) //isometric distortions to be added will be b/w meshes 0-target, 1-target, 2-target, .., (target-1)-target, where for all cases except the last i need to compose maps
				{
					vector< int > composedMap;
					//compose the map b/w src and target
					for (int match = 0; match < N; match++) //for each of N matches in gverts[ path->maps[src] ]->map
					{
						int i = gverts[ path->maps[src] ]->map[2*match], ci = gverts[ path->maps[src] ]->map[2*match + 1]; //match i-ci where i from meshes[src], ci from meshes[src+1]
						for (int a = src+1; a < target; a++) //propagate ci info to the map at right
						{
							vector< int > map = (a < target-1 ? gverts[ path->maps[a] ]->map : tryMap); //in the last iteration, path->maps[a] undefined and therefore needs to be replaced w/ tryMap
							//look for the match ci-cii in the map at right, namely gverts[ path->maps[a] ]->map
							int cii = -1;
							for (int match2 = 0; match2 < N; match2++) //for each of N matches in gverts[ path->maps[a] ]->map
								if (map[2*match2] == ci)
								{
									cii = map[2*match2 + 1];
									break;
								}
if (cii == -1) {cout << "cii cannot remain -1 during map composition\t" << target << endl; exit(0);} //INF constant may be too small to pass the first < test below whose result may make cii = -1; just increase INF or do some /= target normalization
							ci = cii; //propagate new ci info to the map at right
						}
						//i-ci as the (match+1)'th match of composedMap b/w meshes[src]-meshes[target]
						composedMap.push_back(i);
						composedMap.push_back(ci);
					} //end of match
					//isometric distortion of the composed map b/w meshes[src]-meshes[target]
					distortionSum += meshes[src]->getCorrespListDistortion(meshes[target], composedMap);
				} //end of src
				//isometric distortion of the map b/w meshes[target-1]-meshes[target] which is not composed 'cos it's simply ready
				distortionSum += meshes[target-1]->getCorrespListDistortion(meshes[target], tryMap);
				//distortionSum /= target; //target += done above so scale distortionSum in [0,1] which is required for large sequences as my INF initialization above may fail after accumulations; not needed 'cos i currently use m == 0 trick in that test


//finally there is the overall distortion of this current path (to which tryMap is appended)
distortionSum += path->overallDistortion;


				//current distortionSum based on the current tryMap
				allDistortionSums[p] = distortionSum;
				correspondingPaths[p] = p;
			} //end of p

			//decide the topP best-fit paths for the current tryMap; correspondingPaths[]'ll be overwritten w.r.t. the sorted allDistortionSums; no!! to save time i don't update correspondingPaths at all
			for (int i = 0; i < P; i++) //when P sorted paths made into correspondingPaths after overwrites below, no need to sort the remaining part [P..nPathsToAMesh-1] to save time; so, i < P instead of i < nPathsToAMesh
			{
				//get the min from the right side of [i] 'cos left of [i] already sorted and fixed
				float minDisto = INF;
				int minDistoIdx = -1;
				for (int j = i; j < nPathsToAMesh; j++)
					if (j == i || allDistortionSums[j] < minDisto) //j==i to make sure the first iteration passes this test 'cos sometimes INF may not be big enough
					{
						minDisto = allDistortionSums[j];
						minDistoIdx = j;
					}
				//replace [i] & [minDistoIdx]
				float tmp = allDistortionSums[i];
				allDistortionSums[i] = minDisto;
				allDistortionSums[minDistoIdx] = tmp;
				//replace corresponding paths
				int tmp2 = correspondingPaths[i];
				correspondingPaths[i] = correspondingPaths[minDistoIdx];
				correspondingPaths[minDistoIdx] = tmp2;
				//pathsTo[target][bestFitIdx] = pathsTo[target-1][bestFitIdx] + tryMap as bestFitIdx and tryMap fit very good
				GPath* bestPath = new GPath();
				for (int j = 0; j < (int) pathsTo[target-1][ correspondingPaths[i] ]->maps.size(); j++)
					bestPath->maps.push_back(pathsTo[target-1][ correspondingPaths[i] ]->maps[j]); //copy the generator path
				bestPath->maps.push_back((target-1)*topK + m); //append the tryMap which is retrieved w/ (target-1)*topK + m above


//finally there is the overall distortion of this newly created bestPath to be used in the next iteration (on p) during distortionSum += above; the new overall distortion is all C(bestPath->maps.size+1, 2) pairwise distortions = 
//all pairwise distortions including the new mesh target (0-target, 1-target, 2-target, .., (target-1)-target) which is stored in distortionSumMin plus all pairwise distortions on the generating path which is stored in path->overallDistortion again via dynamic programming
bestPath->overallDistortion = //(minDisto + pathsTo[target-1][ correspondingPaths[minDistoIdx] ]->overallDistortion);// / 2; //distortionSumMin in [0,1] plus path->overallDistortion in [0,1] so div by 2 to ensure bestPath->overallDistortion in [0,1] (redundant but anyway); no normalization anymore 'cos INF test is stable currently thanks to j==i check
							  minDisto; //path->overallDistortion already added above


				pathsTo[target].push_back(bestPath); //N! additions'll be done here for sure 'cos p runs until (int) pathsTo[target-1].size() = N!
//if (target==3)cout << "gverts" << (target-1)*topK + m << " appended to pathsTo[" << target-1 << "][" << i << "] = " << correspondingPaths[i] << " w/ cost\t" << minDisto << "\t" << (int) pathsTo[target].size()-1 << endl;
			} //end of i
//if (target==3){for(int ma=0;ma<pathsTo[target-1][10]->maps.size();ma++)cout<<pathsTo[target-1][10]->maps[ma] << " --> ";cout<<endl;}

			delete [] allDistortionSums;
			delete [] correspondingPaths; //AClass* obj = new AClass(); requires delete obj but arrays allocated with new[] must be deallocated with delete[]
		}
	else //initial call b/w meshes[0] and [1] cannot use pathsTo[] so write here a modified version of the block just-above
		//path is a vector/collection of maps where maps[i] represents gverts[ maps[i] ] as a path vertex/map; for the initial 0-1 case here, each of N! paths will consist of 1-size maps, i.e. path'll have only 1 vertex/map on it
		for (int p = 0; p < nPathsToAMesh; p++)
		{
			//(p+1)'th path of N! pathes from 0 to target=1: path (none existing) + m (bestFitIdx)
			GPath* bestPath = new GPath();
			//no generator path to copy
			if (p < topK)
			{
				bestPath->maps.push_back(p); //append the best fit map, no need to call it best fit actually 'cos there's no existing path that p is tried to be fit
				bestPath->overallDistortion = gverts[p]->distortion; //bestPath consists of 1 map and hence the overall distortion is the corresponding gverts[p].distortion
			}
			else
			{
				bestPath->maps.push_back(topK-1); //append dummy map (the last one of the first layer) to guarantee pathsTo[target-1].size = nPathsToAMesh for the upcoming iteration
				bestPath->overallDistortion = 0.99f; //dummy distortion to prevent the selection of this dummy map forever (not that dummy to prevent crashes in map compositions but due to .99 it'll never be in the final output path)
			}
			pathsTo[target].push_back(bestPath);
		}
}*/

void Graph::findTheOptimalPath(bool evalAllPairwiseDistortions, bool interpCoarseToDense)
{
	//evaluate the overall distortion implied by each of N! paths from meshes[0] to meshes[nMeshes-1] select the one yielding the min summation of all implied/composed pairwise maps as the optimal; if evalAllPairwiseDistortions not wanted, then simply return the path pathsTo[nMeshes-1] w/ min iverallDistortion value

	cout << "optimal path sought in a set of paths from the first mesh to the last..\n";

	//heavily adapted from computePaths() but summation uses more terms here; all C(nMeshes, 2) pairs are composed here but in computePaths() i composed just nMeshes-1 of them (0-target, 1-target, .., (target-1)-target) 'cos the others wouldn't have affected the distortionSumMin test w/ their fixed contributions
	//best/optimal path minimizes the sum of distortions of all of the composed pairwise maps on it
	float distortionSumMin = INF, overallGrdTruthDistortion = 0.0f, tmp; //ground-truth stuff and symmetric flip counters for evaluation only
	int bestPathIdx, nAdds, nsf;
	overallSymFlips = 0;
	for (int p = 0; p < (int) pathsTo[nMeshes-1].size(); p++) //for each of N! existing path b/w the first meshes[0] to the last meshes[nMeshes-1]
	{
		//path is a vector/collection of maps where maps[i] represents gverts[ maps[i] ] as a path vertex/map
		GPath* path = pathsTo[nMeshes-1][p]; //(p+1)'th path from meshes[0] to nMeshes-1, so of size nMeshes-1 for sure: when nMeshes = 5, 1 map b/w 0-1, 1 b/w 1-2, 1 b/w 2-3, 1 b/w 3-4 hence nMeshes-1 = 4 maps b/w 0-4
if ((int) path->maps.size() != nMeshes-1) {cout << "# maps in path wrong\n"; exit(0);}
		if (evalAllPairwiseDistortions) //if isometric and/or ground-truth evaluation demanded, then i need explicit composition of each of C(nMeshes, 2) pairwise maps
		{
			//sum of isometric distortions imposed by the path of maps: gverts[ path->maps[0] ]->map, gverts[ path->maps[1] ]->map, .., gverts[ path->maps[nMeshes-2] ]->map
			float distortionSum = 0.0f, ogtd = 0.0f; //based on the current tryMap
			nAdds = nsf = 0;
			for (int target = nMeshes-1; target > 0; target--)
				/////////////// 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1] ////////////////////////
				for (int src = 0; src < target; src++) //isometric distortions to be added will be b/w meshes 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1]
				{
					vector< int > composedMap;
					//compose the map b/w src and target
					for (int match = 0; match < N; match++) //for each of N matches in gverts[ path->maps[src] ]->map
					{
						int i = gverts[ path->maps[src] ]->map[2*match], ci = gverts[ path->maps[src] ]->map[2*match + 1]; //match i-ci where i from meshes[src], ci from meshes[src+1]
						for (int a = src+1; a < target; a++) //propagate ci info to the map at right
						{
							//path->maps[a] defined for all a values here unlike computePaths2() which requires a replacement w/ tryMap at the last iteration
							//look for the match ci-cii in the map at right, namely gverts[ path->maps[a] ]->map
							int cii = -1;
							for (int match2 = 0; match2 < N; match2++) //for each of N matches in gverts[ path->maps[a] ]->map
								if (gverts[ path->maps[a] ]->map[2*match2] == ci)
								{
									cii = gverts[ path->maps[a] ]->map[2*match2 + 1];
									break;
								}
//if (cii == -1) {cout << "cii cannot remain -1 during map composition\n"; exit(0);}
							ci = cii; //propagate new ci info to the map at right
						}
						//i-ci as the (match+1)'th match of composedMap b/w meshes[src]-meshes[target]
						composedMap.push_back(i);
						composedMap.push_back(ci);
					} //end of match
					//isometric distortion of the composed map b/w meshes[src]-meshes[target]
					distortionSum += meshes[src]->getCorrespListDistortion(meshes[target], composedMap);
					nAdds++;
					tmp = (grdTruth ? getGrdTruthDistortion(src, target, composedMap) : 0.0f);
					ogtd += tmp;
					nsf += (symFlipped(tmp) ? 1 : 0);
				} //end of src
			/////////////// 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1] ends ////////////////////////
			//current distortionSum based on the current tryMap
			if (p == 0 || distortionSum < distortionSumMin)
			{
				distortionSumMin = distortionSum;
				bestPathIdx = p; //idx used in the retrieval of this path
				overallGrdTruthDistortion = ogtd;
				overallSymFlips = nsf;
			}
		}
		else //just the optimal path is demanded so no composition stuff
		{			
			if (p == 0 || path->overallDistortion < distortionSumMin) //distortionSumMin naming not that good but anyway
			{
				distortionSumMin = path->overallDistortion;
				bestPathIdx = p; //idx used in the retrieval of this path
			}
		}
	}
	
	//evaluate the pairwise maps b/w consecutive meshes which does not necessarily be the shortest from 0 to nMeshes-1 but necessarily be the approximately optimal one, i.e. the one that minimizes the sum of all composed pairwise maps, not just b/w the consecutives
	cout << "optimal path (" << bestPathIdx << ") from the first mesh to the last:\n";
	float totalDistortion = 0.0f, optPathGrdTruthDistortion = 0.0f;
	int optPathSymFlips = 0;
	for (int i = 0; i < (int) pathsTo[nMeshes-1][bestPathIdx]->maps.size(); i++)
	{
		int optPathMap = pathsTo[nMeshes-1][bestPathIdx]->maps[i];
		mapsOnTheOptimalPath.push_back(optPathMap);
		//info
		cout << optPathMap << "\t";
		totalDistortion += gverts[optPathMap]->distortion;
		gverts[optPathMap]->grdTruthDistortion = (grdTruth ? getGrdTruthDistortion(gverts[optPathMap]->sidx, gverts[optPathMap]->tidx, gverts[optPathMap]->map) : 0.0f);
		optPathGrdTruthDistortion += gverts[optPathMap]->grdTruthDistortion;
		optPathSymFlips += (symFlipped(gverts[optPathMap]->grdTruthDistortion) ? 1 : 0);		
	}
	//this total or average is not necessarily smaller than shortest path solution in solve() 'cos here our claim is to get the optimal path which minimizes the total or average of all C(nMeshes,2) pairwise distortions, i.e. the ones on the path plus the implied/composed ones
	cout << "\ntotal & average isometric distortion on this optimal path: " << totalDistortion << " & " << totalDistortion / (int) mapsOnTheOptimalPath.size() << endl;
	if (evalAllPairwiseDistortions)
		//distortionSumMin is actually minimized and can be compared w/ solve().distortionSumMin and we expect distortionSumMin here to be better/smaller
		cout << "\ntotal & average isometric distortion over ALL (optimal path + compositions) pairwise distortions: " << distortionSumMin << " & " << distortionSumMin / nAdds << " (for C(" << nMeshes << ", 2) = " << nAdds << " pairwise maps)\n\n";
	//valuable grd-truth results for evaluation; similar to just-above, only overallGrdTruthDistortion's are comparable b/w here and solve()
	cout << "\ntotal & average ground-truth distortion on this optimal path: " << optPathGrdTruthDistortion << " & " << optPathGrdTruthDistortion / (int) mapsOnTheOptimalPath.size() << endl;
	if (evalAllPairwiseDistortions)
		 cout << "\ntotal & average ground-truth distortion over ALL (optimal path + compositions) pairwise distortions: " << overallGrdTruthDistortion << " & " << overallGrdTruthDistortion / nAdds << " (for C(" << nMeshes << ", 2) = " << nAdds << " pairwise maps)\n\n";
	//fprint the result
	for (int i = 0; i < (int) mapsOnTheOptimalPath.size(); i++)
		printMap(i, mapsOnTheOptimalPath[i]);
	cout << "done!\nmaps on the optimal path written to \\output folder\n\n";
	if (evalAllPairwiseDistortions)
		 cout << ">>>>>> SUMMARY (optimal path isometric distortion) & (over ALL isometric distortions): " << totalDistortion / (int) mapsOnTheOptimalPath.size() << " & " << distortionSumMin / nAdds << "\n"
			  << ">>>>>> SUMMARY (optimal path ground-truth distortion) & (over ALL ground-truth distortions): " << optPathGrdTruthDistortion / (int) mapsOnTheOptimalPath.size() << " & " << overallGrdTruthDistortion / nAdds << "\n"
			  << ">>>>>> SUMMARY (optimal path # symmetric flips) & (over ALL # symmetric flips): " << optPathSymFlips << " / " << (int) mapsOnTheOptimalPath.size() << " & " << overallSymFlips << " / " << nAdds << "\n\n";
	else
		 cout << ">>>>>> SUMMARY (optimal path isometric distortion) & (optimal path ground-truth distortion): " << totalDistortion / (int) mapsOnTheOptimalPath.size() << " & " << overallGrdTruthDistortion / nAdds << "\n\n";
	//to be used for triple runs only, i.e. if MESH_TRIPLETS_ONLY is defined
	overallIsoDistOfTriplet = distortionSumMin / nAdds;
	overallGrdTruthDistOfTriplet = overallGrdTruthDistortion / nAdds;

	if (interpCoarseToDense) //explicitly compose all maps implied by optPath; heavily adapted from evalAllPairwiseDistortions part above
	{
		cout << "interpolation of good coarse matches into 80x80 dense map..\n";
		GPath* optPath = pathsTo[nMeshes-1][bestPathIdx];
		float distortionSum = 0.0f, ogtd = 0.0f;
		nAdds = nsf = 0;
		for (int target = nMeshes-1; target > 0; target--)
			/////////////// 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1] ////////////////////////
			for (int src = 0; src < target; src++) //isometric distortions to be added will be b/w meshes 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1]
			{
				vector< int > composedMap, goodMatches;
				//compose the map b/w src and target
				for (int match = 0; match < N; match++) //for each of N matches in gverts[ path->maps[src] ]->map
				{
					int i = gverts[ optPath->maps[src] ]->map[2*match], ci = gverts[ optPath->maps[src] ]->map[2*match + 1]; //match i-ci where i from meshes[src], ci from meshes[src+1]
					for (int a = src+1; a < target; a++) //propagate ci info to the map at right
					{
						//path->maps[a] defined for all a values here unlike computePaths2() which requires a replacement w/ tryMap at the last iteration
						//look for the match ci-cii in the map at right, namely gverts[ path->maps[a] ]->map
						int cii = -1;
						for (int match2 = 0; match2 < N; match2++) //for each of N matches in gverts[ path->maps[a] ]->map
							if (gverts[ optPath->maps[a] ]->map[2*match2] == ci)
							{
								cii = gverts[ optPath->maps[a] ]->map[2*match2 + 1];
								break;
							}
//if (cii == -1) {cout << "cii cannot remain -1 during map composition\n"; exit(0);}
						ci = cii; //propagate new ci info to the map at right
					}
					//i-ci as the (match+1)'th match of composedMap b/w meshes[src]-meshes[target]
					composedMap.push_back(i);
					composedMap.push_back(ci);
				} //end of match
				
				//interpolate coarse composedMap into a dense map via perfect matching based on known goodMatches in composedMap; for goodMatches use the matches in composedMap before the first significant jump (of individual distortions plot)
				////////////////////////////// goodMatches[] creation //////////////////////////////////
				meshes[src]->getCorrespListDistortion(meshes[target], composedMap); //meshes[src].v.diso's set
				//fill goodMatches in sorted order from composedMap w.r.t. diso's
				for (int c = 0; c < N; c++)
					meshes[src]->verts[ composedMap[2*c] ]->processed = false; //true means added to goodMatches[]
				float jump = 0; //sum of first 2 distortion differences
				for (int c = 0; c < N; c++)
				{
					int nextMinCostMatch; //next match whose diso is current minimum (idx to composedMap)
					float minCost = INF;
					for (int cc = 0; cc < N; cc++)
						if (! meshes[src]->verts[ composedMap[2*cc] ]->processed && meshes[src]->verts[ composedMap[2*cc] ]->diso < minCost)
						{
							minCost = meshes[src]->verts[ composedMap[2*cc] ]->diso;
							nextMinCostMatch = 2*cc;
						}
					goodMatches.push_back(composedMap[nextMinCostMatch]);
					goodMatches.push_back(composedMap[nextMinCostMatch+1]);
					meshes[src]->verts[ composedMap[nextMinCostMatch] ]->processed = true;
					jump += (c == 1 || c == 2 ? meshes[src]->verts[ goodMatches[2*c] ]->diso - meshes[src]->verts[ goodMatches[2*c-2] ]->diso : 0.0f);
				}
				//get rid of the bad matches in goodMatches[] which are the ones before the first big jump
				composedMap.clear(); //naming sucks but goods of the sorted goodMatches[] are kept in this to be transferred to goodMatches[] later after the loop
				for (int c = 0; c < (int) goodMatches.size(); c += 2)
					if (c <= 2 || meshes[src]->verts[ goodMatches[c] ]->diso - meshes[src]->verts[ goodMatches[c-2] ]->diso < jump)
					{
						composedMap.push_back(goodMatches[c]);
						composedMap.push_back(goodMatches[c+1]);
					}
					else
						break;
				goodMatches = composedMap; //transfer
				////////////////////////////// goodMatches[] creation ends //////////////////////////////////

				////////////////////////////// dense matching //////////////////////////////////
				int denseN = 80; //# of dense evenly-spaced samples
				meshes[src]->computeEvenlySpacedSamples(denseN, -1);
				meshes[target]->computeEvenlySpacedSamples(denseN, -1);
				vector< int > denseMap = meshes[src]->denseMatchWith(meshes[target], denseN, goodMatches);				
				////////////////////////////// dense matching ends //////////////////////////////////

				//isometric and grd-truth distortion of the dense map b/w meshes[src]-meshes[target]
				float isoDisto = meshes[src]->getCorrespListDistortion(meshes[target], denseMap), grdDisto = getGrdTruthDistortion(src, target, denseMap);
				distortionSum += isoDisto;
				nAdds++;
				tmp = (grdTruth ? grdDisto : 0.0f);
				ogtd += tmp;
				nsf += (symFlipped(tmp) ? 1 : 0);
				if (src == target - 1) //consecutive shapes in collection, so this is on the path from the first mesh to the last; it is also optimal for sure 'cos i interpolate only the optimal coarse maps here
				{
					char fName[250], dName[250];
					sprintf(dName, "");//sprintf(dName, "SCAPE\\");
					sprintf(fName, "output\\DYNAMIC_PROGRAMMING_ALGO\\coarse\\interpolated dense\\%s%d-%d.dat", dName, src, target);
					FILE* fPtr = fopen(fName, "w");
					fprintf(fPtr, "Optimal map of size: %d\nMax geodesics on two surfaces: %f\t%f\n", denseN, meshes[src]->maxGeoDist, meshes[target]->maxGeoDist);
					for (int i = 0; i < (int) denseMap.size(); i += 2)
						fprintf(fPtr, "%d\t%d\n", denseMap[i], denseMap[i+1]);
					fprintf(fPtr, "%f\n", isoDisto); //isometric distortion
					fprintf(fPtr, "%f\n", grdDisto); //ground-truth distortion
					fclose(fPtr);
				}
			} //end of src
		/////////////// 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1] ends ////////////////////////
		if (nAdds > 0)
			cout << "\nresults on interpolated dense maps:\n>>>>>> SUMMARY for interp-dense (over ALL isometric distortions): " << distortionSum / nAdds << "\n"
				 << ">>>>>> SUMMARY for interp-dense (over ALL ground-truth distortions): " << ogtd / nAdds << "\n"
				 << ">>>>>> SUMMARY for interp-dense (over ALL # symmetric flips): " << nsf << " / " << nAdds << "\ndone!\n\n";
	}
}
void Graph::computeOptimalMultipleCorrespondence(int version, bool interpCoarseToDense)
{
	//w/ dynamic programming, finds the optimal path from the first mesh (meshes[0]) to the last (meshes[nMeshes-1]), which implies intermediate paths and hence maps as well, where optimality is defined as minimization of sum of all pairwise distortions, not just the distortions on the final path b/w the first and the last mesh (solve() does this)

	cout << "\n\n................DYNAMIC PROGRAMMING OVER THE GRAPH................\n";

	cout << "paths are created dynamically..\n";
	//initial pathes from 0 to 0 is null and the rest to be filled during computePaths(); pathsTo[i] keeps N! (a vector of) paths from meshes[0] to meshes[i]; N! is actually topK 'cos i use topK mappings at all times
	pathsTo = new vector< GPath* >[nMeshes];
	pathsTo[0].push_back(NULL); //pathes from 0 to 0 is null and will never be accessed; place holder just to make the upcoming idxs intuitive
	//compute paths over the graph that will be updated via dynamic programming as i increases
	for (int i = 1; i < nMeshes; i++)
		if (version == 3)
			cout << "version3 is commented out as it give the same result as v2\n";//computePaths3(i);
		else if (version == 1)
			computePaths(i);
		else
			computePaths2(i);
	cout << "done!\n\n";

	//now i've a set of N! paths from meshes[0] to meshes[nMeshes-1]; evaluate the overall distortion implied by each and select the shortest one
	findTheOptimalPath(true, interpCoarseToDense);
}

void Graph::bruteForceMultipleCorrespondence()
{
	//w/ exhaustive search, finds the optimal path from the first mesh (meshes[0]) to the last (meshes[nMeshes-1]), which implies intermediate paths and hence maps as well, where optimality is defined as minimization of sum of all pairwise distortions, not just the distortions on the final path b/w the first and the last mesh (solve() does this)

	cout << "\n\n................EXHAUSTIVE SEARCH OVER THE GRAPH................\n";
	float overallIsometricDistortion = INF, overallGrdTruthDistortion;
	int* bestPathOfMaps = new int[nMeshes-1];
	int nAdds, nsf;
	bool first = true;
	if (nMeshes == 3)
		for (int l0 = 0; l0 < topK; l0++) //for each N! (= topK in implementation) pairwise maps b/w meshes[0] and [1], i.e. layer0
		{
			int* pathOfMaps = new int[nMeshes-1]; //this path of maps keeps idxs to the maps/gverts so that i can compose the implied maps in the inner loop
			pathOfMaps[0] = l0; //vector< int > l0Map = gverts[l0]->map; effect
			for (int l1 = 0; l1 < topK; l1++) //for each N! (= topK in implementation) pairwise maps b/w meshes[1] and [2], i.e. layer1
			{
				pathOfMaps[1] = topK + l1; //vector< int > l1Map = gverts[topK + l1]->map; effect
				//sum of isometric distortions imposed by the path of maps: l0Map, l1Map = gverts[ pathOfMaps[0] ]->map, gverts[ pathOfMaps[1] ]->map, .., gverts[ pathOfMaps[nMeshes-2] ]->map
				float distortionSum = 0.0f, ogtd = 0.0f, tmp; //based on the current tryMap
				nAdds = nsf = 0;
				for (int target = nMeshes-1; target > 0; target--)
					/////////////// 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1] ////////////////////////
					for (int src = 0; src < target; src++) //isometric distortions to be added will be b/w meshes 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1]
					{
						vector< int > composedMap;
						//compose the map b/w src and target
						for (int match = 0; match < N; match++) //for each of N matches in gverts[ path->maps[src] ]->map
						{
							int i = gverts[ pathOfMaps[src] ]->map[2*match], ci = gverts[ pathOfMaps[src] ]->map[2*match + 1]; //match i-ci where i from meshes[src], ci from meshes[src+1]
							for (int a = src+1; a < target; a++) //propagate ci info to the map at right
							{
								//pathOfMaps[a] defined for all a values here unlike computePaths2() which requires a replacement w/ tryMap at the last iteration
								//look for the match ci-cii in the map at right, namely gverts[ path->maps[a] ]->map
								int cii = -1;
								for (int match2 = 0; match2 < N; match2++) //for each of N matches in gverts[ path->maps[a] ]->map
									if (gverts[ pathOfMaps[a] ]->map[2*match2] == ci)
									{
										cii = gverts[ pathOfMaps[a] ]->map[2*match2 + 1];
										break;
									}
if (cii == -1) {cout << "cii cannot remain -1 during map composition\n"; exit(0);}
								ci = cii; //propagate new ci info to the map at right
							}
							//i-ci as the (match+1)'th match of composedMap b/w meshes[src]-meshes[target]
							composedMap.push_back(i);
							composedMap.push_back(ci);
						} //end of match
						//isometric distortion of the composed map b/w meshes[src]-meshes[target]
						distortionSum += meshes[src]->getCorrespListDistortion(meshes[target], composedMap);
						nAdds++;
						tmp = (grdTruth ? getGrdTruthDistortion(src, target, composedMap) : 0.0f);
						ogtd += tmp;
						nsf += (symFlipped(tmp) ? 1 : 0);
					} //end of src
				/////////////// 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1] ends ////////////////////////
				//current distortionSum based on the current pathOfMaps
				if (first || distortionSum < overallIsometricDistortion) //make sure the very first call passes here
				{
					overallIsometricDistortion = distortionSum;						
					for (int i = 0; i < nMeshes-1; i++) //bestPathOfMaps = pathOfMaps; sucks due to delete [] below; so make explicit copy w/ a loop
						bestPathOfMaps[i] = pathOfMaps[i];
					overallGrdTruthDistortion = ogtd;
					overallSymFlips = nsf;
					first = false;
				}
			} //end of l1
			delete [] pathOfMaps;
		}
	else if (nMeshes == 4)
		for (int l0 = 0; l0 < topK; l0++) //for each N! (= topK in implementation) pairwise maps b/w meshes[0] and [1], i.e. layer0
		{
			int* pathOfMaps = new int[nMeshes-1]; //this path of maps keeps idxs to the maps/gverts so that i can compose the implied maps in the inner loop
			pathOfMaps[0] = l0; //vector< int > l0Map = gverts[l0]->map; effect
			for (int l1 = 0; l1 < topK; l1++) //for each N! (= topK in implementation) pairwise maps b/w meshes[1] and [2], i.e. layer1
			{
				pathOfMaps[1] = topK + l1; //vector< int > l1Map = gverts[topK + l1]->map; effect
				for (int l2 = 0; l2 < topK; l2++) //for each N! (= topK in implementation) pairwise maps b/w meshes[2] and [3], i.e. layer2
				{
					pathOfMaps[2] = 2*topK + l2; //vector< int > l2Map = gverts[2*topK + l2]->map; effect
					//sum of isometric distortions imposed by the path of maps: l0Map, l1Map, l2Map	= gverts[ pathOfMaps[0] ]->map, gverts[ pathOfMaps[1] ]->map, .., gverts[ pathOfMaps[nMeshes-2] ]->map
					float distortionSum = 0.0f, ogtd = 0.0f, tmp; //based on the current tryMap
					nAdds = nsf = 0;
					for (int target = nMeshes-1; target > 0; target--)
						/////////////// 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1] ////////////////////////
						for (int src = 0; src < target; src++) //isometric distortions to be added will be b/w meshes 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1]
						{
							vector< int > composedMap;
							//compose the map b/w src and target
							for (int match = 0; match < N; match++) //for each of N matches in gverts[ path->maps[src] ]->map
							{
								int i = gverts[ pathOfMaps[src] ]->map[2*match], ci = gverts[ pathOfMaps[src] ]->map[2*match + 1]; //match i-ci where i from meshes[src], ci from meshes[src+1]
								for (int a = src+1; a < target; a++) //propagate ci info to the map at right
								{
									//pathOfMaps[a] defined for all a values here unlike computePaths2() which requires a replacement w/ tryMap at the last iteration
									//look for the match ci-cii in the map at right, namely gverts[ path->maps[a] ]->map
									int cii = -1;
									for (int match2 = 0; match2 < N; match2++) //for each of N matches in gverts[ path->maps[a] ]->map
										if (gverts[ pathOfMaps[a] ]->map[2*match2] == ci)
										{
											cii = gverts[ pathOfMaps[a] ]->map[2*match2 + 1];
											break;
										}
if (cii == -1) {cout << "cii cannot remain -1 during map composition\n"; exit(0);}
									ci = cii; //propagate new ci info to the map at right
								}
								//i-ci as the (match+1)'th match of composedMap b/w meshes[src]-meshes[target]
								composedMap.push_back(i);
								composedMap.push_back(ci);
							} //end of match
							//isometric distortion of the composed map b/w meshes[src]-meshes[target]
							distortionSum += meshes[src]->getCorrespListDistortion(meshes[target], composedMap);
							nAdds++;
							tmp = (grdTruth ? getGrdTruthDistortion(src, target, composedMap) : 0.0f);
							ogtd += tmp;
							nsf += (symFlipped(tmp) ? 1 : 0);
						} //end of src
					/////////////// 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1] ends ////////////////////////
					//current distortionSum based on the current pathOfMaps
					if (first || distortionSum < overallIsometricDistortion) //make sure the very first call passes here
					{
						overallIsometricDistortion = distortionSum;						
						for (int i = 0; i < nMeshes-1; i++) //bestPathOfMaps = pathOfMaps; sucks due to delete [] below; so make explicit copy w/ a loop
							bestPathOfMaps[i] = pathOfMaps[i];
						overallGrdTruthDistortion = ogtd;
						overallSymFlips = nsf;
						first = false;
					}
				} //end of l2
			}
			delete [] pathOfMaps;
		}
	else if (nMeshes == 5)
		for (int l0 = 0; l0 < topK; l0++) //for each N! (= topK in implementation) pairwise maps b/w meshes[0] and [1], i.e. layer0
		{
			int* pathOfMaps = new int[nMeshes-1]; //this path of maps keeps idxs to the maps/gverts so that i can compose the implied maps in the inner loop
			pathOfMaps[0] = l0; //vector< int > l0Map = gverts[l0]->map; effect
			for (int l1 = 0; l1 < topK; l1++) //for each N! (= topK in implementation) pairwise maps b/w meshes[1] and [2], i.e. layer1
			{
				pathOfMaps[1] = topK + l1; //vector< int > l1Map = gverts[topK + l1]->map; effect
				for (int l2 = 0; l2 < topK; l2++) //for each N! (= topK in implementation) pairwise maps b/w meshes[2] and [3], i.e. layer2
				{
					pathOfMaps[2] = 2*topK + l2; //vector< int > l2Map = gverts[2*topK + l2]->map; effect
					for (int l3 = 0; l3 < topK; l3++) //for each N! (= topK in implementation) pairwise maps b/w meshes[3] and [4], i.e. layer3
					{
						pathOfMaps[3] = 3*topK + l3; //vector< int > l3Map = gverts[3*topK + l3]->map; effect
						//sum of isometric distortions imposed by the path of maps: l0Map, l1Map, l2Map, l3Map	= gverts[ pathOfMaps[0] ]->map, gverts[ pathOfMaps[1] ]->map, .., gverts[ pathOfMaps[nMeshes-2] ]->map
						float distortionSum = 0.0f, ogtd = 0.0f, tmp; //based on the current tryMap
						nAdds = nsf = 0;
						for (int target = nMeshes-1; target > 0; target--)
							/////////////// 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1] ////////////////////////
							for (int src = 0; src < target; src++) //isometric distortions to be added will be b/w meshes 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1]
							{
								vector< int > composedMap;
								//compose the map b/w src and target
								for (int match = 0; match < N; match++) //for each of N matches in gverts[ path->maps[src] ]->map
								{
									int i = gverts[ pathOfMaps[src] ]->map[2*match], ci = gverts[ pathOfMaps[src] ]->map[2*match + 1]; //match i-ci where i from meshes[src], ci from meshes[src+1]
									for (int a = src+1; a < target; a++) //propagate ci info to the map at right
									{
										//pathOfMaps[a] defined for all a values here unlike computePaths2() which requires a replacement w/ tryMap at the last iteration
										//look for the match ci-cii in the map at right, namely gverts[ path->maps[a] ]->map
										int cii = -1;
										for (int match2 = 0; match2 < N; match2++) //for each of N matches in gverts[ path->maps[a] ]->map
											if (gverts[ pathOfMaps[a] ]->map[2*match2] == ci)
											{
												cii = gverts[ pathOfMaps[a] ]->map[2*match2 + 1];
												break;
											}
if (cii == -1) {cout << "cii cannot remain -1 during map composition\n"; exit(0);}
										ci = cii; //propagate new ci info to the map at right
									}
									//i-ci as the (match+1)'th match of composedMap b/w meshes[src]-meshes[target]
									composedMap.push_back(i);
									composedMap.push_back(ci);
								} //end of match
								//isometric distortion of the composed map b/w meshes[src]-meshes[target]
								distortionSum += meshes[src]->getCorrespListDistortion(meshes[target], composedMap);
								nAdds++;
								tmp = (grdTruth ? getGrdTruthDistortion(src, target, composedMap) : 0.0f);
								ogtd += tmp;
								nsf += (symFlipped(tmp) ? 1 : 0);
							} //end of src
						/////////////// 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1] ends ////////////////////////
						//current distortionSum based on the current pathOfMaps
						if (first || distortionSum < overallIsometricDistortion) //make sure the very first call passes here
						{
							overallIsometricDistortion = distortionSum;						
							for (int i = 0; i < nMeshes-1; i++) //bestPathOfMaps = pathOfMaps; sucks due to delete [] below; so make explicit copy w/ a loop
								bestPathOfMaps[i] = pathOfMaps[i];
							overallGrdTruthDistortion = ogtd;
							overallSymFlips = nsf;
							first = false;
						}
					} //end of l3
				} //end of l2
			}
			delete [] pathOfMaps;
		}
	else if (nMeshes == 6)
		for (int l0 = 0; l0 < topK; l0++) //for each N! (= topK in implementation) pairwise maps b/w meshes[0] and [1], i.e. layer0
		{
			int* pathOfMaps = new int[nMeshes-1]; //this path of maps keeps idxs to the maps/gverts so that i can compose the implied maps in the inner loop
			pathOfMaps[0] = l0; //vector< int > l0Map = gverts[l0]->map; effect
			for (int l1 = 0; l1 < topK; l1++) //for each N! (= topK in implementation) pairwise maps b/w meshes[1] and [2], i.e. layer1
			{
				pathOfMaps[1] = topK + l1; //vector< int > l1Map = gverts[topK + l1]->map; effect
				for (int l2 = 0; l2 < topK; l2++) //for each N! (= topK in implementation) pairwise maps b/w meshes[2] and [3], i.e. layer2
				{
					pathOfMaps[2] = 2*topK + l2; //vector< int > l2Map = gverts[2*topK + l2]->map; effect
					for (int l3 = 0; l3 < topK; l3++) //for each N! (= topK in implementation) pairwise maps b/w meshes[3] and [4], i.e. layer3
					{
						pathOfMaps[3] = 3*topK + l3; //vector< int > l3Map = gverts[3*topK + l3]->map; effect
						for (int l4 = 0; l4 < topK; l4++) //for each N! (= topK in implementation) pairwise maps b/w meshes[4] and [5], i.e. layer4
						{
							pathOfMaps[4] = 4*topK + l4; //vector< int > l4Map = gverts[4*topK + l4]->map; effect
							//sum of isometric distortions imposed by the path of maps: l0Map, l1Map, l2Map, l3Map, l4Map = gverts[ pathOfMaps[0] ]->map, gverts[ pathOfMaps[1] ]->map, .., gverts[ pathOfMaps[nMeshes-2] ]->map
							float distortionSum = 0.0f, ogtd = 0.0f, tmp; //based on the current tryMap
							nAdds = nsf = 0;
							for (int target = nMeshes-1; target > 0; target--)
								/////////////// 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1] ////////////////////////
								for (int src = 0; src < target; src++) //isometric distortions to be added will be b/w meshes 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1]
								{
									vector< int > composedMap;
									//compose the map b/w src and target
									for (int match = 0; match < N; match++) //for each of N matches in gverts[ path->maps[src] ]->map
									{
										int i = gverts[ pathOfMaps[src] ]->map[2*match], ci = gverts[ pathOfMaps[src] ]->map[2*match + 1]; //match i-ci where i from meshes[src], ci from meshes[src+1]
										for (int a = src+1; a < target; a++) //propagate ci info to the map at right
										{
											//pathOfMaps[a] defined for all a values here unlike computePaths2() which requires a replacement w/ tryMap at the last iteration
											//look for the match ci-cii in the map at right, namely gverts[ path->maps[a] ]->map
											int cii = -1;
											for (int match2 = 0; match2 < N; match2++) //for each of N matches in gverts[ path->maps[a] ]->map
												if (gverts[ pathOfMaps[a] ]->map[2*match2] == ci)
												{
													cii = gverts[ pathOfMaps[a] ]->map[2*match2 + 1];
													break;
												}
if (cii == -1) {cout << "cii cannot remain -1 during map composition\n"; exit(0);}
											ci = cii; //propagate new ci info to the map at right
										}
										//i-ci as the (match+1)'th match of composedMap b/w meshes[src]-meshes[target]
										composedMap.push_back(i);
										composedMap.push_back(ci);
									} //end of match
									//isometric distortion of the composed map b/w meshes[src]-meshes[target]
									distortionSum += meshes[src]->getCorrespListDistortion(meshes[target], composedMap);
									nAdds++;
									tmp = (grdTruth ? getGrdTruthDistortion(src, target, composedMap) : 0.0f);
									ogtd += tmp;
									nsf += (symFlipped(tmp) ? 1 : 0);
								} //end of src
							/////////////// 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1] ends ////////////////////////
							//current distortionSum based on the current pathOfMaps
							if (first || distortionSum < overallIsometricDistortion) //make sure the very first call passes here
							{
								overallIsometricDistortion = distortionSum;						
								for (int i = 0; i < nMeshes-1; i++) //bestPathOfMaps = pathOfMaps; sucks due to delete [] below; so make explicit copy w/ a loop
									bestPathOfMaps[i] = pathOfMaps[i];
								overallGrdTruthDistortion = ogtd;
								overallSymFlips = nsf;
								first = false;
							}
						} //end of l4
					} //end of l3
				} //end of l2
			}
			delete [] pathOfMaps;
		}
	else if (nMeshes == 7)
		for (int l0 = 0; l0 < topK; l0++) //for each N! (= topK in implementation) pairwise maps b/w meshes[0] and [1], i.e. layer0
		{
			int* pathOfMaps = new int[nMeshes-1]; //this path of maps keeps idxs to the maps/gverts so that i can compose the implied maps in the inner loop
			pathOfMaps[0] = l0; //vector< int > l0Map = gverts[l0]->map; effect
			for (int l1 = 0; l1 < topK; l1++) //for each N! (= topK in implementation) pairwise maps b/w meshes[1] and [2], i.e. layer1
			{
				pathOfMaps[1] = topK + l1; //vector< int > l1Map = gverts[topK + l1]->map; effect
				for (int l2 = 0; l2 < topK; l2++) //for each N! (= topK in implementation) pairwise maps b/w meshes[2] and [3], i.e. layer2
				{
					pathOfMaps[2] = 2*topK + l2; //vector< int > l2Map = gverts[2*topK + l2]->map; effect
					for (int l3 = 0; l3 < topK; l3++) //for each N! (= topK in implementation) pairwise maps b/w meshes[3] and [4], i.e. layer3
					{
						pathOfMaps[3] = 3*topK + l3; //vector< int > l3Map = gverts[3*topK + l3]->map; effect
						for (int l4 = 0; l4 < topK; l4++) //for each N! (= topK in implementation) pairwise maps b/w meshes[4] and [5], i.e. layer4
						{
							pathOfMaps[4] = 4*topK + l4; //vector< int > l4Map = gverts[4*topK + l4]->map; effect
							for (int l5 = 0; l5 < topK; l5++) //for each N! (= topK in implementation) pairwise maps b/w meshes[5] and [6], i.e. layer5
							{
								pathOfMaps[5] = 5*topK + l5; //vector< int > l5Map = gverts[5*topK + l5]->map; effect
								//sum of isometric distortions imposed by the path of maps: l0Map, l1Map, l2Map, l3Map, l4Map, l5Map = gverts[ pathOfMaps[0] ]->map, gverts[ pathOfMaps[1] ]->map, .., gverts[ pathOfMaps[nMeshes-2] ]->map
								float distortionSum = 0.0f, ogtd = 0.0f, tmp; //based on the current tryMap
								nAdds = nsf = 0;
								for (int target = nMeshes-1; target > 0; target--)
									/////////////// 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1] ////////////////////////
									for (int src = 0; src < target; src++) //isometric distortions to be added will be b/w meshes 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1]
									{
										vector< int > composedMap;
										//compose the map b/w src and target
										for (int match = 0; match < N; match++) //for each of N matches in gverts[ path->maps[src] ]->map
										{
											int i = gverts[ pathOfMaps[src] ]->map[2*match], ci = gverts[ pathOfMaps[src] ]->map[2*match + 1]; //match i-ci where i from meshes[src], ci from meshes[src+1]
											for (int a = src+1; a < target; a++) //propagate ci info to the map at right
											{
												//pathOfMaps[a] defined for all a values here unlike computePaths2() which requires a replacement w/ tryMap at the last iteration
												//look for the match ci-cii in the map at right, namely gverts[ path->maps[a] ]->map
												int cii = -1;
												for (int match2 = 0; match2 < N; match2++) //for each of N matches in gverts[ path->maps[a] ]->map
													if (gverts[ pathOfMaps[a] ]->map[2*match2] == ci)
													{
														cii = gverts[ pathOfMaps[a] ]->map[2*match2 + 1];
														break;
													}
if (cii == -1) {cout << "cii cannot remain -1 during map composition\n"; exit(0);}
												ci = cii; //propagate new ci info to the map at right
											}
											//i-ci as the (match+1)'th match of composedMap b/w meshes[src]-meshes[target]
											composedMap.push_back(i);
											composedMap.push_back(ci);
										} //end of match
										//isometric distortion of the composed map b/w meshes[src]-meshes[target]
										distortionSum += meshes[src]->getCorrespListDistortion(meshes[target], composedMap);
										nAdds++;
										tmp = (grdTruth ? getGrdTruthDistortion(src, target, composedMap) : 0.0f);
										ogtd += tmp;
										nsf += (symFlipped(tmp) ? 1 : 0);
									} //end of src
								/////////////// 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1] ends ////////////////////////
								//current distortionSum based on the current pathOfMaps
								if (first || distortionSum < overallIsometricDistortion) //make sure the very first call passes here
								{
									overallIsometricDistortion = distortionSum;						
									for (int i = 0; i < nMeshes-1; i++) //bestPathOfMaps = pathOfMaps; sucks due to delete [] below; so make explicit copy w/ a loop
										bestPathOfMaps[i] = pathOfMaps[i];
									overallGrdTruthDistortion = ogtd;
									overallSymFlips = nsf;
									first = false;
								}
							} //end of l5
						} //end of l4
					} //end of l3
				} //end of l2
			}
			delete [] pathOfMaps;
		}
	else if (nMeshes == 8)
		for (int l0 = 0; l0 < topK; l0++) //for each N! (= topK in implementation) pairwise maps b/w meshes[0] and [1], i.e. layer0
		{
			int* pathOfMaps = new int[nMeshes-1]; //this path of maps keeps idxs to the maps/gverts so that i can compose the implied maps in the inner loop
			pathOfMaps[0] = l0; //vector< int > l0Map = gverts[l0]->map; effect
			for (int l1 = 0; l1 < topK; l1++) //for each N! (= topK in implementation) pairwise maps b/w meshes[1] and [2], i.e. layer1
			{
				pathOfMaps[1] = topK + l1; //vector< int > l1Map = gverts[topK + l1]->map; effect
				for (int l2 = 0; l2 < topK; l2++) //for each N! (= topK in implementation) pairwise maps b/w meshes[2] and [3], i.e. layer2
				{
					pathOfMaps[2] = 2*topK + l2; //vector< int > l2Map = gverts[2*topK + l2]->map; effect
					for (int l3 = 0; l3 < topK; l3++) //for each N! (= topK in implementation) pairwise maps b/w meshes[3] and [4], i.e. layer3
					{
						pathOfMaps[3] = 3*topK + l3; //vector< int > l3Map = gverts[3*topK + l3]->map; effect
						for (int l4 = 0; l4 < topK; l4++) //for each N! (= topK in implementation) pairwise maps b/w meshes[4] and [5], i.e. layer4
						{
							pathOfMaps[4] = 4*topK + l4; //vector< int > l4Map = gverts[4*topK + l4]->map; effect
							for (int l5 = 0; l5 < topK; l5++) //for each N! (= topK in implementation) pairwise maps b/w meshes[5] and [6], i.e. layer5
							{
								pathOfMaps[5] = 5*topK + l5; //vector< int > l5Map = gverts[5*topK + l5]->map; effect
								for (int l6 = 0; l6 < topK; l6++) //for each N! (= topK in implementation) pairwise maps b/w meshes[6] and [7], i.e. layer6
								{
									pathOfMaps[6] = 6*topK + l6; //vector< int > l6Map = gverts[6*topK + l6]->map; effect
									//sum of isometric distortions imposed by the path of maps: l0Map, l1Map, l2Map, l3Map, l4Map, l5Map, l6Map = gverts[ pathOfMaps[0] ]->map, gverts[ pathOfMaps[1] ]->map, .., gverts[ pathOfMaps[nMeshes-2] ]->map
									float distortionSum = 0.0f, ogtd = 0.0f, tmp; //based on the current tryMap
									nAdds = nsf = 0;
									for (int target = nMeshes-1; target > 0; target--)
										/////////////// 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1] ////////////////////////
										for (int src = 0; src < target; src++) //isometric distortions to be added will be b/w meshes 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1]
										{
											vector< int > composedMap;
											//compose the map b/w src and target
											for (int match = 0; match < N; match++) //for each of N matches in gverts[ path->maps[src] ]->map
											{
												int i = gverts[ pathOfMaps[src] ]->map[2*match], ci = gverts[ pathOfMaps[src] ]->map[2*match + 1]; //match i-ci where i from meshes[src], ci from meshes[src+1]
												for (int a = src+1; a < target; a++) //propagate ci info to the map at right
												{
													//pathOfMaps[a] defined for all a values here unlike computePaths2() which requires a replacement w/ tryMap at the last iteration
													//look for the match ci-cii in the map at right, namely gverts[ path->maps[a] ]->map
													int cii = -1;
													for (int match2 = 0; match2 < N; match2++) //for each of N matches in gverts[ path->maps[a] ]->map
														if (gverts[ pathOfMaps[a] ]->map[2*match2] == ci)
														{
															cii = gverts[ pathOfMaps[a] ]->map[2*match2 + 1];
															break;
														}
if (cii == -1) {cout << "cii cannot remain -1 during map composition\n"; exit(0);}
													ci = cii; //propagate new ci info to the map at right
												}
												//i-ci as the (match+1)'th match of composedMap b/w meshes[src]-meshes[target]
												composedMap.push_back(i);
												composedMap.push_back(ci);
											} //end of match
											//isometric distortion of the composed map b/w meshes[src]-meshes[target]
											distortionSum += meshes[src]->getCorrespListDistortion(meshes[target], composedMap);
											nAdds++;
											tmp = (grdTruth ? getGrdTruthDistortion(src, target, composedMap) : 0.0f);
											ogtd += tmp;
											nsf += (symFlipped(tmp) ? 1 : 0);
										} //end of src
									/////////////// 0-target, 1-target, 2-target, .., (target-1)-target for target in [nMeshes-1, 1] ends ////////////////////////
									//current distortionSum based on the current pathOfMaps
									if (first || distortionSum < overallIsometricDistortion) //make sure the very first call passes here
									{
										overallIsometricDistortion = distortionSum;						
										for (int i = 0; i < nMeshes-1; i++) //bestPathOfMaps = pathOfMaps; sucks due to delete [] below; so make explicit copy w/ a loop
											bestPathOfMaps[i] = pathOfMaps[i];
										overallGrdTruthDistortion = ogtd;
										overallSymFlips = nsf;
										first = false;
									}
								} //end of l6
							} //end of l5
						} //end of l4
					} //end of l3
				} //end of l2
			}
			delete [] pathOfMaps;
		}

	//evaluate the pairwise maps b/w consecutive meshes which does not necessarily be the shortest from 0 to nMeshes-1 but necessarily be the optimal one, i.e. the one that minimizes the sum of all composed pairwise maps, not just b/w the consecutives
	cout << "optimal path from the first mesh to the last:\n";
	float totalDistortion = 0.0f, optPathGrdTruthDistortion = 0.0f;
	int optPathSymFlips = 0;
	for (int i = 0; i < nMeshes-1; i++)
	{
		int optPathMap = bestPathOfMaps[i];
		//info
		cout << optPathMap << "\t";
		totalDistortion += gverts[optPathMap]->distortion;
		gverts[optPathMap]->grdTruthDistortion = (grdTruth ? getGrdTruthDistortion(gverts[optPathMap]->sidx, gverts[optPathMap]->tidx, gverts[optPathMap]->map) : 0.0f);
		optPathGrdTruthDistortion += gverts[optPathMap]->grdTruthDistortion;
		optPathSymFlips += (symFlipped(gverts[optPathMap]->grdTruthDistortion) ? 1 : 0);
		//fprint the result
		printMap(i, bestPathOfMaps[i]);
	}

	cout << "\n\n>>>>>> SUMMARY (optimal path isometric distortion) & (over ALL isometric distortions): " << totalDistortion / (nMeshes-1) << " & " << overallIsometricDistortion / nAdds << "\n"
		 << ">>>>>> SUMMARY (optimal path ground-truth distortion) & (over ALL ground-truth distortions): " << optPathGrdTruthDistortion / (nMeshes-1) << " & " << overallGrdTruthDistortion / nAdds << "\n"
		 << ">>>>>> SUMMARY (optimal path # symmetric flips) & (over ALL # symmetric flips): " << optPathSymFlips << " / " << (nMeshes-1) << " & " << overallSymFlips << " / " << nAdds << "\n\n";
}

void Graph::printMap(int fIdx, int mIdx, bool toCmdPrompt)
{
	//fprints and optionally prints map info

	fprintf(outFiles[fIdx], "Optimal map of size: %d\nMax geodesics on two surfaces: %f\t%f\n", N, meshes[ gverts[mIdx]->sidx ]->maxGeoDist, meshes[ gverts[mIdx]->tidx ]->maxGeoDist);
	for (int i = 0; i < (int) gverts[mIdx]->map.size(); i += 2)
	{
		fprintf(outFiles[fIdx], "%d\t%d\n", gverts[mIdx]->map[i], gverts[mIdx]->map[i+1]);
		if (toCmdPrompt)
			cout << gverts[mIdx]->map[i] << " - " << gverts[mIdx]->map[i+1] << endl;
	}
	fprintf(outFiles[fIdx], "%f\n", gverts[mIdx]->distortion); //= minDisto
	if (grdTruth)
		fprintf(outFiles[fIdx], "%f\n", gverts[mIdx]->grdTruthDistortion);
	if (toCmdPrompt)
	{
		if (grdTruth)
			cout << "isometric & ground-truth distortion: " << gverts[mIdx]->distortion << " & " << gverts[mIdx]->grdTruthDistortion << endl;
		else
			cout << "isometric distortion: " << gverts[mIdx]->distortion << endl;
	}
	fclose(outFiles[fIdx]);
}

inline float Graph::getGrdTruthDistortion(int m1, int m2, vector< int > tmpSafeCorrespTry)
{
	//returns ground-truth distortion assuming tha ground-truth matches are i-i

	float grdTruthCost = 0.0f;
	for (int s = 0; s < (int) tmpSafeCorrespTry.size() / 2; s++)
		if ((m1 + m2) % 2 == 0) //m1 & m2 both even or odd hence comparing meshes of the same class (even david, odd victoria) in which case grdTruthOnMesh1 is wrong and i use i-i ground-truths directly
			grdTruthCost += meshes[m1]->verts[ tmpSafeCorrespTry[2*s] ]->dSpanning[ tmpSafeCorrespTry[2*s + 1] ]; //this simple version cannot handle manual ground-truth settings b/w hybrid classes such as even vs. odd (or odd vs. even) in TOSCA_DAVIDVICTORIA case
		else
			grdTruthCost += meshes[m1]->verts[ tmpSafeCorrespTry[2*s] ]->dSpanning[ meshes[m2]->verts[ tmpSafeCorrespTry[2*s + 1] ]->grdTruthOnMesh1 ];
	grdTruthCost /= ((int) tmpSafeCorrespTry.size() / 2);

	return grdTruthCost;
}

inline bool Graph::symFlipped(float grdTruthDistortion)
{
//	if (((meshes[0]->theInput == SCAPE || meshes[0]->theInput == TOSCA_DAVID ||  meshes[0]->theInput == TOSCA_DAVIDVICTORIA || meshes[0]->theInput == TOSCA_VICTORIA || meshes[0]->theInput == TOSCA_GORILLA) && grdTruthDistortion > 0.18f) || //for human shapes .18 grd-truth distortion implies a sym flip, e.g. just hands flipped, or .5f for both hands and legs
//		((meshes[0]->theInput == TOSCA_CAT || meshes[0]->theInput == TOSCA_CENTAUR || meshes[0]->theInput == TOSCA_DOG || meshes[0]->theInput == TOSCA_WOLF) && grdTruthDistortion > 0.1f) || (meshes[0]->theInput == TOSCA_HORSE && grdTruthDistortion > 0.19f)) //see symmetricFlip deciders.txt for more
	if (grdTruthDistortion > 0.18f) //in this public release i can't decide this input-specific value automatically; so pick the most common choice which is assumed to fit all models from humans to 4-leggeds to centaurs, a safe assumption i might say; this is for evaluation purposes only and you can always interpret the resulting grdTruthDistortion yourself
		return true;
	return false;
}
