Context Navigation

source: cpp/frams/model/similarity/simil_model.cpp @ 969

Last change on this file since 969 was 877, checked in by Maciej Komosinski, 6 years ago
Introduced a shared function to avoid code duplication; more informative error message; frees memory in case of error
Property svn:eol-style set to `native`
File size: 68.7 KB

Rev	Line
[349]	1	// This file is a part of Framsticks SDK. http://www.framsticks.com/
[869]	2	// Copyright (C) 1999-2019 Maciej Komosinski and Szymon Ulatowski.
[349]	3	// See LICENSE.txt for details.
	4
[869]	5	// implementation of the ModelSimil class.
	6
[349]	7	#include "SVD/matrix_tools.h"
[869]	8	#include "hungarian/hungarian.h"
[349]	9	#include "simil_model.h"
	10	#include "simil_match.h"
	11	#include "frams/model/modelparts.h"
	12	#include "frams/util/list.h"
	13	#include "common/nonstd.h"
	14	#include <frams/vm/classes/genoobj.h>
[492]	15	#ifdef EMSCRIPTEN
[606]	16	#include <cstdlib>
[492]	17	#else
[606]	18	#include <stdlib.h>
[492]	19	#endif
[349]	20	#include <math.h>
	21	#include <string>
	22	#include <limits>
	23	#include <assert.h>
	24	#include <vector>
	25	#include <algorithm>
[666]	26	#include <cstdlib> //required for std::qsort in macos xcode
[349]	27
	28	#define DB(x) //define as x if you want to print debug information
	29
	30	const int ModelSimil::iNOFactors = 4;
	31	//depth of the fuzzy neighbourhood
	32	int fuzDepth = 0;
	33
	34	#define FIELDSTRUCT ModelSimil
	35
	36	static ParamEntry MSparam_tab[] = {
[869]	37	{ "Creature: Similarity", 1, 8, "ModelSimilarity", "Evaluates morphological dissimilarity. More information:\nhttp://www.framsticks.com/bib/Komosinski-et-al-2001\nhttp://www.framsticks.com/bib/Komosinski-and-Kubiak-2011\nhttp://www.framsticks.com/bib/Komosinski-2016\nhttps://doi.org/10.1007/978-3-030-16692-2_8", },
	38	{ "simil_method", 0, 0, "Similarity algorithm", "d 0 1 0 ~New (flexible criteria order, optimal matching)~Old (vertex degree order, greedy matching)", FIELD(matching_method), "",},
[606]	39	{ "simil_parts", 0, 0, "Weight of parts count", "f 0 100 0", FIELD(m_adFactors[0]), "Differing number of parts is also handled by the 'part degree' similarity component.", },
	40	{ "simil_partdeg", 0, 0, "Weight of parts' degree", "f 0 100 1", FIELD(m_adFactors[1]), "", },
	41	{ "simil_neuro", 0, 0, "Weight of neurons count", "f 0 100 0.1", FIELD(m_adFactors[2]), "", },
	42	{ "simil_partgeom", 0, 0, "Weight of parts' geometric distances", "f 0 100 0", FIELD(m_adFactors[3]), "", },
	43	{ "simil_fixedZaxis", 0, 0, "Fix 'z' (vertical) axis?", "d 0 1 0", FIELD(fixedZaxis), "", },
[818]	44	{ "simil_weightedMDS", 0, 0, "Should weighted MDS be used?", "d 0 1 0", FIELD(wMDS), "If activated, weighted MDS with vertex (i.e., Part) degrees as weights is used for 3D alignment of body structure.", },
[606]	45	{ "evaluateDistance", 0, PARAM_DONTSAVE \| PARAM_USERHIDDEN, "evaluate model dissimilarity", "p f(oGeno,oGeno)", PROCEDURE(p_evaldistance), "Calculates dissimilarity between two models created from Geno objects.", },
	46	{ 0, },
[349]	47	};
	48
	49	#undef FIELDSTRUCT
	50
	51	//////////////////////////////////////////////////////////////////////
	52	// Construction/Destruction
	53	//////////////////////////////////////////////////////////////////////
	54
	55	/** Constructor. Sets default weights. Initializes other fields with zeros.
	56	*/
[356]	57	ModelSimil::ModelSimil() : localpar(MSparam_tab, this), m_iDV(0), m_iDD(0), m_iDN(0), m_dDG(0.0)
[349]	58	{
[606]	59	localpar.setDefault();
[349]	60
[606]	61	m_Gen[0] = NULL;
	62	m_Gen[1] = NULL;
	63	m_Mod[0] = NULL;
	64	m_Mod[1] = NULL;
	65	m_aDegrees[0] = NULL;
	66	m_aDegrees[1] = NULL;
	67	m_aPositions[0] = NULL;
	68	m_aPositions[1] = NULL;
	69	m_fuzzyNeighb[0] = NULL;
	70	m_fuzzyNeighb[1] = NULL;
	71	m_Neighbours[0] = NULL;
	72	m_Neighbours[1] = NULL;
	73	m_pMatching = NULL;
[349]	74
[606]	75	//Determines whether "fuzzy vertex degree" should be used.
	76	//Currently "fuzzy vertex degree" is inactive.
[872]	77	isFuzzy = false;
[606]	78	fuzzyDepth = 10;
[869]	79
[817]	80	//Determines whether weighted MDS should be used.
	81	wMDS = 0;
[869]	82	//Determines whether best matching should be saved using hungarian similarity measure.
[872]	83	saveMatching = false;
[349]	84	}
	85
[869]	86	double ModelSimil::EvaluateDistanceGreedy(const Geno G0, const Geno G1)
[349]	87	{
[877]	88	double dResult = 0.0;
[349]	89
[606]	90	m_Gen[0] = G0;
	91	m_Gen[1] = G1;
[349]	92
[606]	93	// create models of objects to compare
[877]	94	m_Mod[0] = newModel(m_Gen[0]);
	95	m_Mod[1] = newModel(m_Gen[1]);
[349]	96
[877]	97	if (m_Mod[0] == NULL \|\| m_Mod[1] == NULL)
[872]	98	return 0.0;
[349]	99
[606]	100	// difference in the number of vertices (Parts) - positive
	101	// find object that has less parts (m_iSmaller)
	102	m_iDV = (m_Mod[0]->getPartCount() - m_Mod[1]->getPartCount());
	103	if (m_iDV > 0)
	104	m_iSmaller = 1;
	105	else
	106	{
	107	m_iSmaller = 0;
	108	m_iDV = -m_iDV;
	109	}
[349]	110
[606]	111	// check if index of the smaller organism is a valid index
	112	assert((m_iSmaller == 0) \|\| (m_iSmaller == 1));
	113	// validate difference in the parts number
	114	assert(m_iDV >= 0);
[349]	115
[606]	116	// create Parts matching object
	117	m_pMatching = new SimilMatching(m_Mod[0]->getPartCount(), m_Mod[1]->getPartCount());
	118	// validate matching object
	119	assert(m_pMatching != NULL);
	120	assert(m_pMatching->IsEmpty() == true);
[349]	121
	122
[606]	123	// assign matching function
	124	int (ModelSimil::* pfMatchingFunction) () = NULL;
[349]	125
[606]	126	pfMatchingFunction = &ModelSimil::MatchPartsGeometry;
[349]	127
[606]	128	// match Parts (vertices of creatures)
	129	if ((this->*pfMatchingFunction)() == 0)
	130	{
[872]	131	logPrintf("ModelSimil", "EvaluateDistanceGreedy", LOG_ERROR, "The matching function returned 0");
	132	return 0.0;
[606]	133	}
[349]	134
[606]	135	// after matching function call we must have full matching
	136	assert(m_pMatching->IsFull() == true);
[349]	137
[606]	138	DB(SaveIntermediateFiles();)
[349]	139
[606]	140	// count differences in matched parts
	141	if (CountPartsDistance() == 0)
	142	{
[872]	143	logPrintf("ModelSimil", "EvaluateDistanceGreedy", LOG_ERROR, "CountPartDistance()==0");
	144	return 0.0;
[606]	145	}
[349]	146
[606]	147	// delete degree arrays created in CreatePartInfoTables
	148	SAFEDELETEARRAY(m_aDegrees[0]);
	149	SAFEDELETEARRAY(m_aDegrees[1]);
[349]	150
[606]	151	// and position arrays
	152	SAFEDELETEARRAY(m_aPositions[0]);
	153	SAFEDELETEARRAY(m_aPositions[1]);
[349]	154
[606]	155	// in fuzzy mode delete arrays of neighbourhood and fuzzy neighbourhood
	156	if (isFuzzy)
	157	{
	158	for (int i = 0; i != 2; ++i)
	159	{
	160	for (int j = 0; j != m_Mod[i]->getPartCount(); ++j)
	161	{
	162	delete[] m_Neighbours[i][j];
	163	delete[] m_fuzzyNeighb[i][j];
	164	}
	165	delete[] m_Neighbours[i];
	166	delete[] m_fuzzyNeighb[i];
	167	}
[349]	168
[606]	169	}
[349]	170
[606]	171	// delete created models
	172	SAFEDELETE(m_Mod[0]);
	173	SAFEDELETE(m_Mod[1]);
	174
	175	// delete created matching
	176	SAFEDELETE(m_pMatching);
	177
	178	dResult = m_adFactors[0] * double(m_iDV) +
	179	m_adFactors[1] * double(m_iDD) +
	180	m_adFactors[2] * double(m_iDN) +
	181	m_adFactors[3] * double(m_dDG);
	182
	183	return dResult;
[349]	184	}
	185
	186	ModelSimil::~ModelSimil()
	187	{
[606]	188	// matching should have been deleted earlier
	189	assert(m_pMatching == NULL);
[349]	190	}
	191
	192	/** Creates files matching.txt, org0.txt and org1.txt containing information
	193	* about the matching and both organisms for visualization purpose.
	194	*/
	195	void ModelSimil::SaveIntermediateFiles()
	196	{
[606]	197	assert(m_pMatching->IsFull() == true);
	198	printf("Saving the matching to file 'matching.txt'\n");
	199	FILE *pMatchingFile = NULL;
	200	// try to open the file
	201	pMatchingFile = fopen("matching.txt", "wt");
	202	assert(pMatchingFile != NULL);
[349]	203
[606]	204	int iOrgPart; // original index of a Part
	205	int nBigger; // index of the larger organism
[349]	206
[606]	207	// check which object is bigger
	208	if (m_pMatching->GetObjectSize(0) >= m_pMatching->GetObjectSize(1))
	209	{
	210	nBigger = 0;
	211	}
	212	else
	213	{
	214	nBigger = 1;
	215	}
[349]	216
[606]	217	// print first line - original indices of Parts of the bigger organism
	218	fprintf(pMatchingFile, "[ ");
	219	for (iOrgPart = 0; iOrgPart < m_pMatching->GetObjectSize(nBigger); iOrgPart++)
	220	{
	221	fprintf(pMatchingFile, "%2i ", iOrgPart);
	222	}
	223	fprintf(pMatchingFile, "] : ORG[%i]\n", nBigger);
[349]	224
[606]	225	// print second line - matched original indices of the second organism
	226	fprintf(pMatchingFile, "[ ");
	227	for (iOrgPart = 0; iOrgPart < m_pMatching->GetObjectSize(nBigger); iOrgPart++)
	228	{
	229	int iSorted; // index of the iOrgPart after sorting (as used by matching)
	230	int iSortedMatched; // index of the matched Part (after sorting)
	231	int iOrginalMatched; // index of the matched Part (the original one)
[349]	232
[606]	233	// find the index of iOrgPart after sorting (in m_aDegrees)
	234	for (iSorted = 0; iSorted < m_Mod[nBigger]->getPartCount(); iSorted++)
	235	{
	236	// for each iSorted, an index in the sorted m_aDegrees array
	237	if (m_aDegrees[nBigger][iSorted][0] == iOrgPart)
	238	{
	239	// if the iSorted Part is the one with iOrgPart as the orginal index
	240	// remember the index
	241	break;
	242	}
	243	}
	244	// if the index iSorted was found, then this condition is met
	245	assert(iSorted < m_Mod[nBigger]->getPartCount());
[349]	246
[606]	247	// find the matched sorted index
	248	if (m_pMatching->IsMatched(nBigger, iSorted))
	249	{
	250	// if Part iOrgPart is matched
	251	// then get the matched Part (sorted) index
	252	iSortedMatched = m_pMatching->GetMatchedIndex(nBigger, iSorted);
	253	assert(iSortedMatched >= 0);
	254	// and find its original index
	255	iOrginalMatched = m_aDegrees[1 - nBigger][iSortedMatched][0];
	256	fprintf(pMatchingFile, "%2i ", iOrginalMatched);
	257	}
	258	else
	259	{
	260	// if the Part iOrgPart is not matched
	261	// just print "X"
	262	fprintf(pMatchingFile, " X ");
	263	}
	264	} // for ( iOrgPart )
[349]	265
[606]	266	// now all matched Part indices are printed out, end the line
	267	fprintf(pMatchingFile, "] : ORG[%i]\n", 1 - nBigger);
[349]	268
[606]	269	// close the file
	270	fclose(pMatchingFile);
	271	// END TEMP
[349]	272
[606]	273	// TEMP
	274	// print out the 2D positions of Parts of both of the organisms
	275	// to files "org0.txt" and "org1.txt" using the original indices of Parts
	276	int iModel; // index of a model (an organism)
	277	FILE *pModelFile;
	278	for (iModel = 0; iModel < 2; iModel++)
	279	{
	280	// for each iModel, a model of a compared organism
	281	// write its (only 2D) positions to a file "org<iModel>.txt"
	282	// construct the model filename "org<iModel>.txt"
	283	std::string sModelFilename("org");
	284	// char *szModelIndex = "0"; // the index of the model (iModel) in the character form
	285	char szModelIndex[2];
	286	sprintf(szModelIndex, "%i", iModel);
	287	sModelFilename += szModelIndex;
	288	sModelFilename += ".txt";
	289	// open the file for writing
	290	pModelFile = fopen(sModelFilename.c_str(), "wt"); //FOPEN_WRITE
	291	assert(pModelFile != NULL);
	292	// write the 2D positions of iModel to the file
	293	int iOriginalPart; // an original index of a Part
	294	for (iOriginalPart = 0; iOriginalPart < m_Mod[iModel]->getPartCount(); iOriginalPart++)
	295	{
	296	// for each iOriginalPart, a Part of the organism iModel
	297	// get the 2D coordinates of the Part
	298	double dPartX = m_aPositions[iModel][iOriginalPart].x;
	299	double dPartY = m_aPositions[iModel][iOriginalPart].y;
	300	// print the line: <iOriginalPart> <dPartX> <dPartY>
	301	fprintf(pModelFile, "%i %.4lf %.4lf\n", iOriginalPart, dPartX, dPartY);
	302	}
	303	// close the file
	304	fclose(pModelFile);
	305	}
[349]	306	}
	307
	308	/** Comparison function required for qsort() call. Used while sorting groups of
[606]	309	Parts with respect to degree. Compares two TDN structures
	310	with respect to their [1] field (degree). Highest degree goes first.
	311	@param pElem1 Pointer to the TDN structure of some Part.
	312	@param pElem2 Pointer to the TDN structure of some Part.
	313	@return (-1) - pElem1 should go first, 0 - equal, (1) - pElem2 should go first.
	314	*/
[349]	315	int ModelSimil::CompareDegrees(const void pElem1, const void pElem2)
	316	{
[606]	317	int tdn1 = (int )pElem1;
	318	int tdn2 = (int )pElem2;
[349]	319
[606]	320	if (tdn1[1] > tdn2[1])
	321	{
	322	// when degree - tdn1[1] - is BIGGER
	323	return -1;
	324	}
	325	else
	326	if (tdn1[1] < tdn2[1])
	327	{
[869]	328	// when degree - tdn2[1] - is BIGGER
	329	return 1;
[606]	330	}
	331	else
	332	{
	333	return 0;
	334	}
[349]	335	}
	336
[869]	337	/** Comparison function required for qsort() call. Used while sorting groups of
	338	Parts with respect to fuzzy vertex degree. Compares two TDN structures
	339	with respect to their [4] field ( fuzzy vertex degree). Highest degree goes first.
	340	@param pElem1 Pointer to the TDN structure of some Part.
	341	@param pElem2 Pointer to the TDN structure of some Part.
	342	@return (-1) - pElem1 should go first, 0 - equal, (1) - pElem2 should go first.
	343	*/
	344	int ModelSimil::CompareFuzzyDegrees(const void pElem1, const void pElem2)
	345	{
	346	int tdn1 = (int )pElem1;
	347	int tdn2 = (int )pElem2;
	348
	349	if (tdn1[4] > tdn2[4])
	350	{
	351	// when degree - tdn1[4] - is BIGGER
	352	return -1;
	353	}
	354	else
	355	if (tdn1[4] < tdn2[4])
	356	{
	357	// when degree - tdn2[4] - is BIGGER
	358	return 1;
	359	}
	360	else
	361	{
	362	return 0;
	363	}
	364	}
	365
[349]	366	/** Comparison function required for qsort() call. Used while sorting groups of Parts with
[606]	367	the same degree. Firstly, compare NIt. Secondly, compare Neu. If both are equal -
	368	compare Parts' original index (they are never equal). So this sorting assures
	369	that the order obtained is deterministic.
	370	@param pElem1 Pointer to the TDN structure of some Part.
	371	@param pElem2 Pointer to the TDN structure of some Part.
	372	@return (-1) - pElem1 should go first, 0 - equal, (1) - pElem2 should go first.
	373	*/
[349]	374	int ModelSimil::CompareConnsNo(const void pElem1, const void pElem2)
	375	{
[606]	376	// pointers to TDN arrays
	377	int tdn1, tdn2;
	378	// definitions of elements being compared
	379	tdn1 = (int *)pElem1;
	380	tdn2 = (int *)pElem2;
[349]	381
[606]	382	// comparison according to Neural Connections (to jest TDN[2])
[782]	383	if (tdn1[NEURO_CONNS] > tdn2[NEURO_CONNS])
[606]	384	{
	385	// when number of NConn Elem1 is BIGGER
	386	return -1;
	387	}
	388	else
[782]	389	if (tdn1[NEURO_CONNS] < tdn2[NEURO_CONNS])
[606]	390	{
[869]	391	// when number of NConn Elem1 is SMALLER
	392	return 1;
[606]	393	}
	394	else
	395	{
	396	// when numbers of NConn are EQUAL
	397	// compare Neu numbers (TDN[3])
	398	if (tdn1[NEURONS] > tdn2[NEURONS])
	399	{
	400	// when number of Neu is BIGGER for Elem1
	401	return -1;
	402	}
	403	else
	404	if (tdn1[NEURONS] < tdn2[NEURONS])
	405	{
[869]	406	// when number of Neu is SMALLER for Elem1
	407	return 1;
[606]	408	}
	409	else
	410	{
	411	// when numbers of Nconn and Neu are equal we check original indices
	412	// of Parts being compared
[349]	413
[606]	414	// comparison according to OrgIndex
	415	if (tdn1[ORIG_IND] > tdn2[ORIG_IND])
	416	{
	417	// when the number of NIt Deg1 id BIGGER
	418	return -1;
	419	}
	420	else
	421	if (tdn1[ORIG_IND] < tdn2[ORIG_IND])
	422	{
[869]	423	// when the number of NIt Deg1 id SMALLER
	424	return 1;
[606]	425	}
	426	else
	427	{
	428	// impossible, indices are alway different
	429	return 0;
	430	}
	431	}
	432	}
[349]	433	}
	434
[606]	435	/** Returns number of factors involved in final distance computation.
	436	These factors include differences in numbers of parts, degrees,
	437	number of neurons.
	438	*/
[349]	439	int ModelSimil::GetNOFactors()
	440	{
[606]	441	return ModelSimil::iNOFactors;
[349]	442	}
	443
	444	/** Prints the array of degrees for the given organism. Debug method.
	445	*/
	446	void ModelSimil::_PrintDegrees(int i)
	447	{
[606]	448	int j;
	449	printf("Organizm %i :", i);
	450	printf("\n ");
	451	for (j = 0; j < m_Mod[i]->getPartCount(); j++)
	452	printf("%3i ", j);
	453	printf("\nInd: ");
	454	for (j = 0; j < m_Mod[i]->getPartCount(); j++)
	455	printf("%3i ", (int)m_aDegrees[i][j][0]);
	456	printf("\nDeg: ");
	457	for (j = 0; j < m_Mod[i]->getPartCount(); j++)
	458	printf("%3i ", (int)m_aDegrees[i][j][1]);
	459	printf("\nNIt: ");
	460	for (j = 0; j < m_Mod[i]->getPartCount(); j++)
	461	printf("%3i ", (int)m_aDegrees[i][j][2]);
	462	printf("\nNeu: ");
	463	for (j = 0; j < m_Mod[i]->getPartCount(); j++)
	464	printf("%3i ", (int)m_aDegrees[i][j][3]);
	465	printf("\n");
[349]	466	}
	467
	468	/** Prints one array of ints. Debug method.
[606]	469	@param array Base pointer of the array.
	470	@param base First index of the array's element.
	471	@param size Number of elements to print.
	472	*/
[349]	473	void ModelSimil::_PrintArray(int *array, int base, int size)
	474	{
[606]	475	int i;
	476	for (i = base; i < base + size; i++)
	477	{
	478	printf("%i ", array[i]);
	479	}
	480	printf("\n");
[349]	481	}
	482
	483	void ModelSimil::_PrintArrayDouble(double *array, int base, int size)
	484	{
[606]	485	int i;
	486	for (i = base; i < base + size; i++)
	487	{
	488	printf("%f ", array[i]);
	489	}
	490	printf("\n");
[349]	491	}
	492
	493	/** Prints one array of parts neighbourhood.
[606]	494	@param index of organism
	495	*/
[349]	496	void ModelSimil::_PrintNeighbourhood(int o)
	497	{
[606]	498	assert(m_Neighbours[o] != 0);
	499	printf("Neighbourhhod of organism %i\n", o);
	500	int size = m_Mod[o]->getPartCount();
	501	for (int i = 0; i < size; i++)
	502	{
	503	for (int j = 0; j < size; j++)
	504	{
	505	printf("%i ", m_Neighbours[o][i][j]);
	506	}
	507	printf("\n");
	508	}
[349]	509	}
	510
[869]	511	/** Prints one array of parts fuzzy neighbourhood.
	512	@param index of organism
	513	*/
	514	void ModelSimil::_PrintFuzzyNeighbourhood(int o)
	515	{
	516	assert(m_fuzzyNeighb[o] != NULL);
	517	printf("Fuzzy neighbourhhod of organism %i\n", o);
	518	int size = m_Mod[o]->getPartCount();
	519	for (int i = 0; i < size; i++)
	520	{
	521	for (int j = 0; j < fuzzyDepth; j++)
	522	{
	523	printf("%f ", m_fuzzyNeighb[o][i][j]);
	524	}
	525	printf("\n");
	526	}
	527	}
	528
[877]	529	Model* ModelSimil::newModel(const Geno *g)
	530	{
	531	if (g == NULL)
	532	{
	533	logPrintf("ModelSimil", "newModel", LOG_ERROR, "NULL genotype pointer");
	534	return NULL;
	535	}
	536	Model m = new Model(g);
	537	if (!m->isValid())
	538	{
	539	logPrintf("ModelSimil", "newModel", LOG_ERROR, "Invalid model for the genotype of '%s'", g->getName().c_str());
	540	delete m;
	541	return NULL;
	542	}
	543	return m;
	544	}
	545
	546
[349]	547	/** Creates arrays holding information about organisms' Parts (m_aDegrees) andm_Neigh
[606]	548	fills them with initial data (original indices and zeros).
	549	Assumptions:
	550	- Models (m_Mod) are created and available.
	551	*/
[877]	552	int ModelSimil::ModelSimil::CreatePartInfoTables()
[349]	553	{
[606]	554	// check assumptions about models
	555	assert((m_Mod[0] != NULL) && (m_Mod[1] != NULL));
	556	assert(m_Mod[0]->isValid() && m_Mod[1]->isValid());
[349]	557
[606]	558	int i, j, partCount;
	559	// utwórz tablice na informacje o stopniach wierzchołków i liczbie neuroitems
	560	for (i = 0; i < 2; i++)
	561	{
	562	partCount = m_Mod[i]->getPartCount();
	563	// utworz i wypelnij tablice dla Parts wartosciami poczatkowymi
	564	m_aDegrees[i] = new TDN[partCount];
[349]	565
[606]	566	if (isFuzzy)
	567	{
	568	m_Neighbours[i] = new int*[partCount];
	569	m_fuzzyNeighb[i] = new float*[partCount];
	570	}
[349]	571
[872]	572	if (m_aDegrees[i] != NULL && ((!isFuzzy) \|\| (m_Neighbours[i] != NULL && m_fuzzyNeighb[i] != NULL)))
[606]	573	{
	574	// wypelnij tablice zgodnie z sensem TDN[0] - orginalny index
	575	// TDN[1], TDN[2], TDN[3] - zerami
	576	DB(printf("m_aDegrees[%i]: %p\n", i, m_aDegrees[i]);)
	577	for (j = 0; j < partCount; j++)
	578	{
[869]	579	m_aDegrees[i][j][0] = j;
	580	m_aDegrees[i][j][1] = 0;
	581	m_aDegrees[i][j][2] = 0;
	582	m_aDegrees[i][j][3] = 0;
	583	m_aDegrees[i][j][4] = 0;
[349]	584
[869]	585	// sprawdz, czy nie piszemy po jakims szalonym miejscu pamieci
	586	assert(m_aDegrees[i][j] != NULL);
[349]	587
[869]	588	if (isFuzzy)
[606]	589	{
[869]	590	m_Neighbours[i][j] = new int[partCount];
	591	for (int k = 0; k < partCount; k++)
	592	{
	593	m_Neighbours[i][j][k] = 0;
	594	}
[349]	595
[869]	596	m_fuzzyNeighb[i][j] = new float[fuzzyDepth];
	597	for (int k = 0; k < fuzzyDepth; k++)
	598	{
	599	m_fuzzyNeighb[i][j][k] = 0;
	600	}
	601
	602	assert(m_Neighbours[i][j] != NULL);
	603	assert(m_fuzzyNeighb[i][j] != NULL);
[606]	604	}
[349]	605
[606]	606	}
	607	}
	608	else
	609	{
[872]	610	logPrintf("ModelSimil", "CreatePartInfoTables", LOG_ERROR, "Not enough memory?");
	611	return 0;
[606]	612	}
	613	// utworz tablice dla pozycji 3D Parts (wielkosc tablicy: liczba Parts organizmu)
	614	m_aPositions[i] = new Pt3D[m_Mod[i]->getPartCount()];
	615	assert(m_aPositions[i] != NULL);
	616	// wypelnij tablice OnJoints i Anywhere wartościami początkowymi
	617	// OnJoint
	618	m_aOnJoint[i][0] = 0;
	619	m_aOnJoint[i][1] = 0;
	620	m_aOnJoint[i][2] = 0;
	621	m_aOnJoint[i][3] = 0;
	622	// Anywhere
	623	m_aAnywhere[i][0] = 0;
	624	m_aAnywhere[i][1] = 0;
	625	m_aAnywhere[i][2] = 0;
	626	m_aAnywhere[i][3] = 0;
	627	}
	628	return 1;
[349]	629	}
	630
	631	/** Computes degrees of Parts of both organisms. Fills in the m_aDegrees arrays
[606]	632	with proper information about degrees.
	633	Assumptions:
	634	- Models (m_Mod) are created and available.
	635	- Arrays m_aDegrees are created.
	636	*/
[349]	637	int ModelSimil::CountPartDegrees()
	638	{
[606]	639	// sprawdz zalozenie - o modelach
	640	assert((m_Mod[0] != NULL) && (m_Mod[1] != NULL));
	641	assert(m_Mod[0]->isValid() && m_Mod[1]->isValid());
[349]	642
[606]	643	// sprawdz zalozenie - o tablicach
	644	assert(m_aDegrees[0] != NULL);
	645	assert(m_aDegrees[1] != NULL);
[349]	646
[606]	647	Part P1, P2;
	648	int i, j, i1, i2;
[349]	649
[606]	650	// dla obu stworzen oblicz stopnie wierzcholkow
	651	for (i = 0; i < 2; i++)
	652	{
	653	// dla wszystkich jointow
	654	for (j = 0; j < m_Mod[i]->getJointCount(); j++)
	655	{
	656	// pobierz kolejny Joint
	657	Joint *J = m_Mod[i]->getJoint(j);
	658	// wez jego konce
	659	P1 = J->part1;
	660	P2 = J->part2;
	661	// znajdz ich indeksy w Modelu (indeksy orginalne)
	662	i1 = m_Mod[i]->findPart(P1);
	663	i2 = m_Mod[i]->findPart(P2);
	664	// zwieksz stopien odpowiednich Parts
	665	m_aDegrees[i][i1][DEGREE]++;
	666	m_aDegrees[i][i2][DEGREE]++;
	667	m_aDegrees[i][i1][FUZZ_DEG]++;
	668	m_aDegrees[i][i2][FUZZ_DEG]++;
	669	if (isFuzzy)
	670	{
	671	m_Neighbours[i][i1][i2] = 1;
	672	m_Neighbours[i][i2][i1] = 1;
	673	}
	674	}
	675	// dla elementow nie osadzonych na Parts (OnJoint, Anywhere) -
	676	// stopnie wierzchołka są już ustalone na zero
	677	}
[349]	678
[606]	679	if (isFuzzy)
	680	{
	681	CountFuzzyNeighb();
	682	}
[349]	683
[606]	684	return 1;
[349]	685	}
	686
	687	void ModelSimil::GetNeighbIndexes(int mod, int partInd, std::vector<int> &indexes)
	688	{
[606]	689	indexes.clear();
	690	int i, size = m_Mod[mod]->getPartCount();
[349]	691
[606]	692	for (i = 0; i < size; i++)
	693	{
	694	if (m_Neighbours[mod][partInd][i] > 0)
	695	{
	696	indexes.push_back(i);
	697	}
	698	}
[349]	699	}
	700
	701	int cmpFuzzyRows(const void pa, const void pb)
	702	{
[606]	703	float a = (float)pa;
	704	float b = (float)pb;
[349]	705
	706
[606]	707	for (int i = 1; i < fuzDepth; i++)
	708	{
	709	if (a[0][i] > b[0][i])
	710	{
[349]	711
[606]	712	return -1;
	713	}
	714	if (a[0][i] < b[0][i])
	715	{
[349]	716
[606]	717	return 1;
	718	}
	719	}
[349]	720
[606]	721	return 0;
[349]	722	}
	723
[869]	724	void ModelSimil::FuzzyOrder()
[349]	725	{
[869]	726	int i, depth, partInd, prevPartInd, partCount;
[606]	727	for (int mod = 0; mod < 2; mod++)
	728	{
	729	partCount = m_Mod[mod]->getPartCount();
[869]	730	partInd = m_fuzzyNeighb[mod][partCount - 1][0];
	731	m_aDegrees[mod][partInd][FUZZ_DEG] = 0;
	732
	733	for (i = (partCount - 2); i >= 0; i--)
[606]	734	{
[869]	735	prevPartInd = partInd;
	736	partInd = m_fuzzyNeighb[mod][i][0];
	737	m_aDegrees[mod][partInd][FUZZ_DEG] = m_aDegrees[mod][prevPartInd][FUZZ_DEG];
	738	for (depth = 1; depth < fuzzyDepth; depth++)
[606]	739	{
[869]	740	if (m_fuzzyNeighb[mod][i][depth] != m_fuzzyNeighb[mod][i + 1][depth])
[606]	741	{
[869]	742	m_aDegrees[mod][partInd][FUZZ_DEG]++;
[606]	743	break;
	744	}
	745	}
	746	}
	747	}
[349]	748	}
	749
	750	//sort according to fuzzy degree
	751	void ModelSimil::SortFuzzyNeighb()
	752	{
[606]	753	fuzDepth = fuzzyDepth;
	754	for (int mod = 0; mod < 2; mod++)
	755	{
	756	std::qsort(m_fuzzyNeighb[mod], (size_t)m_Mod[mod]->getPartCount(), sizeof(m_fuzzyNeighb[mod][0]), cmpFuzzyRows);
	757	}
[349]	758	}
	759
	760	//computes fuzzy vertex degree
	761	void ModelSimil::CountFuzzyNeighb()
	762	{
[606]	763	assert(m_aDegrees[0] != NULL);
	764	assert(m_aDegrees[1] != NULL);
[349]	765
[606]	766	assert(m_Neighbours[0] != NULL);
	767	assert(m_Neighbours[1] != NULL);
[349]	768
[606]	769	assert(m_fuzzyNeighb[0] != NULL);
	770	assert(m_fuzzyNeighb[1] != NULL);
[349]	771
[606]	772	std::vector<int> nIndexes;
	773	float newDeg = 0;
[349]	774
[606]	775	for (int mod = 0; mod < 2; mod++)
	776	{
	777	int partCount = m_Mod[mod]->getPartCount();
[349]	778
[606]	779	for (int depth = 0; depth < fuzzyDepth; depth++)
	780	{
	781	//use first column for storing indices
	782	for (int partInd = 0; partInd < partCount; partInd++)
	783	{
	784	if (depth == 0)
	785	{
	786	m_fuzzyNeighb[mod][partInd][depth] = partInd;
	787	}
	788	else if (depth == 1)
	789	{
	790	m_fuzzyNeighb[mod][partInd][depth] = m_aDegrees[mod][partInd][DEGREE];
	791	}
	792	else
	793	{
	794	GetNeighbIndexes(mod, partInd, nIndexes);
[361]	795	for (unsigned int k = 0; k < nIndexes.size(); k++)
[606]	796	{
	797	newDeg += m_fuzzyNeighb[mod][nIndexes.at(k)][depth - 1];
	798	}
	799	newDeg /= nIndexes.size();
	800	m_fuzzyNeighb[mod][partInd][depth] = newDeg;
	801	for (int mod = 0; mod < 2; mod++)
	802	{
	803	int partCount = m_Mod[mod]->getPartCount();
	804	for (int i = partCount - 1; i >= 0; i--)
	805	{
[349]	806
[606]	807	}
	808	}
	809	newDeg = 0;
	810	}
	811	}
	812	}
	813	}
[349]	814
[606]	815	SortFuzzyNeighb();
[869]	816	FuzzyOrder();
[349]	817	}
	818
	819	/** Gets information about Parts' positions in 3D from models into the arrays
[606]	820	m_aPositions.
	821	Assumptions:
	822	- Models (m_Mod) are created and available.
	823	- Arrays m_aPositions are created.
	824	@return 1 if success, otherwise 0.
	825	*/
[349]	826	int ModelSimil::GetPartPositions()
	827	{
[606]	828	// sprawdz zalozenie - o modelach
	829	assert((m_Mod[0] != NULL) && (m_Mod[1] != NULL));
	830	assert(m_Mod[0]->isValid() && m_Mod[1]->isValid());
[349]	831
[606]	832	// sprawdz zalozenie - o tablicach m_aPositions
	833	assert(m_aPositions[0] != NULL);
	834	assert(m_aPositions[1] != NULL);
[349]	835
[606]	836	// dla każdego stworzenia skopiuj informację o polozeniu jego Parts
	837	// do tablic m_aPositions
	838	Part *pPart;
	839	int iMod; // licznik modeli (organizmow)
	840	int iPart; // licznik Parts
	841	for (iMod = 0; iMod < 2; iMod++)
	842	{
	843	// dla każdego z modeli iMod
	844	for (iPart = 0; iPart < m_Mod[iMod]->getPartCount(); iPart++)
	845	{
	846	// dla każdego iPart organizmu iMod
	847	// pobierz tego Part
	848	pPart = m_Mod[iMod]->getPart(iPart);
	849	// zapamietaj jego pozycje 3D w tablicy m_aPositions
	850	m_aPositions[iMod][iPart].x = pPart->p.x;
	851	m_aPositions[iMod][iPart].y = pPart->p.y;
	852	m_aPositions[iMod][iPart].z = pPart->p.z;
	853	}
	854	}
[349]	855
[606]	856	return 1;
[349]	857	}
	858
	859	/** Computes numbers of neurons and neurons' inputs for each Part of each
[606]	860	organisms and fills in the m_aDegrees array.
	861	Assumptions:
	862	- Models (m_Mod) are created and available.
	863	- Arrays m_aDegrees are created.
	864	*/
[349]	865	int ModelSimil::CountPartNeurons()
	866	{
[606]	867	// sprawdz zalozenie - o modelach
	868	assert((m_Mod[0] != NULL) && (m_Mod[1] != NULL));
	869	assert(m_Mod[0]->isValid() && m_Mod[1]->isValid());
[349]	870
[606]	871	// sprawdz zalozenie - o tablicach
	872	assert(m_aDegrees[0] != NULL);
	873	assert(m_aDegrees[1] != NULL);
[349]	874
[606]	875	Part *P1;
	876	Joint *J1;
	877	int i, j, i2, neuro_connections;
[349]	878
[606]	879	// dla obu stworzen oblicz liczbe Neurons + connections dla Parts
	880	// a takze dla OnJoints i Anywhere
	881	for (i = 0; i < 2; i++)
	882	{
	883	for (j = 0; j < m_Mod[i]->getNeuroCount(); j++)
	884	{
	885	// pobierz kolejny Neuron
	886	Neuro *N = m_Mod[i]->getNeuro(j);
	887	// policz liczbe jego wejść i jego samego tez
	888	// czy warto w ogole liczyc polaczenia...? co to da/spowoduje?
	889	neuro_connections = N->getInputCount() + 1;
	890	// wez Part, na ktorym jest Neuron
	891	P1 = N->getPart();
	892	if (P1)
	893	{
	894	// dla neuronow osadzonych na Partach
	895	i2 = m_Mod[i]->findPart(P1); // znajdz indeks Part w Modelu
	896	m_aDegrees[i][i2][2] += neuro_connections; // zwieksz liczbe Connections+Neurons dla tego Part (TDN[2])
	897	m_aDegrees[i][i2][3]++; // zwieksz liczbe Neurons dla tego Part (TDN[3])
	898	}
	899	else
	900	{
	901	// dla neuronow nie osadzonych na partach
	902	J1 = N->getJoint();
	903	if (J1)
	904	{
	905	// dla tych na Jointach
	906	m_aOnJoint[i][2] += neuro_connections; // zwieksz liczbe Connections+Neurons
	907	m_aOnJoint[i][3]++; // zwieksz liczbe Neurons
	908	}
	909	else
	910	{
	911	// dla tych "gdziekolwiek"
	912	m_aAnywhere[i][2] += neuro_connections; // zwieksz liczbe Connections+Neurons
	913	m_aAnywhere[i][3]++; // zwieksz liczbe Neurons
	914	}
	915	}
	916	}
	917	}
	918	return 1;
[349]	919	}
	920
	921	/** Sorts arrays m_aDegrees (for each organism) by Part's degree, and then by
[606]	922	number of neural connections and neurons in groups of Parts with the same
	923	degree.
	924	Assumptions:
	925	- Models (m_Mod) are created and available.
	926	- Arrays m_aDegrees are created.
	927	@saeDegrees, CompareItemNo
	928	*/
[349]	929	int ModelSimil::SortPartInfoTables()
	930	{
[606]	931	// sprawdz zalozenie - o modelach
	932	assert((m_Mod[0] != NULL) && (m_Mod[1] != NULL));
	933	assert(m_Mod[0]->isValid() && m_Mod[1]->isValid());
[349]	934
[606]	935	// sprawdz zalozenie - o tablicach
	936	assert(m_aDegrees[0] != NULL);
	937	assert(m_aDegrees[1] != NULL);
[349]	938
[606]	939	int i;
[869]	940	int(pfDegreeFunction) (const void, const void*) = NULL;
[872]	941	pfDegreeFunction = isFuzzy ? &CompareFuzzyDegrees : &CompareDegrees;
[606]	942	// sortowanie obu tablic wg stopni punktów - TDN[1]
[869]	943	for (i = 0; i < 2; i++)
[606]	944	{
[869]	945	DB(_PrintDegrees(i));
	946	std::qsort(m_aDegrees[i], (size_t)(m_Mod[i]->getPartCount()),
	947	sizeof(TDN), pfDegreeFunction);
	948	DB(_PrintDegrees(i));
[606]	949	}
[349]	950
[606]	951	// sprawdzenie wartosci parametru
	952	DB(i = sizeof(TDN);)
[872]	953	int degreeType = isFuzzy ? FUZZ_DEG : DEGREE;
[349]	954
[606]	955	// sortowanie obu tablic m_aDegrees wedlug liczby neuronów i
	956	// czesci neuronu - ale w obrebie grup o tym samym stopniu
	957	for (i = 0; i < 2; i++)
	958	{
	959	int iPocz = 0;
	960	int iDeg, iNewDeg, iPartCount, j;
	961	// stopien pierwszego punktu w tablicy Degrees odniesienie
	962	iDeg = m_aDegrees[i][0][degreeType];
	963	iPartCount = m_Mod[i]->getPartCount();
	964	// po kolei dla kazdego punktu w organizmie
	965	for (j = 0; j <= iPartCount; j++)
	966	{
	967	// sprawdz stopien punktu (lub nadaj 0 - gdy juz koniec tablicy)
	968	// iNewDeg = (j != iPartCount) ? m_aDegrees[i][j][1] : 0;
	969	// usunieto stara wersje porownania!!! wprowadzono znak porownania <
[349]	970
[606]	971	iNewDeg = (j < iPartCount) ? m_aDegrees[i][j][degreeType] : 0;
	972	// skoro tablice sa posortowane wg stopni, to mamy na pewno taka kolejnosc
	973	assert(iNewDeg <= iDeg);
	974	if (iNewDeg != iDeg)
	975	{
	976	// gdy znaleziono koniec grupy o tym samym stopniu
	977	// sortuj po liczbie neuronow w obrebie grupy
	978	DB(_PrintDegrees(i));
	979	DB(printf("qsort( poczatek=%i, rozmiar=%i, sizeof(TDN)=%i)\n", iPocz, (j - iPocz), sizeof(TDN));)
	980	// wyswietlamy z jedna komorka po zakonczeniu tablicy
	981	DB(_PrintArray(m_aDegrees[i][iPocz], 0, (j - iPocz) * 4);)
[349]	982
[606]	983	std::qsort(m_aDegrees[i][iPocz], (size_t)(j - iPocz),
[869]	984	sizeof(TDN), ModelSimil::CompareConnsNo);
[606]	985	DB(_PrintDegrees(i));
	986	// wyswietlamy z jedna komorka po zakonczeniu tablicy
	987	DB(_PrintArray(m_aDegrees[i][iPocz], 0, (j - iPocz) * 4);)
	988	// rozpocznij nowa grupe
	989	iPocz = j;
	990	iDeg = iNewDeg;
	991	}
	992	}
	993	}
	994	return 1;
[349]	995	}
	996
	997
	998	/** Prints the state of the matching object. Debug method.
	999	*/
	1000	void ModelSimil::_PrintPartsMatching()
	1001	{
[606]	1002	// assure that matching exists
	1003	assert(m_pMatching != NULL);
[349]	1004
[606]	1005	printf("Parts matching:\n");
	1006	m_pMatching->PrintMatching();
[349]	1007	}
	1008
	1009	void ModelSimil::ComputeMatching()
	1010	{
[606]	1011	// uniwersalne liczniki
	1012	int i, j;
[349]	1013
[606]	1014	assert(m_pMatching != NULL);
	1015	assert(m_pMatching->IsEmpty() == true);
[349]	1016
[606]	1017	// rozpoczynamy etap dopasowywania Parts w organizmach
	1018	// czy dopasowano już wszystkie Parts?
	1019	int iCzyDopasowane = 0;
	1020	// koniec grupy aktualnie dopasowywanej w każdym organizmie
	1021	int aiKoniecGrupyDopasowania[2] = { 0, 0 };
	1022	// koniec grupy już w całości dopasowanej
	1023	// (Pomiedzy tymi dwoma indeksami znajduja sie Parts w tablicy
	1024	// m_aDegrees, ktore moga byc dopasowywane (tam nadal moga
	1025	// byc tez dopasowane - ale nie musi to byc w sposob
	1026	// ciagly)
	1027	int aiKoniecPierwszejGrupy[2] = { 0, 0 };
	1028	// Tablica przechowująca odległości poszczególnych Parts z aktualnych
	1029	// grup dopasowania. Rozmiar - prostokąt o bokach równych liczbie elementów w
	1030	// dopasowywanych aktualnie grupach. Pierwszy wymiar - pierwszy organizm.
	1031	// Drugi wymiar - drugi organizm (nie zależy to od tego, który jest mniejszy).
	1032	// Wliczane w rozmiar tablicy są nawet już dopasowane elementy - tablice
	1033	// paiCzyDopasowany pamiętają stan dopasowania tych elementów.
	1034	typedef double *TPDouble;
	1035	double **aadOdleglosciParts;
	1036	// dwie tablice okreslajace Parts, ktore moga byc do siebie dopasowywane
	1037	// rozmiary: [0] - aiRozmiarCalychGrup[1]
	1038	// [1] - aiRozmiarCalychGrup[0]
	1039	std::vector<bool> *apvbCzyMinimalnaOdleglosc[2];
	1040	// rozmiar aktualnie dopasowywanej grupy w odpowiednim organizmie (tylko elementy
	1041	// jeszcze niedopasowane).
	1042	int aiRozmiarGrupy[2];
	1043	// rozmiar aktualnie dopasowywanych grup w odpowiednim organizmie (włączone są
	1044	// w to również elementy już dopasowane).
	1045	int aiRozmiarCalychGrup[2] = { 0, 0 };
[349]	1046
[606]	1047	// utworzenie tablicy rozmiarow
	1048	for (i = 0; i < 2; i++)
	1049	{
	1050	m_aiPartCount[i] = m_Mod[i]->getPartCount();
	1051	}
[349]	1052
[606]	1053	// DOPASOWYWANIE PARTS
	1054	while (!iCzyDopasowane)
	1055	{
	1056	// znajdz konce obu grup aktualnie dopasowywanych w obu organizmach
	1057	for (i = 0; i < 2; i++)
	1058	{
	1059	// czyli poszukaj miejsca zmiany stopnia lub konca tablicy
	1060	for (j = aiKoniecPierwszejGrupy[i] + 1; j < m_aiPartCount[i]; j++)
	1061	{
	1062	if (m_aDegrees[i][j][DEGREE] < m_aDegrees[i][j - 1][DEGREE])
	1063	{
	1064	break;
	1065	}
	1066	}
	1067	aiKoniecGrupyDopasowania[i] = j;
[349]	1068
[606]	1069	// sprawdz poprawnosc tego indeksu
	1070	assert((aiKoniecGrupyDopasowania[i] > 0) &&
	1071	(aiKoniecGrupyDopasowania[i] <= m_aiPartCount[i]));
[349]	1072
[606]	1073	// oblicz rozmiary całych grup - łącznie z dopasowanymi już elementami
	1074	aiRozmiarCalychGrup[i] = aiKoniecGrupyDopasowania[i] -
	1075	aiKoniecPierwszejGrupy[i];
[349]	1076
[606]	1077	// sprawdz teraz rozmiar tej grupy w sensie liczby niedopasowanzch
	1078	// nie moze to byc puste!
	1079	aiRozmiarGrupy[i] = 0;
	1080	for (j = aiKoniecPierwszejGrupy[i]; j < aiKoniecGrupyDopasowania[i]; j++)
	1081	{
	1082	// od poczatku do konca grupy
	1083	if (!m_pMatching->IsMatched(i, j))
	1084	{
	1085	// jesli niedopasowany, to zwieksz licznik
	1086	aiRozmiarGrupy[i]++;
	1087	}
	1088	}
	1089	// grupa nie moze byc pusta!
	1090	assert(aiRozmiarGrupy[i] > 0);
	1091	}
[349]	1092
[606]	1093	// DOPASOWYWANIE PARTS Z GRUP
[349]	1094
[606]	1095	// stworzenie tablicy odległości lokalnych
	1096	// stwórz pierwszy wymiar - wg rozmiaru zerowego organizmu
	1097	aadOdleglosciParts = new TPDouble[aiRozmiarCalychGrup[0]];
	1098	assert(aadOdleglosciParts != NULL);
	1099	// stwórz drugi wymiar - wg rozmiaru drugiego organizmu
	1100	for (i = 0; i < aiRozmiarCalychGrup[0]; i++)
	1101	{
	1102	aadOdleglosciParts[i] = new double[aiRozmiarCalychGrup[1]];
	1103	assert(aadOdleglosciParts[i] != NULL);
	1104	}
[349]	1105
[606]	1106	// stworzenie tablic mozliwosci dopasowania (indykatorow minimalnej odleglosci)
	1107	apvbCzyMinimalnaOdleglosc[0] = new std::vector<bool>(aiRozmiarCalychGrup[1], false);
	1108	apvbCzyMinimalnaOdleglosc[1] = new std::vector<bool>(aiRozmiarCalychGrup[0], false);
	1109	// sprawdz stworzenie tablic
	1110	assert(apvbCzyMinimalnaOdleglosc[0] != NULL);
	1111	assert(apvbCzyMinimalnaOdleglosc[1] != NULL);
[349]	1112
[606]	1113	// wypełnienie elementów macierzy (i,j) odległościami pomiędzy
	1114	// odpowiednimi Parts: (i) w organizmie 0 i (j) w organizmie 1.
	1115	// UWAGA! Uwzględniamy tylko te Parts, które nie są jeszcze dopasowane
	1116	// (reszta to byłaby po prostu strata czasu).
	1117	int iDeg, iNeu; // ilościowe określenie różnic stopnia, liczby neuronów i połączeń Parts
	1118	//int iNIt;
	1119	double dGeo; // ilościowe określenie różnic geometrycznych pomiędzy Parts
	1120	// indeksy konkretnych Parts - indeksy sa ZERO-BASED, choć właściwy dostep
	1121	// do informacji o Part wymaga dodania aiKoniecPierwszejGrupy[]
	1122	// tylko aadOdleglosciParts[][] indeksuje sie bezposrednio zawartoscia iIndex[]
	1123	int iIndex[2];
	1124	int iPartIndex[2] = { -1, -1 }; // at [iModel]: original index of a Part for the specified model (iModel)
[349]	1125
[606]	1126	// debug - wypisz zakres dopasowywanych indeksow
	1127	DB(printf("Organizm 0: grupa: (%2i, %2i)\n", aiKoniecPierwszejGrupy[0],
	1128	aiKoniecGrupyDopasowania[0]);)
	1129	DB(printf("Organizm 1: grupa: (%2i, %2i)\n", aiKoniecPierwszejGrupy[1],
[869]	1130	aiKoniecGrupyDopasowania[1]);)
[349]	1131
[606]	1132	for (i = 0; i < aiRozmiarCalychGrup[0]; i++)
	1133	{
[349]	1134
[869]	1135	// iterujemy i - Parts organizmu 0
	1136	// (indeks podstawowy - aiKoniecPierwszejGrupy[0])
[349]	1137
[869]	1138	if (!(m_pMatching->IsMatched(0, aiKoniecPierwszejGrupy[0] + i)))
[606]	1139	{
[869]	1140	// interesuja nas tylko te niedopasowane jeszcze (i)
	1141	for (j = 0; j < aiRozmiarCalychGrup[1]; j++)
	1142	{
[349]	1143
[869]	1144	// iterujemy j - Parts organizmu 1
	1145	// (indeks podstawowy - aiKoniecPierwszejGrupy[1])
[349]	1146
[869]	1147	if (!(m_pMatching->IsMatched(1, aiKoniecPierwszejGrupy[1] + j)))
	1148	{
	1149	// interesuja nas tylko te niedopasowane jeszcze (j)
	1150	// teraz obliczymy lokalne różnice pomiędzy Parts
	1151	iDeg = abs(m_aDegrees[0][aiKoniecPierwszejGrupy[0] + i][1]
	1152	- m_aDegrees[1][aiKoniecPierwszejGrupy[1] + j][1]);
[349]	1153
[869]	1154	//iNit currently is not a component of distance measure
	1155	//iNIt = abs(m_aDegrees[0][ aiKoniecPierwszejGrupy[0] + i ][2]
	1156	// - m_aDegrees[1][ aiKoniecPierwszejGrupy[1] + j ][2]);
[349]	1157
[869]	1158	iNeu = abs(m_aDegrees[0][aiKoniecPierwszejGrupy[0] + i][3]
	1159	- m_aDegrees[1][aiKoniecPierwszejGrupy[1] + j][3]);
[349]	1160
[869]	1161	// obliczenie także lokalnych różnic geometrycznych pomiędzy Parts
	1162	// find original indices of Parts for both of the models
	1163	iPartIndex[0] = m_aDegrees[0][aiKoniecPierwszejGrupy[0] + i][0];
	1164	iPartIndex[1] = m_aDegrees[1][aiKoniecPierwszejGrupy[1] + j][0];
	1165	// now compute the geometrical distance of these Parts (use m_aPositions
	1166	// which should be computed by SVD)
	1167	Pt3D Part0Pos(m_aPositions[0][iPartIndex[0]]);
	1168	Pt3D Part1Pos(m_aPositions[1][iPartIndex[1]]);
	1169	dGeo = m_adFactors[3] == 0 ? 0 : Part0Pos.distanceTo(Part1Pos); //no need to compute distane when m_dDG weight is 0
[349]	1170
[869]	1171	// tutaj prawdopodobnie należy jeszcze dodać sprawdzanie
	1172	// identyczności pozostałych własności (oprócz geometrii)
	1173	// - żeby móc stwierdzić w ogóle identyczność Parts
[349]	1174
[869]	1175	// i ostateczna odleglosc indukowana przez te roznice
	1176	// (tutaj nie ma różnicy w liczbie wszystkich wierzchołków)
	1177	aadOdleglosciParts[i][j] = m_adFactors[1] * double(iDeg) +
	1178	m_adFactors[2] * double(iNeu) +
	1179	m_adFactors[3] * dGeo;
	1180	DB(printf("Parts Distance (%2i,%2i) = %.3lf\n", aiKoniecPierwszejGrupy[0] + i,
	1181	aiKoniecPierwszejGrupy[1] + j, aadOdleglosciParts[i][j]);)
	1182	DB(printf("Parts geometrical distance = %.20lf\n", dGeo);)
	1183	DB(printf("Pos0: (%.3lf %.3lf %.3lf)\n", Part0Pos.x, Part0Pos.y, Part0Pos.z);)
	1184	DB(printf("Pos1: (%.3lf %.3lf %.3lf)\n", Part1Pos.x, Part1Pos.y, Part1Pos.z);)
	1185	}
[606]	1186	}
	1187	}
	1188	}
[349]	1189
[606]	1190	// tutaj - sprawdzic tylko, czy tablica odleglosci lokalnych jest poprawnie obliczona
[349]	1191
[606]	1192	// WYKORZYSTANIE TABLICY ODLEGLOSCI DO BUDOWY DOPASOWANIA
[349]	1193
[606]	1194	// trzeba raczej iterować aż do zebrania wszystkich możliwych dopasowań w grupie
	1195	// dlatego wprowadzamy dodatkowa zmienna - czy skonczyla sie juz grupa
	1196	bool bCzyKoniecGrupy = false;
	1197	while (!bCzyKoniecGrupy)
	1198	{
	1199	for (i = 0; i < 2; i++)
	1200	{
	1201	// iterujemy (i) po organizmach
	1202	// szukamy najpierw jakiegoś niedopasowanego jeszcze Part w organizmach
[349]	1203
[606]	1204	// zakładamy, że nie ma takiego Part
	1205	iIndex[i] = -1;
[349]	1206
[606]	1207	for (j = 0; j < aiRozmiarCalychGrup[i]; j++)
	1208	{
	1209	// iterujemy (j) - Parts organizmu (i)
	1210	// (indeks podstawowy - aiKoniecPierwszejGrupy[0])
	1211	if (!(m_pMatching->IsMatched(i, aiKoniecPierwszejGrupy[i] + j)))
	1212	{
	1213	// gdy mamy w tej grupie jakis niedopasowany element, to ustawiamy
	1214	// iIndex[i] (chcemy w zasadzie pierwszy taki)
	1215	iIndex[i] = j;
	1216	break;
	1217	}
	1218	}
[349]	1219
[606]	1220	// sprawdzamy, czy w ogole znaleziono taki Part
	1221	if (iIndex[i] < 0)
	1222	{
	1223	// gdy nie znaleziono takiego Part - mamy koniec dopasowywania w
	1224	// tych grupach
	1225	bCzyKoniecGrupy = true;
	1226	}
	1227	// sprawdz poprawnosc indeksu niedopasowanego Part - musi byc w aktualnej grupie
	1228	assert((iIndex[i] >= -1) && (iIndex[i] < aiRozmiarCalychGrup[i]));
	1229	}
[349]	1230
	1231
[606]	1232	// teraz mamy sytuacje:
	1233	// - mamy w iIndex[0] i iIndex[1] indeksy pierwszych niedopasowanych Part
	1234	// w organizmach, albo
	1235	// - nie ma w ogóle już czego dopasowywać (należy przejść do innej grupy)
	1236	if (!bCzyKoniecGrupy)
	1237	{
	1238	// gdy nie ma jeszcze końca żadnej z grup - możemy dopasowywać
	1239	// najpierw wyszukujemy w tablicy minimum odległości od tych
	1240	// wyznaczonych Parts
[349]	1241
[606]	1242	// najpierw wyczyscimy wektory potencjalnych dopasowan
	1243	// dla organizmu 1 (o rozmiarze grupy z 0)
	1244	for (i = 0; i < aiRozmiarCalychGrup[0]; i++)
	1245	{
	1246	apvbCzyMinimalnaOdleglosc[1]->operator[](i) = false;
	1247	}
	1248	// dla organizmu 0 (o rozmiarze grup z 1)
	1249	for (i = 0; i < aiRozmiarCalychGrup[1]; i++)
	1250	{
	1251	apvbCzyMinimalnaOdleglosc[0]->operator[](i) = false;
	1252	}
[349]	1253
[606]	1254	// szukanie minimum dla Part o indeksie iIndex[0] w organizmie 0
	1255	// wsrod niedopasowanych Parts z organizmu 1
	1256	// zakładamy, że nie znaleliśmy jeszcze minimum
	1257	double dMinimum = HUGE_VAL;
	1258	for (i = 0; i < aiRozmiarCalychGrup[1]; i++)
	1259	{
	1260	if (!(m_pMatching->IsMatched(1, aiKoniecPierwszejGrupy[1] + i)))
	1261	{
[349]	1262
[606]	1263	// szukamy minimum obliczonej lokalnej odleglosci tylko wsrod
	1264	// niedopasowanych jeszcze Parts
	1265	if (aadOdleglosciParts[iIndex[0]][i] < dMinimum)
	1266	{
	1267	dMinimum = aadOdleglosciParts[iIndex[0]][i];
	1268	}
[349]	1269
[606]	1270	// przy okazji - sprawdz, czy HUGE_VAL jest rzeczywiscie max dla double
	1271	assert(aadOdleglosciParts[iIndex[0]][i] < HUGE_VAL);
	1272	}
	1273	}
	1274	// sprawdz, czy minimum znaleziono - musi takie byc, bo jest cos niedopasowanego
	1275	assert((dMinimum >= 0.0) && (dMinimum < HUGE_VAL));
[349]	1276
[606]	1277	// teraz zaznaczamy w tablicy te wszystkie Parts, ktore maja
	1278	// rzeczywiscie te minimalna odleglosc do Part iIndex[0] w organizmie 0
	1279	for (i = 0; i < aiRozmiarCalychGrup[1]; i++)
	1280	{
	1281	if (!(m_pMatching->IsMatched(1, aiKoniecPierwszejGrupy[1] + i)))
	1282	{
	1283	if (aadOdleglosciParts[iIndex[0]][i] == dMinimum)
	1284	{
	1285	// jesli w danym miejscu tablicy odleglosci jest faktyczne
	1286	// minimum odleglosci, to mamy potencjalna pare dla iIndex[0]
	1287	apvbCzyMinimalnaOdleglosc[0]->operator[](i) = true;
	1288	}
[349]	1289
[606]	1290	// sprawdz poprawnosc znalezionego wczesniej minimum
	1291	assert(aadOdleglosciParts[iIndex[0]][i] >= dMinimum);
	1292	}
	1293	}
[349]	1294
[606]	1295	// podobnie szukamy minimum dla Part o indeksie iIndex[1] w organizmie 1
	1296	// wsrod niedopasowanych Parts z organizmu 0
	1297	dMinimum = HUGE_VAL;
	1298	for (i = 0; i < aiRozmiarCalychGrup[0]; i++)
	1299	{
	1300	if (!(m_pMatching->IsMatched(0, aiKoniecPierwszejGrupy[0] + i)))
	1301	{
	1302	// szukamy minimum obliczonej lokalnej odleglosci tylko wsrod
	1303	// niedopasowanych jeszcze Parts
	1304	if (aadOdleglosciParts[i][iIndex[1]] < dMinimum)
	1305	{
	1306	dMinimum = aadOdleglosciParts[i][iIndex[1]];
	1307	}
	1308	// przy okazji - sprawdz, czy HUGE_VAL jest rzeczywiscie max dla double
	1309	assert(aadOdleglosciParts[i][iIndex[1]] < HUGE_VAL);
	1310	}
	1311	}
	1312	// sprawdz, czy minimum znaleziono - musi takie byc, bo jest cos niedopasowanego
	1313	assert((dMinimum >= 0.0) && (dMinimum < HUGE_VAL));
[349]	1314
[606]	1315	// teraz zaznaczamy w tablicy te wszystkie Parts, ktore maja
	1316	// rzeczywiscie te minimalne odleglosci do Part iIndex[1] w organizmie 1
	1317	for (i = 0; i < aiRozmiarCalychGrup[0]; i++)
	1318	{
	1319	if (!(m_pMatching->IsMatched(0, aiKoniecPierwszejGrupy[0] + i)))
	1320	{
	1321	if (aadOdleglosciParts[i][iIndex[1]] == dMinimum)
	1322	{
	1323	// jesli w danym miejscu tablicy odleglosci jest faktyczne
	1324	// minimum odleglosci, to mamy potencjalna pare dla iIndex[1]
	1325	apvbCzyMinimalnaOdleglosc[1]->operator[](i) = true;
	1326	}
[349]	1327
[606]	1328	// sprawdz poprawnosc znalezionego wczesniej minimum
	1329	assert(aadOdleglosciParts[i][iIndex[1]] >= dMinimum);
	1330	}
	1331	}
[349]	1332
[606]	1333	// teraz mamy juz poszukane potencjalne grupy dopasowania - musimy
	1334	// zdecydowac, co do czego dopasujemy!
	1335	// szukamy Part iIndex[0] posrod mozliwych do dopasowania dla Part iIndex[1]
	1336	// szukamy takze Part iIndex[1] posrod mozliwych do dopasowania dla Part iIndex[0]
	1337	bool PartZ1NaLiscie0 = apvbCzyMinimalnaOdleglosc[0]->operator[](iIndex[1]);
	1338	bool PartZ0NaLiscie1 = apvbCzyMinimalnaOdleglosc[1]->operator[](iIndex[0]);
[349]	1339
[606]	1340	if (PartZ1NaLiscie0 && PartZ0NaLiscie1)
	1341	{
	1342	// PRZYPADEK 1. Oba Parts maja sie wzajemnie na listach mozliwych
	1343	// dopasowan.
	1344	// AKCJA. Dopasowanie wzajemne do siebie.
[349]	1345
[606]	1346	m_pMatching->Match(0, aiKoniecPierwszejGrupy[0] + iIndex[0],
	1347	1, aiKoniecPierwszejGrupy[1] + iIndex[1]);
[349]	1348
[606]	1349	// zmniejsz liczby niedopasowanych elementow w grupach
	1350	aiRozmiarGrupy[0]--;
	1351	aiRozmiarGrupy[1]--;
	1352	// debug - co zostalo dopasowane
	1353	DB(printf("Przypadek 1.: dopasowane Parts: (%2i, %2i)\n", aiKoniecPierwszejGrupy[0]
	1354	+ iIndex[0], aiKoniecPierwszejGrupy[1] + iIndex[1]);)
[349]	1355
[606]	1356	}// PRZYPADEK 1.
	1357	else
	1358	if (PartZ1NaLiscie0 \|\| PartZ0NaLiscie1)
	1359	{
[869]	1360	// PRZYPADEK 2. Tylko jeden z Parts ma drugiego na swojej liscie
	1361	// mozliwych dopasowan
	1362	// AKCJA. Dopasowanie jednego jest proste (tego, ktory nie ma
	1363	// na swojej liscie drugiego). Dla tego drugiego nalezy powtorzyc
	1364	// duza czesc obliczen (znalezc mu nowa mozliwa pare)
[349]	1365
[869]	1366	// indeks organizmu, ktorego Part nie ma dopasowywanego Part
	1367	// z drugiego organizmu na swojej liscie
	1368	int iBezDrugiego;
[349]	1369
[869]	1370	// okreslenie indeksu organizmu z dopasowywanym Part
	1371	if (!PartZ1NaLiscie0)
	1372	{
	1373	iBezDrugiego = 0;
	1374	}
	1375	else
	1376	{
	1377	iBezDrugiego = 1;
	1378	}
	1379	// sprawdz indeks organizmu
	1380	assert((iBezDrugiego == 0) \|\| (iBezDrugiego == 1));
[349]	1381
[869]	1382	int iDopasowywany = -1;
	1383	// poszukujemy pierwszego z listy dopasowania
	1384	for (i = 0; i < aiRozmiarCalychGrup[1 - iBezDrugiego]; i++)
[606]	1385	{
[869]	1386	if (apvbCzyMinimalnaOdleglosc[iBezDrugiego]->operator[](i))
	1387	{
	1388	iDopasowywany = i;
	1389	break;
	1390	}
[606]	1391	}
[869]	1392	// sprawdz poprawnosc indeksu dopasowywanego (musimy cos znalezc!)
	1393	// nieujemny i w odpowiedniej grupie!
	1394	assert((iDopasowywany >= 0) &&
	1395	(iDopasowywany < aiRozmiarCalychGrup[1 - iBezDrugiego]));
[349]	1396
[869]	1397	// znalezlismy 1. Part z listy dopasowania - dopasowujemy!
	1398	m_pMatching->Match(
	1399	iBezDrugiego,
	1400	aiKoniecPierwszejGrupy[iBezDrugiego] + iIndex[iBezDrugiego],
	1401	1 - iBezDrugiego,
	1402	aiKoniecPierwszejGrupy[1 - iBezDrugiego] + iDopasowywany);
	1403	DB(printf("Przypadek 2.1.: dopasowane Parts dopasowanie bez drugiego: (%2i, %2i)\n", aiKoniecPierwszejGrupy[iBezDrugiego] + iIndex[iBezDrugiego],
	1404	aiKoniecPierwszejGrupy[1 - iBezDrugiego] + iDopasowywany);)
[349]	1405
[869]	1406	// zmniejsz liczby niedopasowanych elementow w grupach
	1407	aiRozmiarGrupy[0]--;
	1408	aiRozmiarGrupy[1]--;
[349]	1409
[869]	1410	// sprawdz, czy jest szansa dopasowania tego Part z drugiej strony
	1411	// (ktora miala mozliwosc dopasowania tego Part z poprzedniego organizmu)
	1412	if ((aiRozmiarGrupy[0] > 0) && (aiRozmiarGrupy[1] > 0))
	1413	{
	1414	// jesli grupy sie jeszcze nie wyczrpaly
	1415	// to jest mozliwosc dopasowania w organizmie
[349]	1416
[869]	1417	int iZDrugim = 1 - iBezDrugiego;
	1418	// sprawdz indeks organizmu
	1419	assert((iZDrugim == 0) \|\| (iZDrugim == 1));
[349]	1420
[869]	1421	// bedziemy szukac minimum wsrod niedopasowanych - musimy wykasowac
	1422	// poprzednie obliczenia minimum
	1423	// dla organizmu 1 (o rozmiarze grupy z 0)
	1424	for (i = 0; i < aiRozmiarCalychGrup[0]; i++)
	1425	{
	1426	apvbCzyMinimalnaOdleglosc[1]->operator[](i) = false;
	1427	}
	1428	// dla organizmu 0 (o rozmiarze grup z 1)
	1429	for (i = 0; i < aiRozmiarCalychGrup[1]; i++)
	1430	{
	1431	apvbCzyMinimalnaOdleglosc[0]->operator[](i) = false;
	1432	}
[349]	1433
[869]	1434	// szukamy na nowo minimum dla Part o indeksie iIndex[ iZDrugim ] w organizmie iZDrugim
	1435	// wsrod niedopasowanych Parts z organizmu 1 - iZDrugim
	1436	dMinimum = HUGE_VAL;
	1437	for (i = 0; i < aiRozmiarCalychGrup[1 - iZDrugim]; i++)
[606]	1438	{
[869]	1439	if (!(m_pMatching->IsMatched(
	1440	1 - iZDrugim,
	1441	aiKoniecPierwszejGrupy[1 - iZDrugim] + i)))
[606]	1442	{
[869]	1443	// szukamy minimum obliczonej lokalnej odleglosci tylko wsrod
	1444	// niedopasowanych jeszcze Parts
	1445	if (iZDrugim == 0)
[606]	1446	{
[869]	1447	// teraz niestety musimy rozpoznac odpowiedni organizm
	1448	// zeby moc indeksowac niesymetryczna tablice
	1449	if (aadOdleglosciParts[iIndex[0]][i] < dMinimum)
	1450	{
	1451	dMinimum = aadOdleglosciParts[iIndex[0]][i];
	1452	}
	1453	// przy okazji - sprawdz, czy HUGE_VAL jest rzeczywiscie max dla double
	1454	assert(aadOdleglosciParts[iIndex[0]][i] < HUGE_VAL);
	1455
[606]	1456	}
[869]	1457	else
[606]	1458	{
[869]	1459	if (aadOdleglosciParts[i][iIndex[1]] < dMinimum)
	1460	{
	1461	dMinimum = aadOdleglosciParts[i][iIndex[1]];
	1462	}
	1463	// przy okazji - sprawdz, czy HUGE_VAL jest rzeczywiscie max dla double
	1464	assert(aadOdleglosciParts[i][iIndex[1]] < HUGE_VAL);
[606]	1465	}
	1466	}
	1467	}
[869]	1468	// sprawdz, czy minimum znaleziono - musi takie byc, bo jest cos niedopasowanego
	1469	assert((dMinimum >= 0.0) && (dMinimum < HUGE_VAL));
[349]	1470
[869]	1471	// teraz zaznaczamy w tablicy te wszystkie Parts, ktore maja
	1472	// rzeczywiscie te minimalne odleglosci do Part w organizmie iZDrugim
	1473	for (i = 0; i < aiRozmiarCalychGrup[1 - iZDrugim]; i++)
[606]	1474	{
[869]	1475	if (!(m_pMatching->IsMatched(
	1476	1 - iZDrugim,
	1477	aiKoniecPierwszejGrupy[1 - iZDrugim] + i)))
[606]	1478	{
[869]	1479	if (iZDrugim == 0)
[606]	1480	{
[869]	1481	// teraz niestety musimy rozpoznac odpowiedni organizm
	1482	// zeby moc indeksowac niesymetryczna tablice
	1483	if (aadOdleglosciParts[iIndex[0]][i] == dMinimum)
	1484	{
	1485	// jesli w danym miejscu tablicy odleglosci jest faktyczne
	1486	// minimum odleglosci, to mamy potencjalna pare dla iIndex[1]
	1487	apvbCzyMinimalnaOdleglosc[0]->operator[](i) = true;
	1488	}
[606]	1489	}
[869]	1490	else
[606]	1491	{
[869]	1492	if (aadOdleglosciParts[i][iIndex[1]] == dMinimum)
	1493	{
	1494	apvbCzyMinimalnaOdleglosc[1]->operator[](i) = true;
	1495	}
[606]	1496	}
	1497	}
	1498	}
[349]	1499
[869]	1500	// a teraz - po znalezieniu potencjalnych elementow do dopasowania
	1501	// dopasowujemy pierwszy z potencjalnych
	1502	iDopasowywany = -1;
	1503	for (i = 0; i < aiRozmiarCalychGrup[1 - iZDrugim]; i++)
[606]	1504	{
[869]	1505	if (apvbCzyMinimalnaOdleglosc[iZDrugim]->operator[](i))
	1506	{
	1507	iDopasowywany = i;
	1508	break;
	1509	}
[606]	1510	}
[869]	1511	// musielismy znalezc jakiegos dopasowywanego
	1512	assert((iDopasowywany >= 0) &&
	1513	(iDopasowywany < aiRozmiarCalychGrup[1 - iZDrugim]));
[349]	1514
[869]	1515	// no to juz mozemy dopasowac
	1516	m_pMatching->Match(
	1517	iZDrugim,
	1518	aiKoniecPierwszejGrupy[iZDrugim] + iIndex[iZDrugim],
	1519	1 - iZDrugim,
	1520	aiKoniecPierwszejGrupy[1 - iZDrugim] + iDopasowywany);
	1521	DB(printf("Przypadek 2.1.: dopasowane Parts dopasowaniebz drugim: (%2i, %2i)\n", aiKoniecPierwszejGrupy[iZDrugim] + iIndex[iZDrugim], aiKoniecPierwszejGrupy[1 - iZDrugim] + iDopasowywany);)
[349]	1522
[869]	1523	//aiKoniecPierwszejGrupy[ 1-iBezDrugiego ] + iDopasowywany ;)
	1524	//aiKoniecPierwszejGrupy[ 1-iBezDrugiego ] + iDopasowywany ;)
	1525	// zmniejsz liczby niedopasowanych elementow w grupach
	1526	aiRozmiarGrupy[0]--;
	1527	aiRozmiarGrupy[1]--;
[349]	1528
[869]	1529	}
	1530	else
	1531	{
	1532	// jedna z grup sie juz wyczerpala
	1533	// wiec nie ma mozliwosci dopasowania tego drugiego Partu
	1534	/// i trzeba poczekac na zmiane grupy
	1535	}
[349]	1536
[869]	1537	DB(printf("Przypadek 2.\n");)
[349]	1538
[606]	1539	}// PRZYPADEK 2.
	1540	else
	1541	{
	1542	// PRZYPADEK 3. Zaden z Parts nie ma na liscie drugiego
	1543	// AKCJA. Niezalezne dopasowanie obu Parts do pierwszych ze swojej listy
[349]	1544
[606]	1545	// najpierw dopasujemy do iIndex[0] w organizmie 0
	1546	int iDopasowywany = -1;
	1547	// poszukujemy pierwszego z listy dopasowania - w organizmie 1
	1548	for (i = 0; i < aiRozmiarCalychGrup[1]; i++)
	1549	{
	1550	if (apvbCzyMinimalnaOdleglosc[0]->operator[](i))
	1551	{
	1552	iDopasowywany = i;
	1553	break;
	1554	}
	1555	}
	1556	// musielismy znalezc jakiegos dopasowywanego
	1557	assert((iDopasowywany >= 0) &&
	1558	(iDopasowywany < aiRozmiarCalychGrup[1]));
[349]	1559
[606]	1560	// teraz wlasnie dopasowujemy
	1561	m_pMatching->Match(
	1562	0,
	1563	aiKoniecPierwszejGrupy[0] + iIndex[0],
	1564	1,
	1565	aiKoniecPierwszejGrupy[1] + iDopasowywany);
[349]	1566
[606]	1567	// zmniejszamy liczbe niedopasowanych Parts
	1568	aiRozmiarGrupy[0]--;
	1569	aiRozmiarGrupy[1]--;
[349]	1570
[606]	1571	// debug - dopasowanie
	1572	DB(printf("Przypadek 3.: dopasowane Parts: (%2i, %2i)\n", aiKoniecPierwszejGrupy[0]
	1573	+ iIndex[0], aiKoniecPierwszejGrupy[1] + iDopasowywany);)
[349]	1574
[606]	1575	// teraz dopasowujemy do iIndex[1] w organizmie 1
	1576	iDopasowywany = -1;
	1577	// poszukujemy pierwszego z listy dopasowania - w organizmie 0
	1578	for (i = 0; i < aiRozmiarCalychGrup[0]; i++)
	1579	{
	1580	if (apvbCzyMinimalnaOdleglosc[1]->operator[](i))
	1581	{
	1582	iDopasowywany = i;
	1583	break;
	1584	}
	1585	}
	1586	// musielismy znalezc jakiegos dopasowywanego
	1587	assert((iDopasowywany >= 0) &&
	1588	(iDopasowywany < aiRozmiarCalychGrup[0]));
[349]	1589
[606]	1590	// no i teraz realizujemy dopasowanie
	1591	m_pMatching->Match(
	1592	0,
	1593	aiKoniecPierwszejGrupy[0] + iDopasowywany,
	1594	1,
	1595	aiKoniecPierwszejGrupy[1] + iIndex[1]);
[349]	1596
[606]	1597	// zmniejszamy liczbe niedopasowanych Parts
	1598	aiRozmiarGrupy[0]--;
	1599	aiRozmiarGrupy[1]--;
[349]	1600
[606]	1601	// debug - dopasowanie
	1602	DB(printf("Przypadek 3.: dopasowane Parts: (%2i, %2i)\n", aiKoniecPierwszejGrupy[0]
	1603	+ iDopasowywany, aiKoniecPierwszejGrupy[1] + iIndex[1]);)
[349]	1604
	1605
[606]	1606	} // PRZYPADEK 3.
[349]	1607
[606]	1608	}// if (! bCzyKoniecGrupy)
	1609	else
	1610	{
[647]	1611	// gdy mamy juz koniec grup - musimy zlikwidowac tablice aadOdleglosciParts
[606]	1612	// bo za chwilke skonczy sie nam petla
	1613	for (i = 0; i < aiRozmiarCalychGrup[0]; i++)
	1614	{
	1615	SAFEDELETEARRAY(aadOdleglosciParts[i]);
	1616	}
	1617	SAFEDELETEARRAY(aadOdleglosciParts);
[349]	1618
[606]	1619	// musimy tez usunac tablice (wektory) mozliwosci dopasowania
	1620	SAFEDELETE(apvbCzyMinimalnaOdleglosc[0]);
	1621	SAFEDELETE(apvbCzyMinimalnaOdleglosc[1]);
	1622	}
	1623	} // while (! bCzyKoniecGrupy)
[349]	1624
[647]	1625	// PO DOPASOWANIU WSZYSTKIEGO Z GRUP (CO NAJMNIEJ JEDNEJ GRUPY W CALOSCI)
[349]	1626
[606]	1627	// gdy rozmiar ktorejkolwiek z grup dopasowania spadl do zera
	1628	// to musimy przesunac KoniecPierwszejGrupy (wszystkie dopasowane)
	1629	for (i = 0; i < 2; i++)
	1630	{
	1631	// grupy nie moga miec mniejszego niz 0 rozmiaru
	1632	assert(aiRozmiarGrupy[i] >= 0);
	1633	if (aiRozmiarGrupy[i] == 0)
	1634	aiKoniecPierwszejGrupy[i] = aiKoniecGrupyDopasowania[i];
	1635	}
[349]	1636
[606]	1637	// sprawdzenie warunku konca dopasowywania - gdy nie
	1638	// ma juz w jakims organizmie co dopasowywac
	1639	if (aiKoniecPierwszejGrupy[0] >= m_aiPartCount[0] \|\|
	1640	aiKoniecPierwszejGrupy[1] >= m_aiPartCount[1])
	1641	{
	1642	iCzyDopasowane = 1;
	1643	}
	1644	} // koniec WHILE - petli dopasowania
[349]	1645	}
	1646
	1647	/** Matches Parts in both organisms so that computation of similarity is possible.
[606]	1648	New algorithm (assures symmetry of the similarity measure) with geometry
	1649	taken into consideration.
	1650	Assumptions:
	1651	- Models (m_Mod) are created and available.
	1652	- Matching (m_pMatching) is created, but empty
	1653	Exit conditions:
	1654	- Matching (m_pMatching) is full
	1655	@return 1 if success, 0 otherwise
	1656	*/
[349]	1657	int ModelSimil::MatchPartsGeometry()
	1658	{
[606]	1659	// zaloz, ze sa modele i sa poprawne
	1660	assert((m_Mod[0] != NULL) && (m_Mod[1] != NULL));
	1661	assert(m_Mod[0]->isValid() && m_Mod[1]->isValid());
[349]	1662
[606]	1663	if (!CreatePartInfoTables())
	1664	return 0;
	1665	if (!CountPartDegrees())
	1666	return 0;
	1667	if (!GetPartPositions())
	1668	{
	1669	return 0;
	1670	}
	1671	if (!CountPartNeurons())
	1672	return 0;
[349]	1673
	1674
[606]	1675	if (m_adFactors[3] > 0)
	1676	{
	1677	if (!ComputePartsPositionsBySVD())
	1678	{
	1679	return 0;
	1680	}
	1681	}
[349]	1682
[606]	1683	DB(printf("Przed sortowaniem:\n");)
	1684	DB(_PrintDegrees(0);)
	1685	DB(_PrintDegrees(1);)
[349]	1686
[606]	1687	if (!SortPartInfoTables())
	1688	return 0;
[349]	1689
[606]	1690	DB(printf("Po sortowaniu:\n");)
	1691	DB(_PrintDegrees(0);)
	1692	DB(_PrintDegrees(1);)
[349]	1693
[606]	1694	if (m_adFactors[3] > 0)
	1695	{
[869]	1696	// tutaj zacznij pętlę po przekształceniach geometrycznych
	1697	const int NO_OF_TRANSFORM = 8; // liczba transformacji geometrycznych (na razie tylko ID i O_YZ)
	1698	// tablice transformacji współrzędnych; nie są to dokładnie tablice tranformacji, ale raczej tablice PRZEJŚĆ
	1699	// pomiędzy transformacjami;
	1700	// wartości orginalne transformacji dOrig uzyskuje się przez:
	1701	// for ( iTrans = 0; iTrans <= TRANS_INDEX; iTrans++ ) dOrig *= dMul[ iTrans ];
	1702	//const char *szTransformNames[NO_OF_TRANSFORM] = { "ID", "S_yz", "S_xz", "S_xy", "R180_z", "R180_y", "R180_z", "S_(0,0,0)" };
	1703	const int dMulX[NO_OF_TRANSFORM] = { 1, -1, -1, 1, -1, 1, -1, -1 };
	1704	const int dMulY[NO_OF_TRANSFORM] = { 1, 1, -1, -1, -1, -1, -1, 1 };
	1705	const int dMulZ[NO_OF_TRANSFORM] = { 1, 1, 1, -1, -1, -1, 1, 1 };
[349]	1706
[361]	1707	#ifdef max
[606]	1708	#undef max //this macro would conflict with line below
[361]	1709	#endif
[869]	1710	double dMinSimValue = std::numeric_limits<double>::max(); // minimum value of similarity
	1711	int iMinSimTransform = -1; // index of the transformation with the minimum similarity
	1712	SimilMatching *pMinSimMatching = NULL; // matching with the minimum value of similarity
[349]	1713
[869]	1714	// remember the original positions of model 0 as computed by SVD in order to restore them later, after
	1715	// all transformations have been computed
	1716	Pt3D *StoredPositions = NULL; // array of positions of Parts, for one (0th) model
	1717	// create the stored positions
	1718	StoredPositions = new Pt3D[m_Mod[0]->getPartCount()];
	1719	assert(StoredPositions != NULL);
	1720	// copy the original positions of model 0 (store them)
	1721	int iPart; // a counter of Parts
[606]	1722	for (iPart = 0; iPart < m_Mod[0]->getPartCount(); iPart++)
	1723	{
[869]	1724	StoredPositions[iPart].x = m_aPositions[0][iPart].x;
	1725	StoredPositions[iPart].y = m_aPositions[0][iPart].y;
	1726	StoredPositions[iPart].z = m_aPositions[0][iPart].z;
[606]	1727	}
[869]	1728	// now the original positions of model 0 are stored
[349]	1729
	1730
[869]	1731	int iTransform; // a counter of geometric transformations
	1732	for (iTransform = 0; iTransform < NO_OF_TRANSFORM; iTransform++)
[606]	1733	{
[869]	1734	// for each geometric transformation to be done
	1735	// entry conditions:
	1736	// - models (m_Mod) exist and are available
	1737	// - matching (m_pMatching) exists and is empty
	1738	// - all properties are created and available (m_aDegrees and m_aPositions)
[349]	1739
[869]	1740	// recompute geometric properties according to the transformation iTransform
	1741	// but only for model 0
	1742	for (iPart = 0; iPart < m_Mod[0]->getPartCount(); iPart++)
	1743	{
	1744	// for each iPart, a part of the model iMod
	1745	m_aPositions[0][iPart].x *= dMulX[iTransform];
	1746	m_aPositions[0][iPart].y *= dMulY[iTransform];
	1747	m_aPositions[0][iPart].z *= dMulZ[iTransform];
	1748	}
	1749	// now the positions are recomputed
	1750	ComputeMatching();
[349]	1751
[869]	1752	// teraz m_pMatching istnieje i jest pełne
	1753	assert(m_pMatching != NULL);
	1754	assert(m_pMatching->IsFull() == true);
[349]	1755
[869]	1756	// wykorzystaj to dopasowanie i geometrię do obliczenia tymczasowej wartości miary
	1757	int iTempRes = CountPartsDistance();
	1758	// załóż sukces
	1759	assert(iTempRes == 1);
	1760	double dCurrentSim = m_adFactors[0] * double(m_iDV) +
	1761	m_adFactors[1] * double(m_iDD) +
	1762	m_adFactors[2] * double(m_iDN) +
	1763	m_adFactors[3] * double(m_dDG);
	1764	// załóż poprawną wartość podobieństwa
	1765	assert(dCurrentSim >= 0.0);
[349]	1766
[869]	1767	// porównaj wartość obliczoną z dotychczasowym minimum
	1768	if (dCurrentSim < dMinSimValue)
	1769	{
	1770	// jeśli uzyskano mniejszą wartość dopasowania,
	1771	// to zapamiętaj to przekształcenie geometryczne i dopasowanie
	1772	dMinSimValue = dCurrentSim;
	1773	iMinSimTransform = iTransform;
	1774	SAFEDELETE(pMinSimMatching);
	1775	pMinSimMatching = new SimilMatching(*m_pMatching);
	1776	assert(pMinSimMatching != NULL);
	1777	}
	1778
	1779	// teraz już można usunąć stare dopasowanie (dla potrzeb następnego przebiegu pętli)
	1780	m_pMatching->Empty();
	1781	} // for ( iTransform )
	1782
	1783	// po pętli przywróć najlepsze dopasowanie
	1784	delete m_pMatching;
	1785	m_pMatching = pMinSimMatching;
	1786
	1787	DB(printf("Matching has been chosen!\n");)
	1788	// print the name of the chosen transformation:
	1789	// printf("Chosen transformation: %s\n", szTransformNames[ iMinSimTransform ] );
	1790
	1791	// i przywróć jednocześnie pozycje Parts modelu 0 dla tego dopasowania
	1792	// - najpierw przywroc Parts pozycje orginalne, po SVD
	1793	for (iPart = 0; iPart < m_Mod[0]->getPartCount(); iPart++)
	1794	{
	1795	m_aPositions[0][iPart].x = StoredPositions[iPart].x;
	1796	m_aPositions[0][iPart].y = StoredPositions[iPart].y;
	1797	m_aPositions[0][iPart].z = StoredPositions[iPart].z;
	1798	}
	1799	// - usun teraz zapamietane pozycje Parts
	1800	delete[] StoredPositions;
	1801	// - a teraz oblicz na nowo wspolrzedne wlasciwego przeksztalcenia dla model 0
	1802	for (iTransform = 0; iTransform <= iMinSimTransform; iTransform++)
[606]	1803	{
[869]	1804	// for each transformation before and INCLUDING iMinTransform
	1805	// do the transformation (only model 0)
	1806	for (iPart = 0; iPart < m_Mod[0]->getPartCount(); iPart++)
	1807	{
	1808	m_aPositions[0][iPart].x *= dMulX[iTransform];
	1809	m_aPositions[0][iPart].y *= dMulY[iTransform];
	1810	m_aPositions[0][iPart].z *= dMulZ[iTransform];
	1811	}
[606]	1812	}
[349]	1813
[606]	1814	}
	1815	else
	1816	{
	1817	ComputeMatching();
	1818	}
	1819	// teraz dopasowanie musi byc pelne - co najmniej w jednym organizmie musza byc
	1820	// wszystkie elementy dopasowane
	1821	assert(m_pMatching->IsFull() == true);
[349]	1822
[606]	1823	// _PrintDegrees(0);
	1824	// _PrintDegrees(1);
[349]	1825
[606]	1826	DB(_PrintPartsMatching();)
[349]	1827
[606]	1828	return 1;
[349]	1829	}
	1830
	1831	void ModelSimil::_PrintSeamnessTable(std::vector<int> *pTable, int iCount)
	1832	{
[606]	1833	int i;
	1834	printf(" ");
	1835	for (i = 0; i < iCount; i++)
	1836	printf("%3i ", i);
	1837	printf("\n ");
	1838	for (i = 0; i < iCount; i++)
	1839	{
[349]	1840
[606]	1841	printf("%3i ", pTable->operator[](i));
	1842	}
	1843	printf("\n");
[349]	1844	}
	1845
	1846	/** Computes elements of similarity (m_iDD, m_iDN, m_dDG) based on underlying matching.
[606]	1847	Assumptions:
	1848	- Matching (m_pMatching) exists and is full.
	1849	- Internal arrays m_aDegrees and m_aPositions exist and are properly filled in
	1850	Exit conditions:
	1851	- Elements of similarity are computed (m)iDD, m_iDN, m_dDG).
	1852	@return 1 if success, otherwise 0.
	1853	*/
[349]	1854	int ModelSimil::CountPartsDistance()
	1855	{
[606]	1856	int i, temp;
[349]	1857
[606]	1858	// assume existence of m_pMatching
	1859	assert(m_pMatching != NULL);
	1860	// musi byc pelne!
	1861	assert(m_pMatching->IsFull() == true);
[349]	1862
[606]	1863	// roznica w stopniach
	1864	m_iDD = 0;
	1865	// roznica w liczbie neuronów
	1866	m_iDN = 0;
	1867	// roznica w odleglosci dopasowanych Parts
	1868	m_dDG = 0.0;
[349]	1869
[606]	1870	int iOrgPart, iOrgMatchedPart; // orginalny indeks Part i jego dopasowanego Part
	1871	int iMatchedPart; // indeks (wg sortowania) dopasowanego Part
[349]	1872
[606]	1873	// wykorzystanie dopasowania do zliczenia m_iDD - roznicy w stopniach
	1874	// i m_iDN - roznicy w liczbie neuronow
	1875	// petla w wiekszej tablicy
	1876	for (i = 0; i < m_aiPartCount[1 - m_iSmaller]; i++)
	1877	{
	1878	// dla kazdego elementu [i] z wiekszego organizmu
	1879	// pobierz jego orginalny indeks w organizmie z tablicy TDN
	1880	iOrgPart = m_aDegrees[1 - m_iSmaller][i][0];
	1881	if (!(m_pMatching->IsMatched(1 - m_iSmaller, i)))
	1882	{
	1883	// gdy nie zostal dopasowany
	1884	// dodaj jego stopien do DD
	1885	m_iDD += m_aDegrees[1 - m_iSmaller][i][1];
	1886	// dodaj liczbe neuronow do DN
	1887	m_iDN += m_aDegrees[1 - m_iSmaller][i][3];
	1888	// i oblicz odleglosc tego Part od srodka organizmu (0,0,0)
	1889	// (uzyj orginalnego indeksu)
	1890	//no need to compute distane when m_dDG weight is 0
	1891	m_dDG += m_adFactors[3] == 0 ? 0 : m_aPositions[1 - m_iSmaller][iOrgPart].length();
	1892	}
	1893	else
	1894	{
	1895	// gdy byl dopasowany
	1896	// pobierz indeks (po sortowaniu) tego dopasowanego Part
	1897	iMatchedPart = m_pMatching->GetMatchedIndex(1 - m_iSmaller, i);
	1898	// pobierz indeks orginalny tego dopasowanego Part
	1899	iOrgMatchedPart = m_aDegrees[m_iSmaller][iMatchedPart][0];
	1900	// dodaj ABS roznicy stopni do DD
	1901	temp = m_aDegrees[1 - m_iSmaller][i][1] -
	1902	m_aDegrees[m_iSmaller][iMatchedPart][1];
	1903	m_iDD += abs(temp);
	1904	// dodaj ABS roznicy neuronow do DN
	1905	temp = m_aDegrees[1 - m_iSmaller][i][3] -
	1906	m_aDegrees[m_iSmaller][iMatchedPart][3];
	1907	m_iDN += abs(temp);
	1908	// pobierz polozenie dopasowanego Part
	1909	Pt3D MatchedPartPos(m_aPositions[m_iSmaller][iOrgMatchedPart]);
	1910	// dodaj euklidesowa odleglosc Parts do sumy odleglosci
	1911	//no need to compute distane when m_dDG weight is 0
	1912	m_dDG += m_adFactors[3] == 0 ? 0 : m_aPositions[1 - m_iSmaller][iOrgPart].distanceTo(MatchedPartPos);
	1913	}
	1914	}
[349]	1915
[606]	1916	// obliczenie i dodanie różnicy w liczbie neuronów OnJoint...
	1917	temp = m_aOnJoint[0][3] - m_aOnJoint[1][3];
	1918	m_iDN += abs(temp);
	1919	DB(printf("OnJoint DN: %i\n", abs(temp));)
	1920	// ... i Anywhere
	1921	temp = m_aAnywhere[0][3] - m_aAnywhere[1][3];
	1922	m_iDN += abs(temp);
	1923	DB(printf("Anywhere DN: %i\n", abs(temp));)
[349]	1924
[606]	1925	return 1;
[349]	1926	}
	1927
	1928	/** Computes new positions of Parts of both of organisms stored in the object.
[606]	1929	Assumptions:
	1930	- models (m_Mod) are created and valid
	1931	- positions (m_aPositions) are created and filled with original positions of Parts
	1932	- the sorting algorithm was not yet run on the array m_aDegrees
	1933	@return true if successful; false otherwise
	1934	*/
[349]	1935	bool ModelSimil::ComputePartsPositionsBySVD()
	1936	{
[606]	1937	bool bResult = true; // the result; assume a success
[349]	1938
[606]	1939	// check assumptions
	1940	// the models
	1941	assert(m_Mod[0] != NULL && m_Mod[0]->isValid());
	1942	assert(m_Mod[1] != NULL && m_Mod[1]->isValid());
	1943	// the position arrays
	1944	assert(m_aPositions[0] != NULL);
	1945	assert(m_aPositions[1] != NULL);
[349]	1946
[606]	1947	int iMod; // a counter of models
	1948	// use SVD to obtain different point of view on organisms
	1949	// and store the new positions (currently the original ones are still stored)
	1950	for (iMod = 0; iMod < 2; iMod++)
	1951	{
	1952	// prepare the vector of errors of approximation for the SVD
	1953	std::vector<double> vEigenvalues;
	1954	int nSize = m_Mod[iMod]->getPartCount();
[349]	1955
[869]	1956	double pDistances = new double[nSize nSize];
[349]	1957
[606]	1958	for (int i = 0; i < nSize; i++)
	1959	{
	1960	pDistances[i] = 0;
	1961	}
[349]	1962
[606]	1963	Model *pModel = m_Mod[iMod]; // use the model of the iMod (current) organism
	1964	int iP1, iP2; // indices of Parts in the model
	1965	Part P1, P2; // pointers to Parts
	1966	Pt3D P1Pos, P2Pos; // positions of parts
	1967	double dDistance; // the distance between Parts
[869]	1968
[817]	1969	double *weights = new double[nSize];
	1970	for (int i = 0; i < nSize; i++)
	1971	{
[869]	1972	if (wMDS == 1)
[817]	1973	weights[i] = 0;
	1974	else
	1975	weights[i] = 1;
	1976	}
[869]	1977
	1978	if (wMDS == 1)
[817]	1979	for (int i = 0; i < pModel->getJointCount(); i++)
	1980	{
	1981	weights[pModel->getJoint(i)->p1_refno]++;
[869]	1982	weights[pModel->getJoint(i)->p2_refno]++;
[817]	1983	}
[869]	1984
[606]	1985	for (iP1 = 0; iP1 < pModel->getPartCount(); iP1++)
	1986	{
	1987	// for each iP1, a Part index in the model of organism iMod
	1988	P1 = pModel->getPart(iP1);
	1989	// get the position of the Part
	1990	P1Pos = P1->p;
[605]	1991	if (fixedZaxis == 1)
[606]	1992	{
	1993	P1Pos.z = 0; //fixed vertical axis, so pretend all points are on the xy plane
	1994	}
	1995	for (iP2 = 0; iP2 < pModel->getPartCount(); iP2++)
	1996	{
	1997	// for each (iP1, iP2), a pair of Parts index in the model
	1998	P2 = pModel->getPart(iP2);
	1999	// get the position of the Part
	2000	P2Pos = P2->p;
[605]	2001	if (fixedZaxis == 1)
[606]	2002	{
	2003	P2Pos.z = 0; //fixed vertical axis, so pretend all points are on the xy plane
	2004	}
	2005	// compute the geometric (Euclidean) distance between the Parts
	2006	dDistance = P1Pos.distanceTo(P2Pos);
	2007	// store the distance
	2008	pDistances[iP1 * nSize + iP2] = dDistance;
	2009	} // for (iP2)
	2010	} // for (iP1)
[349]	2011
[817]	2012	MatrixTools::weightedMDS(vEigenvalues, nSize, pDistances, m_aPositions[iMod], weights);
[605]	2013	if (fixedZaxis == 1) //restore the original vertical coordinate of each Part
[601]	2014	{
[606]	2015	for (int part = 0; part < pModel->getPartCount(); part++)
	2016	{
	2017	m_aPositions[iMod][part].z = pModel->getPart(part)->p.z;
	2018	}
[601]	2019	}
[606]	2020
[817]	2021	delete[] pDistances;
	2022	delete[] weights;
[601]	2023	}
[349]	2024
[606]	2025	return bResult;
[349]	2026	}
	2027
[869]	2028	/** Evaluates distance between two given genotypes. The distance depends strongly
	2029	on weights set and the matching algorithm used.
	2030	@param G0 Pointer to the first of compared genotypes
	2031	@param G1 Pointer to the second of compared genotypes.
	2032	@return Distance between two genotypes.
	2033	@sa m_adFactors, matching_method
	2034	*/
	2035	double ModelSimil::EvaluateDistance(const Geno G0, const Geno G1)
	2036	{
	2037	return matching_method == 0 ? EvaluateDistanceHungarian(G0, G1) : EvaluateDistanceGreedy(G0, G1);
	2038	}
	2039
[349]	2040	void ModelSimil::p_evaldistance(ExtValue args, ExtValue ret)
	2041	{
[606]	2042	Geno *g1 = GenoObj::fromObject(args[1]);
	2043	Geno *g2 = GenoObj::fromObject(args[0]);
	2044	if ((!g1) \|\| (!g2))
	2045	ret->setEmpty();
	2046	else
	2047	ret->setDouble(EvaluateDistance(g1, g2));
[356]	2048	}
[869]	2049
	2050	void ModelSimil::FillPartsDistances(double*& dist, int bigger, int smaller, bool geo)
	2051	{
	2052	for (int i = 0; i < bigger; i++)
	2053	{
	2054	for (int j = 0; j < bigger; j++)
	2055	{
	2056	// assign penalty for unassignment for vertex from bigger model
	2057	if (j >= smaller)
	2058	{
	2059	if (geo)
	2060	dist[ibigger + j] += m_adFactors[3] m_aPositions[1 - m_iSmaller][i].length();
	2061	else
	2062	dist[ibigger + j] = m_adFactors[1] m_aDegrees[1 - m_iSmaller][i][DEGREE] + m_adFactors[2] * m_aDegrees[1 - m_iSmaller][i][NEURONS];
	2063	}
	2064	// compute distance between parts
	2065	else
	2066	{
	2067	if (geo)
	2068	dist[ibigger + j] += m_adFactors[3] m_aPositions[1 - m_iSmaller][i].distanceTo(m_aPositions[m_iSmaller][j]);
	2069	else
	2070	dist[ibigger + j] = m_adFactors[1] abs(m_aDegrees[1 - m_iSmaller][i][DEGREE] - m_aDegrees[m_iSmaller][j][DEGREE])
	2071	+ m_adFactors[2] * abs(m_aDegrees[1 - m_iSmaller][i][NEURONS] - m_aDegrees[m_iSmaller][j][NEURONS]);
	2072	}
	2073
	2074	}
	2075	}
	2076	}
	2077
	2078	double ModelSimil::EvaluateDistanceHungarian(const Geno G0, const Geno G1)
	2079	{
[877]	2080	double dResult = 0.0;
[869]	2081
	2082	m_Gen[0] = G0;
	2083	m_Gen[1] = G1;
	2084
	2085	// create models of objects to compare
[877]	2086	m_Mod[0] = newModel(m_Gen[0]);
	2087	m_Mod[1] = newModel(m_Gen[1]);
[869]	2088
[877]	2089	if (m_Mod[0] == NULL \|\| m_Mod[1] == NULL)
[872]	2090	return 0.0;
[869]	2091
	2092	//Get information about vertex degrees, neurons and 3D location
	2093	if (!CreatePartInfoTables())
[877]	2094	return 0.0;
[869]	2095	if (!CountPartDegrees())
[877]	2096	return 0.0;
[869]	2097	if (!GetPartPositions())
[877]	2098	return 0.0;
[869]	2099	if (!CountPartNeurons())
[877]	2100	return 0.0;
[869]	2101
	2102	m_iSmaller = m_Mod[0]->getPartCount() <= m_Mod[1]->getPartCount() ? 0 : 1;
	2103	int nSmaller = m_Mod[m_iSmaller]->getPartCount();
	2104	int nBigger = m_Mod[1 - m_iSmaller]->getPartCount();
	2105
	2106	double* partsDistances = new double[nBigger*nBigger]();
	2107	FillPartsDistances(partsDistances, nBigger, nSmaller, false);
	2108	int *assignment = new int[nBigger]();
	2109
	2110	HungarianAlgorithm hungarian;
	2111
	2112	if (m_adFactors[3] > 0)
	2113	{
	2114	if (!ComputePartsPositionsBySVD())
	2115	{
[877]	2116	return 0.0;
[869]	2117	}
	2118
	2119	// tutaj zacznij pętlę po przekształceniach geometrycznych
	2120	const int NO_OF_TRANSFORM = 8; // liczba transformacji geometrycznych (na razie tylko ID i O_YZ)
	2121	// tablice transformacji współrzędnych; nie są to dokładnie tablice tranformacji, ale raczej tablice PRZEJŚĆ
	2122	// pomiędzy transformacjami;
	2123	const int dMulX[NO_OF_TRANSFORM] = { 1, -1, -1, 1, -1, 1, -1, -1 };
	2124	const int dMulY[NO_OF_TRANSFORM] = { 1, 1, -1, -1, -1, -1, -1, 1 };
	2125	const int dMulZ[NO_OF_TRANSFORM] = { 1, 1, 1, -1, -1, -1, 1, 1 };
	2126
	2127	std::vector<int> minAssignment(nBigger);
	2128	#ifdef max
	2129	#undef max //this macro would conflict with line below
	2130	#endif
	2131	double dMinSimValue = std::numeric_limits<double>::max(); // minimum value of similarity
	2132
	2133	int iTransform; // a counter of geometric transformations
	2134	for (iTransform = 0; iTransform < NO_OF_TRANSFORM; iTransform++)
	2135	{
	2136	// for each geometric transformation to be done
	2137	// entry conditions:
	2138	// - models (m_Mod) exist and are available
	2139	// - all properties are created and available (m_aDegrees and m_aPositions)
	2140	double* tmpPartsDistances = new double[nBigger*nBigger]();
	2141	std::copy(partsDistances, partsDistances + nBigger * nBigger, tmpPartsDistances);
	2142	// recompute geometric properties according to the transformation iTransform
	2143	// but only for model 0
	2144	for (int iPart = 0; iPart < m_Mod[m_iSmaller]->getPartCount(); iPart++)
	2145	{
	2146	// for each iPart, a part of the model iMod
	2147	m_aPositions[m_iSmaller][iPart].x *= dMulX[iTransform];
	2148	m_aPositions[m_iSmaller][iPart].y *= dMulY[iTransform];
	2149	m_aPositions[m_iSmaller][iPart].z *= dMulZ[iTransform];
	2150	}
	2151	// now the positions are recomputed
	2152
	2153	FillPartsDistances(tmpPartsDistances, nBigger, nSmaller, true);
	2154	std::fill_n(assignment, nBigger, 0);
	2155	double dCurrentSim = hungarian.Solve(tmpPartsDistances, assignment, nBigger, nBigger);
	2156
	2157	delete[] tmpPartsDistances;
	2158	// załóż poprawną wartość podobieństwa
	2159	assert(dCurrentSim >= 0.0);
	2160
	2161	// porównaj wartość obliczoną z dotychczasowym minimum
	2162	if (dCurrentSim < dMinSimValue)
	2163	{
	2164	dMinSimValue = dCurrentSim;
[872]	2165	if (saveMatching)
[869]	2166	{
	2167	minAssignment.clear();
	2168	minAssignment.insert(minAssignment.begin(), assignment, assignment + nBigger);
	2169	}
	2170	}
	2171	}
	2172
	2173	dResult = dMinSimValue;
[872]	2174	if (saveMatching)
[869]	2175	std::copy(minAssignment.begin(), minAssignment.end(), assignment);
	2176	}
	2177	else
	2178	{
	2179	dResult = hungarian.Solve(partsDistances, assignment, nBigger, nBigger);
	2180	}
	2181
	2182	//add difference in anywhere and onJoint neurons
	2183	dResult += m_adFactors[2] * (abs(m_aOnJoint[0][3] - m_aOnJoint[1][3]) + abs(m_aAnywhere[0][3] - m_aAnywhere[1][3]));
	2184	//add difference in part numbers
	2185	dResult += (nBigger - nSmaller) * m_adFactors[0];
	2186
	2187	// delete degree arrays created in CreatePartInfoTables
	2188	SAFEDELETEARRAY(m_aDegrees[0]);
	2189	SAFEDELETEARRAY(m_aDegrees[1]);
	2190
	2191	// and position arrays
	2192	SAFEDELETEARRAY(m_aPositions[0]);
	2193	SAFEDELETEARRAY(m_aPositions[1]);
	2194
	2195	// delete created models
	2196	SAFEDELETE(m_Mod[0]);
	2197	SAFEDELETE(m_Mod[1]);
	2198
	2199	delete[] assignment;
	2200	delete[] partsDistances;
	2201
	2202	return dResult;
	2203	}

Note: See TracBrowser for help on using the repository browser.

Download in other formats: