source: cpp/frams/model/similarity/measure-greedy.cpp @ 1304

// This file is a part of Framsticks SDK.  http://www.framsticks.com/
// Copyright (C) 1999-2021  Maciej Komosinski and Szymon Ulatowski.
// See LICENSE.txt for details.

#include "measure-greedy.h"
#include <cstdlib> //std::qsort()
#include <assert.h>

#define DB(x)  //define as x if you want to print debug information

const int SimilMeasureGreedy::iNOFactors = 4;
int fuzzDepth = 0; //TODO make local, but not crucial because currently "fuzzy vertex degree" is not activated by default

#define FIELDSTRUCT SimilMeasureGreedy

static ParamEntry simil_greedy_paramtab[] = {
                { "Creature: Similarity: Graph greedy", 1, 7, "SimilMeasureGreedy", "Evaluates morphological dissimilarity using the greedy measure. 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://www.framsticks.com/bib/Komosinski-and-Mensfelt-2019", },
                { "simil_greedy_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.", },
                { "simil_greedy_partdeg", 0, 0, "Weight of parts' degree", "f 0 100 1", FIELD(m_adFactors[1]), "", },
                { "simil_greedy_neuro", 0, 0, "Weight of neurons count", "f 0 100 0.1", FIELD(m_adFactors[2]), "", },
                { "simil_greedy_partgeom", 0, 0, "Weight of parts' geometric distances", "f 0 100 0", FIELD(m_adFactors[3]), "", },
                { "simil_greedy_fixedZaxis", 0, 0, "Fix 'z' (vertical) axis?", "d 0 1 0", FIELD(fixedZaxis), "", },
                { "simil_greedy_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.", },
                { "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.", },
                { 0, },
};

#undef FIELDSTRUCT

SimilMeasureGreedy::SimilMeasureGreedy(): localpar(simil_greedy_paramtab, this), m_iDV(0), m_iDD(0), m_iDN(0), m_dDG(0.0)
{       
        localpar.setDefault();
        
        for (int i = 0; i < 2; i++)
        {
                m_aDegrees[0] = nullptr;
                m_fuzzyNeighb[0] = nullptr;
                m_Neighbours[0] = nullptr;
        }
        
        m_pMatching = nullptr;
        tempMatching = nullptr;
        
        //Determines whether "fuzzy vertex degree" should be used.
        //In preliminary experiments in 2017, "fuzzy vertex degree" turned out to be not beneficial.
        isFuzzy = false;
        fuzzyDepth = 10;
        save_matching = false;
}

double SimilMeasureGreedy::distanceForTransformation()
{
        // now the positions are recomputed
        computeMatching();

        // teraz m_pMatching istnieje i jest pełne
        assert(m_pMatching != NULL);
        assert(m_pMatching->isFull() == true);
        // wykorzystaj to dopasowanie i geometrię do obliczenia tymczasowej wartości miary
        
        int iTempRes = countPartsDistance();
        // załóż sukces
        assert(iTempRes == 1);
        double dCurrentSim = m_adFactors[0] * double(m_iDV) +
                m_adFactors[1] * double(m_iDD) +
                m_adFactors[2] * double(m_iDN) +
                m_adFactors[3] * double(m_dDG);
        // załóż poprawną wartość podobieństwa
        assert(dCurrentSim >= 0.0);
        SAFEDELETE(tempMatching);
        return dCurrentSim;
}


double SimilMeasureGreedy::distanceWithoutAlignment()
{
        return distanceForTransformation();
}

void SimilMeasureGreedy::beforeTransformation()
{
        if (m_pMatching != NULL)
                if (!m_pMatching->isEmpty())
                        m_pMatching->empty();
}

/** Computes elements of similarity (m_iDD, m_iDN, m_dDG) based on underlying matching.
        Assumptions:
        - Matching (m_pMatching) exists and is full.
        - Internal arrays m_aDegrees and coordinates exist and are properly filled in
        Exit conditions:
        - Elements of similarity are computed (m)iDD, m_iDN, m_dDG).
        @return 1 if success, otherwise 0.
        */
int SimilMeasureGreedy::countPartsDistance()
{
        int i, temp;

        // assume existence of m_pMatching
        assert(m_pMatching != NULL);
        // musi byc pelne!
        assert(m_pMatching->isFull() == true);

        // roznica w stopniach
        m_iDD = 0;
        // roznica w liczbie neuronów
        m_iDN = 0;
        // roznica w odleglosci dopasowanych Parts
        m_dDG = 0.0;

        int iOrgPart, iOrgMatchedPart; // orginalny indeks Part i jego dopasowanego Part
        int iMatchedPart; // indeks (wg sortowania) dopasowanego Part

        // wykorzystanie dopasowania do zliczenia m_iDD - roznicy w stopniach
        // i m_iDN - roznicy w liczbie neuronow
        // petla w wiekszej tablicy
        for (i = 0; i < m_aiPartCount[1 - m_iSmaller]; i++)
        {
                // dla kazdego elementu [i] z wiekszego organizmu
                // pobierz jego orginalny indeks w organizmie z tablicy TDN
                iOrgPart = m_aDegrees[1 - m_iSmaller][i][0];
                if (!(m_pMatching->isMatched(1 - m_iSmaller, i)))
                {
                        // gdy nie zostal dopasowany
                        // dodaj jego stopien do DD
                        m_iDD += m_aDegrees[1 - m_iSmaller][i][1];
                        // dodaj liczbe neuronow do DN
                        m_iDN += m_aDegrees[1 - m_iSmaller][i][3];
                        // i oblicz odleglosc tego Part od srodka organizmu (0,0,0)
                        // (uzyj orginalnego indeksu)
                        //no need to compute distane when m_dDG weight is 0
                        m_dDG += m_adFactors[3] == 0 ? 0 : coordinates[1 - m_iSmaller][iOrgPart].length();
                }
                else
                {
                        // gdy byl dopasowany
                        // pobierz indeks (po sortowaniu) tego dopasowanego Part
                        iMatchedPart = m_pMatching->getMatchedIndex(1 - m_iSmaller, i);
                        // pobierz indeks orginalny tego dopasowanego Part
                        iOrgMatchedPart = m_aDegrees[m_iSmaller][iMatchedPart][0];
                        // dodaj ABS roznicy stopni do DD
                        temp = m_aDegrees[1 - m_iSmaller][i][1] -
                                m_aDegrees[m_iSmaller][iMatchedPart][1];
                        m_iDD += abs(temp);
                        // dodaj ABS roznicy neuronow do DN
                        temp = m_aDegrees[1 - m_iSmaller][i][3] -
                                m_aDegrees[m_iSmaller][iMatchedPart][3];
                        m_iDN += abs(temp);
                        // pobierz polozenie dopasowanego Part
                        if (m_adFactors[3] > 0) //no need to compute distance when m_dDG weight is 0
                        {
                            Pt3D MatchedPartPos(coordinates[m_iSmaller][iOrgMatchedPart]);
                            // dodaj euklidesowa odleglosc Parts do sumy odleglosci
                            m_dDG += coordinates[1 - m_iSmaller][iOrgPart].distanceTo(MatchedPartPos);
                        }
                }
        }

        // obliczenie i dodanie różnicy w liczbie neuronów OnJoint...
        temp = m_aOnJoint[0][3] - m_aOnJoint[1][3];
        m_iDN += abs(temp);
        DB(printf("OnJoint DN: %i\n", abs(temp));)
                // ... i Anywhere
                temp = m_aAnywhere[0][3] - m_aAnywhere[1][3];
        m_iDN += abs(temp);
        DB(printf("Anywhere DN: %i\n", abs(temp));)

                return 1;
}

void SimilMeasureGreedy::computeMatching()
{
        // uniwersalne liczniki
        int i, j;

        assert(m_pMatching != NULL);
        assert(m_pMatching->isEmpty() == true);

        // rozpoczynamy etap dopasowywania Parts w organizmach
        // czy dopasowano już wszystkie Parts?
        int iCzyDopasowane = 0;
        // koniec grupy aktualnie dopasowywanej w każdym organizmie
        int aiKoniecGrupyDopasowania[2] = { 0, 0 };
        // koniec grupy już w całości dopasowanej
        // (Pomiedzy tymi dwoma indeksami znajduja sie Parts w tablicy
        // m_aDegrees, ktore moga byc dopasowywane (tam nadal moga
        // byc tez dopasowane - ale nie musi to byc w sposob
        // ciagly)
        int aiKoniecPierwszejGrupy[2] = { 0, 0 };
        // Tablica przechowująca odległości poszczególnych Parts z aktualnych
        // grup dopasowania. Rozmiar - prostokąt o bokach równych liczbie elementów w 
        // dopasowywanych aktualnie grupach. Pierwszy wymiar - pierwszy organizm.
        // Drugi wymiar - drugi organizm (nie zależy to od tego, który jest mniejszy).
        // Wliczane w rozmiar tablicy są nawet już dopasowane elementy - tablice
        // paiCzyDopasowany pamiętają stan dopasowania tych elementów.
        typedef double *TPDouble;
        double **aadOdleglosciParts;
        // dwie tablice okreslajace Parts, ktore moga byc do siebie dopasowywane
        // rozmiary: [0] - aiRozmiarCalychGrup[1]
        //                       [1] - aiRozmiarCalychGrup[0]
        std::vector<bool> *apvbCzyMinimalnaOdleglosc[2];
        // rozmiar aktualnie dopasowywanej grupy w odpowiednim organizmie (tylko elementy
        // jeszcze niedopasowane).
        int aiRozmiarGrupy[2];
        // rozmiar aktualnie dopasowywanych grup w odpowiednim organizmie (włączone są
        // w to również elementy już dopasowane).
        int aiRozmiarCalychGrup[2] = { 0, 0 };

        // DOPASOWYWANIE PARTS
        while (!iCzyDopasowane)
        {
                // znajdz konce obu grup aktualnie dopasowywanych w obu organizmach
                for (i = 0; i < 2; i++)
                {
                        // czyli poszukaj miejsca zmiany stopnia lub konca tablicy
                        for (j = aiKoniecPierwszejGrupy[i] + 1; j < m_aiPartCount[i]; j++)
                        {
                                if (m_aDegrees[i][j][DEG] < m_aDegrees[i][j - 1][DEG])
                                {
                                        break;
                                }
                        }
                        aiKoniecGrupyDopasowania[i] = j;

                        // sprawdz poprawnosc tego indeksu
                        assert((aiKoniecGrupyDopasowania[i] > 0) &&
                                (aiKoniecGrupyDopasowania[i] <= m_aiPartCount[i]));

                        // oblicz rozmiary całych grup - łącznie z dopasowanymi już elementami
                        aiRozmiarCalychGrup[i] = aiKoniecGrupyDopasowania[i] -
                                aiKoniecPierwszejGrupy[i];

                        // sprawdz teraz rozmiar tej grupy w sensie liczby niedopasowanzch
                        // nie moze to byc puste!
                        aiRozmiarGrupy[i] = 0;
                        for (j = aiKoniecPierwszejGrupy[i]; j < aiKoniecGrupyDopasowania[i]; j++)
                        {
                                // od poczatku do konca grupy
                                if (!m_pMatching->isMatched(i, j))
                                {
                                        // jesli niedopasowany, to zwieksz licznik
                                        aiRozmiarGrupy[i]++;
                                }
                        }
                        // grupa nie moze byc pusta!
                        assert(aiRozmiarGrupy[i] > 0);
                }

                // DOPASOWYWANIE PARTS Z GRUP

                // stworzenie tablicy odległości lokalnych
                // stwórz pierwszy wymiar - wg rozmiaru zerowego organizmu
                aadOdleglosciParts = new TPDouble[aiRozmiarCalychGrup[0]];
                assert(aadOdleglosciParts != NULL);
                // stwórz drugi wymiar - wg rozmiaru drugiego organizmu
                for (i = 0; i < aiRozmiarCalychGrup[0]; i++)
                {
                        aadOdleglosciParts[i] = new double[aiRozmiarCalychGrup[1]];
                        assert(aadOdleglosciParts[i] != NULL);
                }

                // stworzenie tablic mozliwosci dopasowania (indykatorow minimalnej odleglosci)
                apvbCzyMinimalnaOdleglosc[0] = new std::vector<bool>(aiRozmiarCalychGrup[1], false);
                apvbCzyMinimalnaOdleglosc[1] = new std::vector<bool>(aiRozmiarCalychGrup[0], false);
                // sprawdz stworzenie tablic
                assert(apvbCzyMinimalnaOdleglosc[0] != NULL);
                assert(apvbCzyMinimalnaOdleglosc[1] != NULL);

                // wypełnienie elementów macierzy (i,j) odległościami pomiędzy
                // odpowiednimi Parts: (i) w organizmie 0 i (j) w organizmie 1.
                // UWAGA! Uwzględniamy tylko te Parts, które nie są jeszcze dopasowane
                // (reszta to byłaby po prostu strata czasu).
                int iDeg, iNeu; // ilościowe określenie różnic stopnia, liczby neuronów i połączeń Parts
                //int iNIt;
                double dGeo; // ilościowe określenie różnic geometrycznych pomiędzy Parts
                // indeksy konkretnych Parts - indeksy sa ZERO-BASED, choć właściwy dostep
                // do informacji o Part wymaga dodania aiKoniecPierwszejGrupy[]
                // tylko aadOdleglosciParts[][] indeksuje sie bezposrednio zawartoscia iIndex[]
                int iIndex[2];
                int iPartIndex[2] = { -1, -1 }; // at [iModel]: original index of a Part for the specified model (iModel)

                // debug - wypisz zakres dopasowywanych indeksow
                DB(printf("Organizm 0: grupa: (%2i, %2i)\n", aiKoniecPierwszejGrupy[0],
                        aiKoniecGrupyDopasowania[0]);)
                        DB(printf("Organizm 1: grupa: (%2i, %2i)\n", aiKoniecPierwszejGrupy[1],
                                aiKoniecGrupyDopasowania[1]);)

                        for (i = 0; i < aiRozmiarCalychGrup[0]; i++)
                        {

                                // iterujemy i - Parts organizmu 0
                                // (indeks podstawowy - aiKoniecPierwszejGrupy[0])

                                if (!(m_pMatching->isMatched(0, aiKoniecPierwszejGrupy[0] + i)))
                                {
                                        // interesuja nas tylko te niedopasowane jeszcze (i)
                                        for (j = 0; j < aiRozmiarCalychGrup[1]; j++)
                                        {

                                                // iterujemy j - Parts organizmu 1
                                                // (indeks podstawowy - aiKoniecPierwszejGrupy[1])

                                                if (!(m_pMatching->isMatched(1, aiKoniecPierwszejGrupy[1] + j)))
                                                {
                                                        // interesuja nas tylko te niedopasowane jeszcze (j)
                                                        // teraz obliczymy lokalne różnice pomiędzy Parts
                                                        iDeg = abs(m_aDegrees[0][aiKoniecPierwszejGrupy[0] + i][1]
                                                                - m_aDegrees[1][aiKoniecPierwszejGrupy[1] + j][1]);

                                                        //iNit currently is not a component of distance measure
                                                        //iNIt = abs(m_aDegrees[0][ aiKoniecPierwszejGrupy[0] + i ][2]
                                                        //- m_aDegrees[1][ aiKoniecPierwszejGrupy[1] + j ][2]);

                                                        iNeu = abs(m_aDegrees[0][aiKoniecPierwszejGrupy[0] + i][3]
                                                                - m_aDegrees[1][aiKoniecPierwszejGrupy[1] + j][3]);

                                                        // obliczenie także lokalnych różnic geometrycznych pomiędzy Parts
                                                        // find original indices of Parts for both of the models
                                                        iPartIndex[0] = m_aDegrees[0][aiKoniecPierwszejGrupy[0] + i][0];
                                                        iPartIndex[1] = m_aDegrees[1][aiKoniecPierwszejGrupy[1] + j][0];
                                                        // now compute the geometrical distance of these Parts (use coordinates
                                                        // which should be computed by SVD)
                                                        
                                                        if (m_adFactors[3] > 0)
                                                        {
                                                                Pt3D Part0Pos(coordinates[0][iPartIndex[0]]);
                                                                Pt3D Part1Pos(coordinates[1][iPartIndex[1]]);
                                                                dGeo = Part0Pos.distanceTo(Part1Pos);
                                                        }
                                                        else
                                                            dGeo = 0; //no need to compute distance when m_dDG weight is 0

                                                        // tutaj prawdopodobnie należy jeszcze dodać sprawdzanie
                                                        // identyczności pozostałych własności (oprócz geometrii)
                                                        // - żeby móc stwierdzić w ogóle identyczność Parts

                                                        // i ostateczna odleglosc indukowana przez te roznice
                                                        // (tutaj nie ma różnicy w liczbie wszystkich wierzchołków)
                                                        aadOdleglosciParts[i][j] = m_adFactors[1] * double(iDeg) +
                                                                m_adFactors[2] * double(iNeu) +
                                                                m_adFactors[3] * dGeo;
                                                        DB(printf("Parts Distance (%2i,%2i) = %.3lf\n", aiKoniecPierwszejGrupy[0] + i,
                                                                aiKoniecPierwszejGrupy[1] + j, aadOdleglosciParts[i][j]);)
                                                                DB(printf("Parts geometrical distance = %.20lf\n", dGeo);)
                                                                DB(printf("Pos0: (%.3lf %.3lf %.3lf)\n", Part0Pos.x, Part0Pos.y, Part0Pos.z);)
                                                                DB(printf("Pos1: (%.3lf %.3lf %.3lf)\n", Part1Pos.x, Part1Pos.y, Part1Pos.z);)
                                                }
                                        }
                                }
                        }

                // tutaj - sprawdzic tylko, czy tablica odleglosci lokalnych jest poprawnie obliczona

                // WYKORZYSTANIE TABLICY ODLEGLOSCI DO BUDOWY DOPASOWANIA

                // trzeba raczej iterować aż do zebrania wszystkich możliwych dopasowań w grupie
                // dlatego wprowadzamy dodatkowa zmienna - czy skonczyla sie juz grupa
                bool bCzyKoniecGrupy = false;
                while (!bCzyKoniecGrupy)
                {
                        for (i = 0; i < 2; i++)
                        {
                                // iterujemy (i) po organizmach
                                // szukamy najpierw jakiegoś niedopasowanego jeszcze Part w organizmach

                                // zakładamy, że nie ma takiego Part
                                iIndex[i] = -1;

                                for (j = 0; j < aiRozmiarCalychGrup[i]; j++)
                                {
                                        // iterujemy (j) - Parts organizmu (i)
                                        // (indeks podstawowy - aiKoniecPierwszejGrupy[0])
                                        if (!(m_pMatching->isMatched(i, aiKoniecPierwszejGrupy[i] + j)))
                                        {
                                                // gdy mamy w tej grupie jakis niedopasowany element, to ustawiamy
                                                // iIndex[i] (chcemy w zasadzie pierwszy taki)
                                                iIndex[i] = j;
                                                break;
                                        }
                                }

                                // sprawdzamy, czy w ogole znaleziono taki Part
                                if (iIndex[i] < 0)
                                {
                                        // gdy nie znaleziono takiego Part - mamy koniec dopasowywania w
                                        // tych grupach
                                        bCzyKoniecGrupy = true;
                                }
                                // sprawdz poprawnosc indeksu niedopasowanego Part - musi byc w aktualnej grupie
                                assert((iIndex[i] >= -1) && (iIndex[i] < aiRozmiarCalychGrup[i]));
                        }


                        // teraz mamy sytuacje:
                        // - mamy w iIndex[0] i iIndex[1] indeksy pierwszych niedopasowanych Part
                        //              w organizmach, albo
                        // - nie ma w ogóle już czego dopasowywać (należy przejść do innej grupy)
                        if (!bCzyKoniecGrupy)
                        {
                                // gdy nie ma jeszcze końca żadnej z grup - możemy dopasowywać
                                // najpierw wyszukujemy w tablicy minimum odległości od tych
                                // wyznaczonych Parts

                                // najpierw wyczyscimy wektory potencjalnych dopasowan
                                // dla organizmu 1 (o rozmiarze grupy z 0)
                                for (i = 0; i < aiRozmiarCalychGrup[0]; i++)
                                {
                                        apvbCzyMinimalnaOdleglosc[1]->operator[](i) = false;
                                }
                                // dla organizmu 0 (o rozmiarze grup z 1)
                                for (i = 0; i < aiRozmiarCalychGrup[1]; i++)
                                {
                                        apvbCzyMinimalnaOdleglosc[0]->operator[](i) = false;
                                }

                                // szukanie minimum dla Part o indeksie iIndex[0] w organizmie 0
                                // wsrod niedopasowanych Parts z organizmu 1
                                // zakładamy, że nie znaleŸliśmy jeszcze minimum
                                double dMinimum = HUGE_VAL;
                                for (i = 0; i < aiRozmiarCalychGrup[1]; i++)
                                {
                                        if (!(m_pMatching->isMatched(1, aiKoniecPierwszejGrupy[1] + i)))
                                        {

                                                // szukamy minimum obliczonej lokalnej odleglosci tylko wsrod
                                                // niedopasowanych jeszcze Parts
                                                if (aadOdleglosciParts[iIndex[0]][i] < dMinimum)
                                                {
                                                        dMinimum = aadOdleglosciParts[iIndex[0]][i];
                                                }
                                                // przy okazji - sprawdz, czy HUGE_VAL jest rzeczywiscie max dla double
                                                assert(aadOdleglosciParts[iIndex[0]][i] < HUGE_VAL);
                                        }
                                }
                                // sprawdz, czy minimum znaleziono - musi takie byc, bo jest cos niedopasowanego
                                assert((dMinimum >= 0.0) && (dMinimum < HUGE_VAL));

                                // teraz zaznaczamy w tablicy te wszystkie Parts, ktore maja
                                // rzeczywiscie te minimalna odleglosc do Part iIndex[0] w organizmie 0
                                for (i = 0; i < aiRozmiarCalychGrup[1]; i++)
                                {
                                        if (!(m_pMatching->isMatched(1, aiKoniecPierwszejGrupy[1] + i)))
                                        {
                                                if (aadOdleglosciParts[iIndex[0]][i] == dMinimum)
                                                {
                                                        // jesli w danym miejscu tablicy odleglosci jest faktyczne
                                                        // minimum odleglosci, to mamy potencjalna pare dla iIndex[0]
                                                        apvbCzyMinimalnaOdleglosc[0]->operator[](i) = true;
                                                }

                                                // sprawdz poprawnosc znalezionego wczesniej minimum
                                                assert(aadOdleglosciParts[iIndex[0]][i] >= dMinimum);
                                        }
                                }

                                // podobnie szukamy minimum dla Part o indeksie iIndex[1] w organizmie 1
                                // wsrod niedopasowanych Parts z organizmu 0
                                dMinimum = HUGE_VAL;
                                for (i = 0; i < aiRozmiarCalychGrup[0]; i++)
                                {
                                        if (!(m_pMatching->isMatched(0, aiKoniecPierwszejGrupy[0] + i)))
                                        {
                                                // szukamy minimum obliczonej lokalnej odleglosci tylko wsrod
                                                // niedopasowanych jeszcze Parts
                                                if (aadOdleglosciParts[i][iIndex[1]] < dMinimum)
                                                {
                                                        dMinimum = aadOdleglosciParts[i][iIndex[1]];
                                                }
                                                // przy okazji - sprawdz, czy HUGE_VAL jest rzeczywiscie max dla double
                                                assert(aadOdleglosciParts[i][iIndex[1]] < HUGE_VAL);
                                        }
                                }
                                // sprawdz, czy minimum znaleziono - musi takie byc, bo jest cos niedopasowanego
                                assert((dMinimum >= 0.0) && (dMinimum < HUGE_VAL));

                                // teraz zaznaczamy w tablicy te wszystkie Parts, ktore maja
                                // rzeczywiscie te minimalne odleglosci do Part iIndex[1] w organizmie 1
                                for (i = 0; i < aiRozmiarCalychGrup[0]; i++)
                                {
                                        if (!(m_pMatching->isMatched(0, aiKoniecPierwszejGrupy[0] + i)))
                                        {
                                                if (aadOdleglosciParts[i][iIndex[1]] == dMinimum)
                                                {
                                                        // jesli w danym miejscu tablicy odleglosci jest faktyczne
                                                        // minimum odleglosci, to mamy potencjalna pare dla iIndex[1]
                                                        apvbCzyMinimalnaOdleglosc[1]->operator[](i) = true;
                                                }

                                                // sprawdz poprawnosc znalezionego wczesniej minimum
                                                assert(aadOdleglosciParts[i][iIndex[1]] >= dMinimum);
                                        }
                                }

                                // teraz mamy juz poszukane potencjalne grupy dopasowania - musimy
                                // zdecydowac, co do czego dopasujemy!
                                // szukamy Part iIndex[0] posrod mozliwych do dopasowania dla Part iIndex[1]
                                // szukamy takze Part iIndex[1] posrod mozliwych do dopasowania dla Part iIndex[0]
                                bool PartZ1NaLiscie0 = apvbCzyMinimalnaOdleglosc[0]->operator[](iIndex[1]);
                                bool PartZ0NaLiscie1 = apvbCzyMinimalnaOdleglosc[1]->operator[](iIndex[0]);

                                if (PartZ1NaLiscie0 && PartZ0NaLiscie1)
                                {
                                        // PRZYPADEK 1. Oba Parts maja sie wzajemnie na listach mozliwych
                                        // dopasowan.
                                        // AKCJA. Dopasowanie wzajemne do siebie.

                                        m_pMatching->match(0, aiKoniecPierwszejGrupy[0] + iIndex[0],
                                                1, aiKoniecPierwszejGrupy[1] + iIndex[1]);

                                        // zmniejsz liczby niedopasowanych elementow w grupach
                                        aiRozmiarGrupy[0]--;
                                        aiRozmiarGrupy[1]--;
                                        // debug - co zostalo dopasowane
                                        DB(printf("Przypadek 1.: dopasowane Parts: (%2i, %2i)\n", aiKoniecPierwszejGrupy[0]
                                                + iIndex[0], aiKoniecPierwszejGrupy[1] + iIndex[1]);)

                                }// PRZYPADEK 1.
                                else
                                        if (PartZ1NaLiscie0 || PartZ0NaLiscie1)
                                        {
                                                // PRZYPADEK 2. Tylko jeden z Parts ma drugiego na swojej liscie
                                                // mozliwych dopasowan
                                                // AKCJA. Dopasowanie jednego jest proste (tego, ktory nie ma 
                                                // na swojej liscie drugiego). Dla tego drugiego nalezy powtorzyc
                                                // duza czesc obliczen (znalezc mu nowa mozliwa pare)

                                                // indeks organizmu, ktorego Part nie ma dopasowywanego Part
                                                // z drugiego organizmu na swojej liscie
                                                int iBezDrugiego;

                                                // okreslenie indeksu organizmu z dopasowywanym Part
                                                if (!PartZ1NaLiscie0)
                                                {
                                                        iBezDrugiego = 0;
                                                }
                                                else
                                                {
                                                        iBezDrugiego = 1;
                                                }
                                                // sprawdz indeks organizmu
                                                assert((iBezDrugiego == 0) || (iBezDrugiego == 1));

                                                int iDopasowywany = -1;
                                                // poszukujemy pierwszego z listy dopasowania
                                                for (i = 0; i < aiRozmiarCalychGrup[1 - iBezDrugiego]; i++)
                                                {
                                                        if (apvbCzyMinimalnaOdleglosc[iBezDrugiego]->operator[](i))
                                                        {
                                                                iDopasowywany = i;
                                                                break;
                                                        }
                                                }
                                                // sprawdz poprawnosc indeksu dopasowywanego (musimy cos znalezc!)
                                                // nieujemny i w odpowiedniej grupie!
                                                assert((iDopasowywany >= 0) &&
                                                        (iDopasowywany < aiRozmiarCalychGrup[1 - iBezDrugiego]));

                                                // znalezlismy 1. Part z listy dopasowania - dopasowujemy!
                                                m_pMatching->match(
                                                        iBezDrugiego,
                                                        aiKoniecPierwszejGrupy[iBezDrugiego] + iIndex[iBezDrugiego],
                                                        1 - iBezDrugiego,
                                                        aiKoniecPierwszejGrupy[1 - iBezDrugiego] + iDopasowywany);
                                                DB(printf("Przypadek 2.1.: dopasowane Parts dopasowanie bez drugiego: (%2i, %2i)\n", aiKoniecPierwszejGrupy[iBezDrugiego] + iIndex[iBezDrugiego],
                                                        aiKoniecPierwszejGrupy[1 - iBezDrugiego] + iDopasowywany);)

                                                        // zmniejsz liczby niedopasowanych elementow w grupach
                                                        aiRozmiarGrupy[0]--;
                                                aiRozmiarGrupy[1]--;

                                                // sprawdz, czy jest szansa dopasowania tego Part z drugiej strony
                                                // (ktora miala mozliwosc dopasowania tego Part z poprzedniego organizmu)
                                                if ((aiRozmiarGrupy[0] > 0) && (aiRozmiarGrupy[1] > 0))
                                                {
                                                        // jesli grupy sie jeszcze nie wyczrpaly
                                                        // to jest mozliwosc dopasowania w organizmie

                                                        int iZDrugim = 1 - iBezDrugiego;
                                                        // sprawdz indeks organizmu
                                                        assert((iZDrugim == 0) || (iZDrugim == 1));

                                                        // bedziemy szukac minimum wsrod niedopasowanych - musimy wykasowac
                                                        // poprzednie obliczenia minimum
                                                        // dla organizmu 1 (o rozmiarze grupy z 0)
                                                        for (i = 0; i < aiRozmiarCalychGrup[0]; i++)
                                                        {
                                                                apvbCzyMinimalnaOdleglosc[1]->operator[](i) = false;
                                                        }
                                                        // dla organizmu 0 (o rozmiarze grup z 1)
                                                        for (i = 0; i < aiRozmiarCalychGrup[1]; i++)
                                                        {
                                                                apvbCzyMinimalnaOdleglosc[0]->operator[](i) = false;
                                                        }

                                                        // szukamy na nowo minimum dla Part o indeksie iIndex[ iZDrugim ] w organizmie iZDrugim
                                                        // wsrod niedopasowanych Parts z organizmu 1 - iZDrugim
                                                        dMinimum = HUGE_VAL;
                                                        for (i = 0; i < aiRozmiarCalychGrup[1 - iZDrugim]; i++)
                                                        {
                                                                if (!(m_pMatching->isMatched(
                                                                        1 - iZDrugim,
                                                                        aiKoniecPierwszejGrupy[1 - iZDrugim] + i)))
                                                                {
                                                                        // szukamy minimum obliczonej lokalnej odleglosci tylko wsrod
                                                                        // niedopasowanych jeszcze Parts
                                                                        if (iZDrugim == 0)
                                                                        {
                                                                                // teraz niestety musimy rozpoznac odpowiedni organizm
                                                                                // zeby moc indeksowac niesymetryczna tablice
                                                                                if (aadOdleglosciParts[iIndex[0]][i] < dMinimum)
                                                                                {
                                                                                        dMinimum = aadOdleglosciParts[iIndex[0]][i];
                                                                                }
                                                                                // przy okazji - sprawdz, czy HUGE_VAL jest rzeczywiscie max dla double
                                                                                assert(aadOdleglosciParts[iIndex[0]][i] < HUGE_VAL);

                                                                        }
                                                                        else
                                                                        {
                                                                                if (aadOdleglosciParts[i][iIndex[1]] < dMinimum)
                                                                                {
                                                                                        dMinimum = aadOdleglosciParts[i][iIndex[1]];
                                                                                }
                                                                                // przy okazji - sprawdz, czy HUGE_VAL jest rzeczywiscie max dla double
                                                                                assert(aadOdleglosciParts[i][iIndex[1]] < HUGE_VAL);
                                                                        }
                                                                }
                                                        }
                                                        // sprawdz, czy minimum znaleziono - musi takie byc, bo jest cos niedopasowanego
                                                        assert((dMinimum >= 0.0) && (dMinimum < HUGE_VAL));

                                                        // teraz zaznaczamy w tablicy te wszystkie Parts, ktore maja
                                                        // rzeczywiscie te minimalne odleglosci do Part w organizmie iZDrugim
                                                        for (i = 0; i < aiRozmiarCalychGrup[1 - iZDrugim]; i++)
                                                        {
                                                                if (!(m_pMatching->isMatched(
                                                                        1 - iZDrugim,
                                                                        aiKoniecPierwszejGrupy[1 - iZDrugim] + i)))
                                                                {
                                                                        if (iZDrugim == 0)
                                                                        {
                                                                                // teraz niestety musimy rozpoznac odpowiedni organizm
                                                                                // zeby moc indeksowac niesymetryczna tablice
                                                                                if (aadOdleglosciParts[iIndex[0]][i] == dMinimum)
                                                                                {
                                                                                        // jesli w danym miejscu tablicy odleglosci jest faktyczne 
                                                                                        // minimum odleglosci, to mamy potencjalna pare dla iIndex[1]
                                                                                        apvbCzyMinimalnaOdleglosc[0]->operator[](i) = true;
                                                                                }
                                                                        }
                                                                        else
                                                                        {
                                                                                if (aadOdleglosciParts[i][iIndex[1]] == dMinimum)
                                                                                {
                                                                                        apvbCzyMinimalnaOdleglosc[1]->operator[](i) = true;
                                                                                }
                                                                        }
                                                                }
                                                        }

                                                        // a teraz - po znalezieniu potencjalnych elementow do dopasowania
                                                        // dopasowujemy pierwszy z potencjalnych
                                                        iDopasowywany = -1;
                                                        for (i = 0; i < aiRozmiarCalychGrup[1 - iZDrugim]; i++)
                                                        {
                                                                if (apvbCzyMinimalnaOdleglosc[iZDrugim]->operator[](i))
                                                                {
                                                                        iDopasowywany = i;
                                                                        break;
                                                                }
                                                        }
                                                        // musielismy znalezc jakiegos dopasowywanego
                                                        assert((iDopasowywany >= 0) &&
                                                                (iDopasowywany < aiRozmiarCalychGrup[1 - iZDrugim]));

                                                        // no to juz mozemy dopasowac
                                                        m_pMatching->match(
                                                                iZDrugim,
                                                                aiKoniecPierwszejGrupy[iZDrugim] + iIndex[iZDrugim],
                                                                1 - iZDrugim,
                                                                aiKoniecPierwszejGrupy[1 - iZDrugim] + iDopasowywany);
                                                        DB(printf("Przypadek 2.1.: dopasowane Parts dopasowaniebz drugim: (%2i, %2i)\n", aiKoniecPierwszejGrupy[iZDrugim] + iIndex[iZDrugim], aiKoniecPierwszejGrupy[1 - iZDrugim] + iDopasowywany);)

                                                                //aiKoniecPierwszejGrupy[ 1-iBezDrugiego ] + iDopasowywany ;)
                                                                //aiKoniecPierwszejGrupy[ 1-iBezDrugiego ] + iDopasowywany ;)
                                                                // zmniejsz liczby niedopasowanych elementow w grupach
                                                                aiRozmiarGrupy[0]--;
                                                        aiRozmiarGrupy[1]--;

                                                }
                                                else
                                                {
                                                        // jedna z grup sie juz wyczerpala
                                                        // wiec nie ma mozliwosci dopasowania tego drugiego Partu
                                                        /// i trzeba poczekac na zmiane grupy
                                                }

                                                DB(printf("Przypadek 2.\n");)

                                        }// PRZYPADEK 2.
                                        else
                                        {
                                                // PRZYPADEK 3. Zaden z Parts nie ma na liscie drugiego
                                                // AKCJA. Niezalezne dopasowanie obu Parts do pierwszych ze swojej listy

                                                // najpierw dopasujemy do iIndex[0] w organizmie 0
                                                int iDopasowywany = -1;
                                                // poszukujemy pierwszego z listy dopasowania - w organizmie 1
                                                for (i = 0; i < aiRozmiarCalychGrup[1]; i++)
                                                {
                                                        if (apvbCzyMinimalnaOdleglosc[0]->operator[](i))
                                                        {
                                                                iDopasowywany = i;
                                                                break;
                                                        }
                                                }
                                                // musielismy znalezc jakiegos dopasowywanego
                                                assert((iDopasowywany >= 0) &&
                                                        (iDopasowywany < aiRozmiarCalychGrup[1]));

                                                // teraz wlasnie dopasowujemy
                                                m_pMatching->match(
                                                        0,
                                                        aiKoniecPierwszejGrupy[0] + iIndex[0],
                                                        1,
                                                        aiKoniecPierwszejGrupy[1] + iDopasowywany);

                                                // zmniejszamy liczbe niedopasowanych Parts
                                                aiRozmiarGrupy[0]--;
                                                aiRozmiarGrupy[1]--;

                                                // debug - dopasowanie
                                                DB(printf("Przypadek 3.: dopasowane Parts: (%2i, %2i)\n", aiKoniecPierwszejGrupy[0]
                                                        + iIndex[0], aiKoniecPierwszejGrupy[1] + iDopasowywany);)

                                                        // teraz dopasowujemy do iIndex[1] w organizmie 1
                                                        iDopasowywany = -1;
                                                // poszukujemy pierwszego z listy dopasowania - w organizmie 0
                                                for (i = 0; i < aiRozmiarCalychGrup[0]; i++)
                                                {
                                                        if (apvbCzyMinimalnaOdleglosc[1]->operator[](i))
                                                        {
                                                                iDopasowywany = i;
                                                                break;
                                                        }
                                                }
                                                // musielismy znalezc jakiegos dopasowywanego
                                                assert((iDopasowywany >= 0) &&
                                                        (iDopasowywany < aiRozmiarCalychGrup[0]));

                                                // no i teraz realizujemy dopasowanie
                                                m_pMatching->match(
                                                        0,
                                                        aiKoniecPierwszejGrupy[0] + iDopasowywany,
                                                        1,
                                                        aiKoniecPierwszejGrupy[1] + iIndex[1]);

                                                // zmniejszamy liczbe niedopasowanych Parts
                                                aiRozmiarGrupy[0]--;
                                                aiRozmiarGrupy[1]--;

                                                // debug - dopasowanie
                                                DB(printf("Przypadek 3.: dopasowane Parts: (%2i, %2i)\n", aiKoniecPierwszejGrupy[0]
                                                        + iDopasowywany, aiKoniecPierwszejGrupy[1] + iIndex[1]);)


                                        } // PRZYPADEK 3.

                        }// if (! bCzyKoniecGrupy)
                        else
                        {
                                // gdy mamy juz koniec grup - musimy zlikwidowac tablice aadOdleglosciParts
                                // bo za chwilke skonczy sie nam petla
                                for (i = 0; i < aiRozmiarCalychGrup[0]; i++)
                                {
                                        SAFEDELETEARRAY(aadOdleglosciParts[i]);
                                }
                                SAFEDELETEARRAY(aadOdleglosciParts);

                                // musimy tez usunac tablice (wektory) mozliwosci dopasowania
                                SAFEDELETE(apvbCzyMinimalnaOdleglosc[0]);
                                SAFEDELETE(apvbCzyMinimalnaOdleglosc[1]);
                        }
                } // while (! bCzyKoniecGrupy)

                // PO DOPASOWANIU WSZYSTKIEGO Z GRUP (CO NAJMNIEJ JEDNEJ GRUPY W CALOSCI)

                // gdy rozmiar ktorejkolwiek z grup dopasowania spadl do zera
                // to musimy przesunac KoniecPierwszejGrupy (wszystkie dopasowane)
                for (i = 0; i < 2; i++)
                {
                        // grupy nie moga miec mniejszego niz 0 rozmiaru
                        assert(aiRozmiarGrupy[i] >= 0);
                        if (aiRozmiarGrupy[i] == 0)
                                aiKoniecPierwszejGrupy[i] = aiKoniecGrupyDopasowania[i];
                }

                // sprawdzenie warunku konca dopasowywania - gdy nie
                // ma juz w jakims organizmie co dopasowywac
                if (aiKoniecPierwszejGrupy[0] >= m_aiPartCount[0] ||
                        aiKoniecPierwszejGrupy[1] >= m_aiPartCount[1])
                {
                        iCzyDopasowane = 1;
                }
        } // koniec WHILE - petli dopasowania
}

void SimilMeasureGreedy::copyMatching()
{
        SAFEDELETE(tempMatching);
        tempMatching = new SimilMatching(*m_pMatching);
}

void SimilMeasureGreedy::prepareData()
{       
        // create Parts matching object
        m_pMatching = new SimilMatching(models[0]->getPartCount(), models[1]->getPartCount());

        // utworzenie tablicy rozmiarow
        for (int i = 0; i < 2; i++)
        {
                m_aiPartCount[i] = models[i]->getPartCount();
        }
        
        // difference in the number of vertices (Parts) - positive
        // find object that has less parts (m_iSmaller)
        m_iDV = (models[0]->getPartCount() - models[1]->getPartCount());
        if (m_iDV < 0)
                m_iDV = -m_iDV;

        // check if index of the smaller organism is a valid index
        assert((m_iSmaller == 0) || (m_iSmaller == 1));
        // validate difference in the parts number
        assert(m_iDV >= 0);
        
        if (!createPartInfoTables())
                logPrintf("GreedyMeasure", "PrepareData", LOG_ERROR, "Unable to create part info tables.");
        if (!countPartDegrees())
                logPrintf("GreedyMeasure", "PrepareData", LOG_ERROR, "Unable to count part degrees.");
        if (!countPartNeurons())
                logPrintf("GreedyMeasure", "PrepareData", LOG_ERROR, "Unable to count part neurons.");

        DB(printf("Przed sortowaniem:\n");)
        DB(_PrintDegrees(0);)
        DB(_PrintDegrees(1);)

        if (!sortPartInfoTables())
                logPrintf("GreedyMeasure", "PrepareData", LOG_ERROR, "Unable to sort part info tables.");

        DB(printf("Po sortowaniu:\n");)
        DB(_PrintDegrees(0);)
        DB(_PrintDegrees(1);)

        if (m_adFactors[3] == 0)
            with_alignment = false;
}

void SimilMeasureGreedy::cleanData()
{
        // delete degree arrays created in CreatePartInfoTables 
        SAFEDELETEARRAY(m_aDegrees[0]);
        SAFEDELETEARRAY(m_aDegrees[1]);
        
        // in fuzzy mode delete arrays of neighbourhood and fuzzy neighbourhood
        if (isFuzzy)
        {
                for (int i = 0; i != 2; ++i)
                {
                        for (int j = 0; j != models[i]->getPartCount(); ++j)
                        {
                                delete[] m_Neighbours[i][j];
                                delete[] m_fuzzyNeighb[i][j];
                        }
                        delete[] m_Neighbours[i];
                        delete[] m_fuzzyNeighb[i];
                }

        }
        
        // delete created matching
        SAFEDELETE(m_pMatching);
        
        with_alignment = true; //restore default value
}

/** Creates arrays holding information about organisms' Parts (m_aDegrees) andm_Neigh
        fills them with initial data (original indices and zeros).
        Assumptions:
        - Models (models) are created and available.
        */
int SimilMeasureGreedy::createPartInfoTables()
{
        int i, j, partCount;
        // utwórz tablice na informacje o stopniach wierzchołków i liczbie neuroitems
        for (i = 0; i < 2; i++)
        {
                partCount = models[i]->getPartCount();
                // utworz i wypelnij tablice dla Parts wartosciami poczatkowymi
                m_aDegrees[i] = new TDN[partCount];

                if (isFuzzy)
                {
                        m_Neighbours[i] = new int*[partCount];
                        m_fuzzyNeighb[i] = new float*[partCount];
                }

                if (m_aDegrees[i] != NULL && ((!isFuzzy) || (m_Neighbours[i] != NULL && m_fuzzyNeighb[i] != NULL)))
                {
                        // wypelnij tablice zgodnie z sensem TDN[0] - orginalny index
                        // TDN[1], TDN[2], TDN[3] - zerami
                        DB(printf("m_aDegrees[%i]: %p\n", i, m_aDegrees[i]);)
                                for (j = 0; j < partCount; j++)
                                {
                                        m_aDegrees[i][j][0] = j;
                                        m_aDegrees[i][j][1] = 0;
                                        m_aDegrees[i][j][2] = 0;
                                        m_aDegrees[i][j][3] = 0;
                                        m_aDegrees[i][j][4] = 0;

                                        // sprawdz, czy nie piszemy po jakims szalonym miejscu pamieci
                                        assert(m_aDegrees[i][j] != NULL);

                                        if (isFuzzy)
                                        {
                                                m_Neighbours[i][j] = new int[partCount];
                                                for (int k = 0; k < partCount; k++)
                                                {
                                                        m_Neighbours[i][j][k] = 0;
                                                }

                                                m_fuzzyNeighb[i][j] = new float[fuzzyDepth];
                                                for (int k = 0; k < fuzzyDepth; k++)
                                                {
                                                        m_fuzzyNeighb[i][j][k] = 0;
                                                }

                                                assert(m_Neighbours[i][j] != NULL);
                                                assert(m_fuzzyNeighb[i][j] != NULL);
                                        }

                                }
                }
                else
                {
                        logPrintf("ModelSimil", "CreatePartInfoTables", LOG_ERROR, "Not enough memory?");
                        return 0;
                }
                // wypelnij tablice OnJoints i Anywhere wartościami początkowymi
                // OnJoint
                m_aOnJoint[i][0] = 0;
                m_aOnJoint[i][1] = 0;
                m_aOnJoint[i][2] = 0;
                m_aOnJoint[i][3] = 0;
                // Anywhere
                m_aAnywhere[i][0] = 0;
                m_aAnywhere[i][1] = 0;
                m_aAnywhere[i][2] = 0;
                m_aAnywhere[i][3] = 0;
        }
        return 1;
}

/** Computes degrees of Parts of both organisms. Fills in the m_aDegrees arrays
        with proper information about degrees.
        Assumptions:
        - Models (models) are created and available.
        - Arrays m_aDegrees are created.
        */
int SimilMeasureGreedy::countPartDegrees()
{
        Part *P1, *P2;
        int i, j, i1, i2;

        // dla obu stworzen oblicz stopnie wierzcholkow
        for (i = 0; i < 2; i++)
        {
                // dla wszystkich jointow
                for (j = 0; j < models[i]->getJointCount(); j++)
                {
                        // pobierz kolejny Joint
                        Joint *J = models[i]->getJoint(j);
                        // wez jego konce
                        P1 = J->part1;
                        P2 = J->part2;
                        // znajdz ich indeksy w Modelu (indeksy orginalne)
                        i1 = models[i]->findPart(P1);
                        i2 = models[i]->findPart(P2);
                        // zwieksz stopien odpowiednich Parts
                        m_aDegrees[i][i1][DEG]++;
                        m_aDegrees[i][i2][DEG]++;
                        m_aDegrees[i][i1][FUZ_DEG]++;
                        m_aDegrees[i][i2][FUZ_DEG]++;
                        if (isFuzzy)
                        {
                                m_Neighbours[i][i1][i2] = 1;
                                m_Neighbours[i][i2][i1] = 1;
                        }
                }
                // dla elementow nie osadzonych na Parts (OnJoint, Anywhere) -
                // stopnie wierzchołka są już ustalone na zero
        }

        if (isFuzzy)
        {
                countFuzzyNeighb();
        }

        return 1;
}

void SimilMeasureGreedy::getNeighbIndexes(int mod, int partInd, std::vector<int> &indexes)
{
        indexes.clear();
        int i, size = models[mod]->getPartCount();

        for (i = 0; i < size; i++)
        {
                if (m_Neighbours[mod][partInd][i] > 0)
                {
                        indexes.push_back(i);
                }
        }
}

int compareFuzzyRows(const void *pa, const void *pb)
{
        float **a = (float**)pa;
        float **b = (float**)pb;

        for (int i = 1; i < fuzzDepth; i++)
        {
                if (a[0][i] > b[0][i])
                {
                        return -1;
                }
                if (a[0][i] < b[0][i])
                {
                        return 1;
                }
        }

        return 0;
}

void SimilMeasureGreedy::fuzzyOrder()
{
        int i, depth, partInd, prevPartInd, partCount;
        for (int mod = 0; mod < 2; mod++)
        {
                partCount = models[mod]->getPartCount();
                partInd = m_fuzzyNeighb[mod][partCount - 1][0];
                m_aDegrees[mod][partInd][FUZ_DEG] = 0;

                for (i = (partCount - 2); i >= 0; i--)
                {
                        prevPartInd = partInd;
                        partInd = m_fuzzyNeighb[mod][i][0];
                        m_aDegrees[mod][partInd][FUZ_DEG] = m_aDegrees[mod][prevPartInd][FUZ_DEG];
                        for (depth = 1; depth < fuzzyDepth; depth++)
                        {
                                if (m_fuzzyNeighb[mod][i][depth] != m_fuzzyNeighb[mod][i + 1][depth])
                                {
                                        m_aDegrees[mod][partInd][FUZ_DEG]++;
                                        break;
                                }
                        }
                }
        }
}

//sort according to fuzzy degree
void SimilMeasureGreedy::sortFuzzyNeighb()
{
        fuzzDepth = fuzzyDepth;
        for (int mod = 0; mod < 2; mod++)
        {
                std::qsort(m_fuzzyNeighb[mod], (size_t)models[mod]->getPartCount(), sizeof(m_fuzzyNeighb[mod][0]), compareFuzzyRows);
        }
}

//computes fuzzy vertex degree
void SimilMeasureGreedy::countFuzzyNeighb()
{
        std::vector<int> nIndexes;
        float newDeg = 0;

        for (int mod = 0; mod < 2; mod++)
        {
                int partCount = models[mod]->getPartCount();

                for (int depth = 0; depth < fuzzyDepth; depth++)
                {
                        //use first column for storing indices
                        for (int partInd = 0; partInd < partCount; partInd++)
                        {
                                if (depth == 0)
                                {
                                        m_fuzzyNeighb[mod][partInd][depth] = partInd;
                                }
                                else if (depth == 1)
                                {
                                        m_fuzzyNeighb[mod][partInd][depth] = m_aDegrees[mod][partInd][DEG];
                                }
                                else
                                {
                                        getNeighbIndexes(mod, partInd, nIndexes);
                                        for (unsigned int k = 0; k < nIndexes.size(); k++)
                                        {
                                                newDeg += m_fuzzyNeighb[mod][nIndexes.at(k)][depth - 1];
                                        }
                                        newDeg /= nIndexes.size();
                                        m_fuzzyNeighb[mod][partInd][depth] = newDeg;
                                        for (int mod = 0; mod < 2; mod++)
                                        {
                                                int partCount = models[mod]->getPartCount();
                                                for (int i = partCount - 1; i >= 0; i--)
                                                {

                                                }
                                        }
                                        newDeg = 0;
                                }
                        }
                }
        }

        sortFuzzyNeighb();
        fuzzyOrder();
}

/** Computes numbers of neurons and neurons' inputs for each Part of each
        organisms and fills in the m_aDegrees array.
        Assumptions:
        - Models (models) are created and available.
        - Arrays m_aDegrees are created.
        */
int SimilMeasureGreedy::countPartNeurons()
{
        // sprawdz zalozenie - o modelach
        assert((models[0] != NULL) && (models[1] != NULL));
        assert(models[0]->isValid() && models[1]->isValid());

        // sprawdz zalozenie - o tablicach
        assert(m_aDegrees[0] != NULL);
        assert(m_aDegrees[1] != NULL);
        
        Part *P1;
        Joint *J1;
        int i, j, i2, neuro_connections;

        // dla obu stworzen oblicz liczbe Neurons + connections dla Parts
        // a takze dla OnJoints i Anywhere
        for (i = 0; i < 2; i++)
        {
                for (j = 0; j < models[i]->getNeuroCount(); j++)
                {
                        // pobierz kolejny Neuron
                        Neuro *N = models[i]->getNeuro(j);
                        // policz liczbe jego wejść i jego samego tez
                        // czy warto w ogole liczyc polaczenia...? co to da/spowoduje?
                        neuro_connections = N->getInputCount() + 1;
                        // wez Part, na ktorym jest Neuron
                        P1 = N->getPart();
                        if (P1)
                        {
                                // dla neuronow osadzonych na Partach
                                i2 = models[i]->findPart(P1); // znajdz indeks Part w Modelu
                                m_aDegrees[i][i2][2] += neuro_connections; // zwieksz liczbe Connections+Neurons dla tego Part (TDN[2])
                                m_aDegrees[i][i2][3]++; // zwieksz liczbe Neurons dla tego Part (TDN[3])
                        }
                        else
                        {
                                // dla neuronow nie osadzonych na partach
                                J1 = N->getJoint();
                                if (J1)
                                {
                                        // dla tych na Jointach
                                        m_aOnJoint[i][2] += neuro_connections; // zwieksz liczbe Connections+Neurons
                                        m_aOnJoint[i][3]++; // zwieksz liczbe Neurons
                                }
                                else
                                {
                                        // dla tych "gdziekolwiek"
                                        m_aAnywhere[i][2] += neuro_connections; // zwieksz liczbe Connections+Neurons
                                        m_aAnywhere[i][3]++; // zwieksz liczbe Neurons
                                }
                        }
                }
        }
        return 1;
}

/** Sorts arrays m_aDegrees (for each organism) by Part's degree, and then by
        number of neural connections and neurons in groups of Parts with the same
        degree.
        Assumptions:
        - Models (models) are created and available.
        - Arrays m_aDegrees are created.
        @saeDegrees, CompareItemNo
        */
int SimilMeasureGreedy::sortPartInfoTables()
{
        int i;
        int(*pfDegreeFunction) (const void*, const void*) = NULL;
        pfDegreeFunction = isFuzzy ? &compareFuzzyDegrees : &compareDegrees;
        // sortowanie obu tablic wg stopni punktów - TDN[1]
        for (i = 0; i < 2; i++)
        {
                DB(_PrintDegrees(i));
                std::qsort(m_aDegrees[i], (size_t)(models[i]->getPartCount()),
                        sizeof(TDN), pfDegreeFunction);
                DB(_PrintDegrees(i));
        }

        // sprawdzenie wartosci parametru
        DB(i = sizeof(TDN);)
                int degreeType = isFuzzy ? FUZ_DEG : DEG;

        // sortowanie obu tablic m_aDegrees wedlug liczby neuronów i
        // czesci neuronu - ale w obrebie grup o tym samym stopniu
        for (i = 0; i < 2; i++)
        {
                int iPocz = 0;
                int iDeg, iNewDeg, iPartCount, j;
                // stopien pierwszego punktu w tablicy Degrees  odniesienie
                iDeg = m_aDegrees[i][0][degreeType];
                iPartCount = models[i]->getPartCount();
                // po kolei dla kazdego punktu w organizmie
                for (j = 0; j <= iPartCount; j++)
                {
                        // sprawdz stopien punktu (lub nadaj 0 - gdy juz koniec tablicy)
                        //                      iNewDeg = (j != iPartCount) ? m_aDegrees[i][j][1] : 0;
                        // usunieto stara wersje porownania!!! wprowadzono znak porownania <

                        iNewDeg = (j < iPartCount) ? m_aDegrees[i][j][degreeType] : 0;
                        // skoro tablice sa posortowane wg stopni, to mamy na pewno taka kolejnosc
                        assert(iNewDeg <= iDeg);
                        if (iNewDeg != iDeg)
                        {
                                // gdy znaleziono koniec grupy o tym samym stopniu
                                // sortuj po liczbie neuronow w obrebie grupy
                                DB(_PrintDegrees(i));
                                DB(printf("qsort( poczatek=%i, rozmiar=%i, sizeof(TDN)=%i)\n", iPocz, (j - iPocz), sizeof(TDN));)
                                        // wyswietlamy z jedna komorka po zakonczeniu tablicy
                                        DB(_PrintArray(m_aDegrees[i][iPocz], 0, (j - iPocz) * 4);)

                                        std::qsort(m_aDegrees[i][iPocz], (size_t)(j - iPocz),
                                                sizeof(TDN), SimilMeasureGreedy::compareConnsNo);
                                DB(_PrintDegrees(i));
                                // wyswietlamy z jedna komorka po zakonczeniu tablicy
                                DB(_PrintArray(m_aDegrees[i][iPocz], 0, (j - iPocz) * 4);)
                                        // rozpocznij nowa grupe
                                        iPocz = j;
                                iDeg = iNewDeg;
                        }
                }
        }
        return 1;
}


/** Prints the state of the matching object. Debug method.
 */
void SimilMeasureGreedy::_PrintPartsMatching()
{
        // assure that matching exists
        assert(m_pMatching != NULL);

        printf("Parts matching:\n");
        m_pMatching->printMatching();
}

/** Comparison function required for qsort() call. Used while sorting groups of
        Parts with respect to degree. Compares two TDN structures
        with respect to their [1] field (degree). Highest degree goes first.
        @param pElem1 Pointer to the TDN structure of some Part.
        @param pElem2 Pointer to the TDN structure of some Part.
        @return (-1) - pElem1 should go first, 0 - equal, (1) - pElem2 should go first.
        */
int SimilMeasureGreedy::compareDegrees(const void *pElem1, const void *pElem2)
{
        int *tdn1 = (int *)pElem1;
        int *tdn2 = (int *)pElem2;

        if (tdn1[1] > tdn2[1])
        {
                // when degree - tdn1[1] - is BIGGER
                return -1;
        }
        else
                if (tdn1[1] < tdn2[1])
                {
                        // when degree - tdn2[1] - is BIGGER
                        return 1;
                }
                else
                {
                        return 0;
                }
}

/** Comparison function required for qsort() call. Used while sorting groups of
        Parts with respect to fuzzy vertex degree. Compares two TDN structures
        with respect to their [4] field ( fuzzy vertex degree). Highest degree goes first.
        @param pElem1 Pointer to the TDN structure of some Part.
        @param pElem2 Pointer to the TDN structure of some Part.
        @return (-1) - pElem1 should go first, 0 - equal, (1) - pElem2 should go first.
        */
int SimilMeasureGreedy::compareFuzzyDegrees(const void *pElem1, const void *pElem2)
{
        int *tdn1 = (int *)pElem1;
        int *tdn2 = (int *)pElem2;

        if (tdn1[4] > tdn2[4])
        {
                // when degree - tdn1[4] - is BIGGER
                return -1;
        }
        else
                if (tdn1[4] < tdn2[4])
                {
                        // when degree - tdn2[4] - is BIGGER
                        return 1;
                }
                else
                {
                        return 0;
                }
}

/** Comparison function required for qsort() call. Used while sorting groups of Parts with
                the same degree. Firstly, compare NIt. Secondly, compare Neu. If both are equal -
                compare Parts' original index (they are never equal). So this sorting assures
                that the order obtained is deterministic.
                @param pElem1 Pointer to the TDN structure of some Part.
                @param pElem2 Pointer to the TDN structure of some Part.
                @return (-1) - pElem1 should go first, 0 - equal, (1) - pElem2 should go first.
                */
int SimilMeasureGreedy::compareConnsNo(const void *pElem1, const void *pElem2)
{
        // pointers to TDN arrays 
        int *tdn1, *tdn2;
        // definitions of elements being compared
        tdn1 = (int *)pElem1;
        tdn2 = (int *)pElem2;

        // comparison according to Neural Connections (to jest TDN[2])
        if (tdn1[NEUR_CONNS] > tdn2[NEUR_CONNS])
        {
                // when number of NConn Elem1 is BIGGER
                return -1;
        }
        else
                if (tdn1[NEUR_CONNS] < tdn2[NEUR_CONNS])
                {
                        // when number of NConn Elem1 is SMALLER
                        return 1;
                }
                else
                {
                        // when numbers of NConn are EQUAL
                        // compare Neu numbers (TDN[3])
                        if (tdn1[NEUROS] > tdn2[NEUROS])
                        {
                                // when number of Neu is BIGGER for Elem1 
                                return -1;
                        }
                        else
                                if (tdn1[NEUROS] < tdn2[NEUROS])
                                {
                                        // when number of Neu is SMALLER for Elem1 
                                        return 1;
                                }
                                else
                                {
                                        // when numbers of Nconn and Neu are equal we check original indices 
                                        // of Parts being compared

                                        // comparison according to OrgIndex
                                        if (tdn1[ORG_IND] > tdn2[ORG_IND])
                                        {
                                                // when the number of NIt Deg1 id BIGGER
                                                return -1;
                                        }
                                        else
                                                if (tdn1[ORG_IND] < tdn2[ORG_IND])
                                                {
                                                        // when the number of NIt Deg1 id SMALLER
                                                        return 1;
                                                }
                                                else
                                                {
                                                        // impossible, indices are alway different
                                                        return 0;
                                                }
                                }
                }
}

/** Returns number of factors involved in final distance computation.
                These factors include differences in numbers of parts, degrees,
                number of neurons.
                */
int SimilMeasureGreedy::getNOFactors()
{
        return SimilMeasureGreedy::iNOFactors;
}

/** Prints the array of degrees for the given organism. Debug method.
 */
void SimilMeasureGreedy::_PrintDegrees(int i)
{
        int j;
        printf("Organizm %i :", i);
        printf("\n      ");
        for (j = 0; j < models[i]->getPartCount(); j++)
                printf("%3i ", j);
        printf("\nInd:  ");
        for (j = 0; j < models[i]->getPartCount(); j++)
                printf("%3i ", (int)m_aDegrees[i][j][0]);
        printf("\nDeg:  ");
        for (j = 0; j < models[i]->getPartCount(); j++)
                printf("%3i ", (int)m_aDegrees[i][j][1]);
        printf("\nNIt:  ");
        for (j = 0; j < models[i]->getPartCount(); j++)
                printf("%3i ", (int)m_aDegrees[i][j][2]);
        printf("\nNeu:  ");
        for (j = 0; j < models[i]->getPartCount(); j++)
                printf("%3i ", (int)m_aDegrees[i][j][3]);
        printf("\n");
}

/** Prints one array of ints. Debug method.
        @param array Base pointer of the array.
        @param base First index of the array's element.
        @param size Number of elements to print.
        */
void SimilMeasureGreedy::_PrintArray(int *array, int base, int size)
{
        int i;
        for (i = base; i < base + size; i++)
        {
                printf("%i ", array[i]);
        }
        printf("\n");
}

void SimilMeasureGreedy::_PrintArrayDouble(double *array, int base, int size)
{
        int i;
        for (i = base; i < base + size; i++)
        {
                printf("%f ", array[i]);
        }
        printf("\n");
}

/** Prints one array of parts neighbourhood.
        @param index of organism
        */
void SimilMeasureGreedy::_PrintNeighbourhood(int o)
{
        assert(m_Neighbours[o] != 0);
        printf("Neighbourhhod of organism %i\n", o);
        int size = models[o]->getPartCount();
        for (int i = 0; i < size; i++)
        {
                for (int j = 0; j < size; j++)
                {
                        printf("%i ", m_Neighbours[o][i][j]);
                }
                printf("\n");
        }
}

/** Prints one array of parts fuzzy neighbourhood.
        @param index of organism
        */
void SimilMeasureGreedy::_PrintFuzzyNeighbourhood(int o)
{
        assert(m_fuzzyNeighb[o] != NULL);
        printf("Fuzzy neighbourhhod of organism %i\n", o);
        int size = models[o]->getPartCount();
        for (int i = 0; i < size; i++)
        {
                for (int j = 0; j < fuzzyDepth; j++)
                {
                        printf("%f ", m_fuzzyNeighb[o][i][j]);
                }
                printf("\n");
        }
}

int SimilMeasureGreedy::setParams(std::vector<double> params)
{
        int i = 0;
        for (i = 0; i < SimilMeasureGreedy::iNOFactors; i++)
                m_adFactors[i] = params.at(i);
        fixedZaxis = params.at(i);
        return 0;
}