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

Last change on this file since 1052 was 1048, checked in by Maciej Komosinski, 3 years ago

SimilMeasure? -> SimilMeasureBase?; introduced a new class (SimilMeasure?) that allows scripts to access all similarity measures; a number of minor fixes and improvements

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