source: cpp/frams/model/similarity/measure-mds-based.cpp @ 1044

Last change on this file since 1044 was 1044, checked in by oriona, 3 months ago

Similarity measures code refactored. Distribution-based similarity measure added.

File size: 5.8 KB
Line 
1// This file is a part of Framsticks SDK.  http://www.framsticks.com/
2// Copyright (C) 1999-2020  Maciej Komosinski and Szymon Ulatowski.
3// See LICENSE.txt for details.
4
5#include "measure-mds-based.h"
6#include "SVD/matrix_tools.h"
7
8SimilMeasureMDSBased::SimilMeasureMDSBased()
9{
10        coordinates[0] = nullptr;
11        coordinates[1] = nullptr;
12       
13        with_alignment = true;
14        save_matching = false;
15        fixedZaxis = 0;
16        wMDS = 1;
17}
18
19double SimilMeasureMDSBased::getDistance()
20{
21        double dResult = -1;
22        m_iSmaller = models[0]->getPartCount() <= models[1]->getPartCount() ? 0 : 1;
23
24        prepareData();
25       
26        if (with_alignment)
27        {
28                if (!getPartPositions())
29                {
30                        logPrintf("SimilMDSBasedMeasure", "getDistance", LOG_ERROR, "Cannot compute part positions");
31                        return -1;
32                }
33               
34                if (!computePartsPositionsByMDS())
35                {
36                        logPrintf("SimilMDSBasedMeasure", "getDistance", LOG_ERROR, "Cannot perform MDS");
37                        return -1;
38                }
39
40                // tutaj zacznij pętlę po przekształceniach  geometrycznych
41                const int NO_OF_TRANSFORM = 8; // liczba transformacji geometrycznych (na razie tylko ID i O_YZ)
42                // tablice transformacji współrzędnych; nie są to dokładnie tablice tranformacji, ale raczej tablice PRZEJŚĆ
43                // pomiędzy transformacjami;
44                const int dMulX[NO_OF_TRANSFORM] = { 1, -1, -1, 1, -1, 1, -1, -1 };
45                const int dMulY[NO_OF_TRANSFORM] = { 1, 1, -1, -1, -1, -1, -1, 1 };
46                const int dMulZ[NO_OF_TRANSFORM] = { 1, 1, 1, -1, -1, -1, 1, 1 };
47
48                #ifdef max
49                #undef max //this macro would conflict with line below
50                #endif
51                double dMinSimValue = std::numeric_limits<double>::max(); // minimum value of similarity
52               
53                int iTransform; // a counter of geometric transformations
54                for (iTransform = 0; iTransform < NO_OF_TRANSFORM; iTransform++)
55                {
56                        beforeTransformation();
57                       
58                        for (int iPart = 0; iPart < models[m_iSmaller]->getPartCount(); iPart++)
59                        {
60                                // for each iPart, a part of the model iMod
61                                coordinates[m_iSmaller][iPart].x *= dMulX[iTransform];
62                                coordinates[m_iSmaller][iPart].y *= dMulY[iTransform];
63                                coordinates[m_iSmaller][iPart].z *= dMulZ[iTransform];
64                        }
65                        // now the positions are recomputed
66                        double dCurrentSim = distanceForTransformation();
67
68                        // załóż poprawną wartość podobieństwa
69                        assert(dCurrentSim >= 0.0);
70
71                        // porównaj wartość obliczoną z dotychczasowym minimum
72                        if (dCurrentSim < dMinSimValue)
73                        {
74                        //printf("lesser sim\n");
75                                dMinSimValue = dCurrentSim;
76                                if (save_matching)
77                                        copyMatching();
78                        }
79                }
80               
81                dResult = dMinSimValue;
82               
83                SAFEDELETEARRAY(coordinates[0]);
84                SAFEDELETEARRAY(coordinates[1]);
85        }
86        else
87        {
88                dResult = distanceWithoutAlignment();
89        }
90
91        cleanData();
92
93        return dResult;
94}
95
96
97/** Gets information about Parts' positions in 3D from models into the arrays
98                coordinates.
99                Assumptions:
100                - Models (models) are created and available.
101                - Arrays coordinates are created.
102                @return 1 if success, otherwise 0.
103                */
104int SimilMeasureMDSBased::getPartPositions()
105{       
106        // for each model copy its part coordinates
107        Part *pPart;
108        for (int iMod = 0; iMod < 2; iMod++)
109        {
110                // utworz tablice dla pozycji 3D Parts (wielkosc tablicy: liczba Parts organizmu)
111                coordinates[iMod] = new Pt3D[models[iMod]->getPartCount()];
112                assert(coordinates[iMod] != NULL);
113                // dla każdego z modeli iMod
114                for (int iPart = 0; iPart < models[iMod]->getPartCount(); iPart++)
115                {
116                        // dla każdego iPart organizmu iMod
117                        // pobierz tego Part
118                        pPart = models[iMod]->getPart(iPart);
119                        // zapamietaj jego pozycje 3D w tablicy coordinates
120                        coordinates[iMod][iPart].x = pPart->p.x;
121                        coordinates[iMod][iPart].y = pPart->p.y;
122                        coordinates[iMod][iPart].z = pPart->p.z;
123                }
124        }
125
126        return 1;
127}
128
129bool SimilMeasureMDSBased::computePartsPositionsByMDS()
130{
131        bool bResult = true;
132       
133        // use MDS to obtain different point of view on organisms
134        // and store the new positions (currently the original ones are still stored)
135        for (int iMod = 0; iMod < 2; iMod++)
136        {
137                // prepare the vector of errors of approximation for the SVD
138                std::vector<double> vEigenvalues;
139                int nSize = models[iMod]->getPartCount();
140
141                double *pDistances = new double[nSize * nSize];
142
143                for (int i = 0; i < nSize; i++)
144                {
145                        pDistances[i] = 0;
146                }
147
148                Model *pModel = models[iMod]; // use the model of the iMod (current) organism
149                int iP1, iP2; // indices of Parts in the model
150                Part *P1, *P2; // pointers to Parts
151                Pt3D P1Pos, P2Pos; // positions of parts
152                double dDistance; // the distance between Parts
153
154                double *weights = new double[nSize];
155                for (int i = 0; i < nSize; i++)
156                {
157                        if (wMDS == 1)
158                                weights[i] = 0;
159                        else
160                                weights[i] = 1;
161                }
162
163                if (wMDS == 1)
164                        for (int i = 0; i < pModel->getJointCount(); i++)
165                        {
166                                weights[pModel->getJoint(i)->p1_refno]++;
167                                weights[pModel->getJoint(i)->p2_refno]++;
168                        }
169
170                for (iP1 = 0; iP1 < pModel->getPartCount(); iP1++)
171                {
172                        // for each iP1, a Part index in the model of organism iMod
173                        P1 = pModel->getPart(iP1);
174                        // get the position of the Part
175                        P1Pos = P1->p;
176                        if (fixedZaxis == 1)
177                        {
178                                P1Pos.z = 0; //fixed vertical axis, so pretend all points are on the xy plane
179                        }
180                        for (iP2 = 0; iP2 < pModel->getPartCount(); iP2++)
181                        {
182                                // for each (iP1, iP2), a pair of Parts index in the model
183                                P2 = pModel->getPart(iP2);
184                                // get the position of the Part
185                                P2Pos = P2->p;
186                                if (fixedZaxis == 1)
187                                {
188                                        P2Pos.z = 0; //fixed vertical axis, so pretend all points are on the xy plane
189                                }
190                                // compute the geometric (Euclidean) distance between the Parts
191                                dDistance = P1Pos.distanceTo(P2Pos);
192                                // store the distance
193                                pDistances[iP1 * nSize + iP2] = dDistance;
194                        } // for (iP2)
195                } // for (iP1)
196               
197                MatrixTools::weightedMDS(vEigenvalues, nSize, pDistances, coordinates[iMod], weights);
198               
199                if (fixedZaxis == 1) //restore the original vertical coordinate of each Part
200                {
201                        for (int part = 0; part < pModel->getPartCount(); part++)
202                        {
203                                coordinates[iMod][part].z = pModel->getPart(part)->p.z;
204                        }
205                }
206               
207                delete[] pDistances;
208                delete[] weights;
209        }
210
211        return bResult;
212}
213
Note: See TracBrowser for help on using the repository browser.