source: cpp/frams/model/similarity/measure-distribution.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: 7.5 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-distribution.h"
6#include <common/nonstd_math.h>
7#include <limits>
8#include "EMD/emd.c"
9#include <iostream>
10
11#define PI 3.1415926535
12
13#define FIELDSTRUCT SimilMeasureDistribution
14
15static ParamEntry simil_distribution_paramtab[] = {
16                { "Similarity: distribution", 1, 4, "SimilMeasureDistribution", "Evaluates morphological dissimilarity using distribution measure.", },
17                { "simil_density", 0, 0, "Density of surface sampling", "f 1 100 10", FIELD(density), "", },
18                { "simil_bin_num", 0, 0, "Number of bins", "d 1 1000 128", FIELD(bin_num), "", },
19                { "simil_samples_num", 0, 0, "Number of samples", "d 1 1048576 1048576", FIELD(samples_num), "", },
20                { "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.", },
21                { 0, },
22};
23
24#undef FIELDSTRUCT
25
26SimilMeasureDistribution::SimilMeasureDistribution() : localpar(simil_distribution_paramtab, this)
27{
28        localpar.setDefault();
29        SimilMeasureDistribution::distribution_fun = &SimilMeasureDistribution::D2; //D1 and D2 are the best descriptors
30}
31
32double SimilMeasureDistribution::getDistance()
33{
34        double dist = 0;
35        for (int i = 0; i < 2; i++)
36        {
37                funs[i] = new std::pair<double, float>[bin_num]();
38                for (int j = 0; j < bin_num; j++)
39                        funs[i][j] = std::make_pair(0, 0);
40        }
41       
42        for (int i = 0; i < 2; i++)
43                sst_models[i] = new SolidsShapeTypeModel((*models[i]));
44       
45        SimilMeasureDistribution::calculateFuns();     
46        dist = SimilMeasureDistribution::compareFuns();
47       
48        for (int i = 0; i < 2; i++)
49        {
50                SAFEDELETE(sst_models[i]);
51                SAFEDELETEARRAY(funs[i]);
52        }
53        return dist;
54}
55
56int SimilMeasureDistribution::setParams(std::vector<double> params)
57{
58        for (unsigned int i = 0; i < params.size(); i++)
59                if (params.at(i) <= 0)
60                {
61                        logPrintf("SimilDistributionMeasure", "setParams", LOG_ERROR, "Param values should be larger than 0.");
62                        return -1;
63                }
64       
65        density = params.at(0);
66        bin_num = params.at(1);
67        samples_num = params.at(2);
68       
69        return 0;
70}
71
72void SimilMeasureDistribution::calculateFun(std::pair<double, float> *fun, Model &sampled)
73{
74        int size = sampled.getPartCount();
75        int samples_taken = samples_num;
76
77        //Check if total number of points pairs is smaller than samples number
78        //but first, prevent exceeding int limits
79        if (size < (int) sqrt((double) std::numeric_limits<int>::max()))
80                samples_taken = min(samples_num, size*size);
81
82        //Get sampled distribution
83        std::vector<double> dist_vect = (this->*SimilMeasureDistribution::distribution_fun)(samples_taken, &sampled);
84
85        auto result = std::minmax_element(dist_vect.begin(), dist_vect.end());
86        double min = *result.first;
87        double max = *result.second;
88
89        //Create histogram
90        float hist[bin_num] = {};       
91        int ind = 0;
92
93        for (unsigned int j = 0; j < dist_vect.size(); j++)
94        {
95                ind = (int) std::floor((dist_vect.at(j)-min)*1/(max-min)*bin_num);
96                if (ind <= (bin_num-1))
97                        hist[ind] += 1;
98                else if (ind == bin_num)
99                        hist[bin_num-1] += 1;
100        }
101
102        //Create pairs
103        for (int j = 0; j < bin_num; j++)
104        {
105                fun[j] = std::make_pair(min+(max-min)/bin_num*(j+0.5), hist[j]);
106        }
107
108        //Normalize
109        float total_mass = 0;
110        for (int j = 0; j < bin_num; j++)
111    {
112        total_mass += fun[j].second;
113    }
114
115        for (int j = 0; j < bin_num; j++)
116                fun[j].second /= total_mass;
117}
118
119void SimilMeasureDistribution::calculateFuns()
120{
121        for (int i = 0; i < 2; i++)
122        {
123                Model sampled = SimilMeasureDistribution::sampleSurface(&sst_models[i]->getModel(), density);
124                SimilMeasureDistribution::calculateFun(funs[i], sampled);
125        }
126}
127
128double SimilMeasureDistribution::compareFuns()
129{
130        return SimilMeasureDistribution::EMD(funs[0], funs[1]);
131}
132
133vector<double> SimilMeasureDistribution::D1(int samples_taken, Model *sampled)
134{
135        vector<double> dist_vect;
136        int size = sampled->getPartCount();
137        double x = 0;
138        double y = 0;
139        double z = 0;
140       
141        for (int i = 0; i < size; i++)
142        {
143                Pt3D pos = sampled->getPart(i)->p;
144                x += pos.x;
145                y += pos.y;
146                z += pos.z;
147        }
148       
149        x = x/size;
150        y = y/size;
151        z = z/size;
152       
153        Pt3D centroid = {x, y, z};
154       
155        for (int i = 0; i < samples_taken; i++)
156        {
157                int p1 = rndUint(size);
158                double dist = sampled->getPart(p1)->p.distanceTo(centroid);
159                if (dist > 0)
160                        dist_vect.push_back(dist);
161        }
162       
163        return dist_vect;
164}
165
166vector<double> SimilMeasureDistribution::D2(int samples_taken, Model *sampled)
167{
168        vector<double> dist_vect;
169        int size = sampled->getPartCount();
170        for (int i = 0; i < samples_taken; i++)
171        {
172                int p1 = rndUint(size);
173                int p2 = rndUint(size);
174                double dist = sampled->getPart(p1)->p.distanceTo(sampled->getPart(p2)->p);
175                if (dist > 0)
176                        dist_vect.push_back(dist);
177        }
178       
179        return dist_vect;
180}
181
182vector<double> SimilMeasureDistribution::D3(int samples_taken, Model *sampled)
183{
184        vector<double> dist_vect;
185        int size = sampled->getPartCount();
186        for (int i = 0; i < samples_taken; i++)
187        {
188                int p1 = rndUint(size);
189                int p2 = rndUint(size);
190                int p3 = rndUint(size);
191               
192                Pt3D v(sampled->getPart(p2)->p);
193                Pt3D w(sampled->getPart(p3)->p);
194                v -= sampled->getPart(p1)->p;
195                w -= sampled->getPart(p1)->p;
196                Pt3D cross_prod(0);
197                cross_prod.vectorProduct(v, w);
198               
199                double dist = 0.5 * cross_prod.length();
200                if (dist > 0)
201                        dist_vect.push_back(dist);
202        }
203       
204        return dist_vect;
205}
206
207vector<double> SimilMeasureDistribution::D4(int samples_taken, Model *sampled)
208{
209        vector<double> dist_vect;
210        int size = sampled->getPartCount();
211        for (int i = 0; i < samples_taken; i++)
212        {
213                int a = rndUint(size);
214                int b = rndUint(size);
215                int c = rndUint(size);
216                int d = rndUint(size);
217               
218                Pt3D ad(sampled->getPart(a)->p);
219                Pt3D bd(sampled->getPart(b)->p);
220                Pt3D cd(sampled->getPart(c)->p);
221               
222                ad -= sampled->getPart(d)->p;
223                bd -= sampled->getPart(d)->p;
224                cd -= sampled->getPart(d)->p;
225               
226                Pt3D cross_prod(0);
227                cross_prod.vectorProduct(bd, cd);
228                cross_prod.entrywiseProduct(ad);
229               
230                double dist = cross_prod.length()/6;
231                if (dist > 0)
232                        dist_vect.push_back(dist);
233        }
234       
235        return dist_vect;
236}
237
238vector<double> SimilMeasureDistribution::A3(int samples_taken, Model *sampled)
239{
240        vector<double> dist_vect;
241        int size = sampled->getPartCount();
242        for (int i = 0; i < samples_taken; i++)
243        {
244                int p1 = rndUint(size);
245                int p2 = rndUint(size);
246                int p3 = rndUint(size);
247                double a = sampled->getPart(p1)->p.distanceTo(sampled->getPart(p3)->p);
248                double b = sampled->getPart(p3)->p.distanceTo(sampled->getPart(p2)->p);
249                double c = sampled->getPart(p1)->p.distanceTo(sampled->getPart(p2)->p);
250                double beta = acos((a*a + b*b - c*c)/(2*a*b));
251       
252                if (!std::isnan(beta))
253                        dist_vect.push_back(beta);
254        }
255       
256        return dist_vect;
257}
258
259
260float dist(feature_t* F1, feature_t* F2)
261{
262        return abs((*F1)-(*F2));
263}
264
265
266void SimilMeasureDistribution::fillPointsWeights(std::pair<double, float> *fun, feature_t *points, float *weights)
267{
268        for (int j = 0; j < bin_num; j++)
269        {
270                points[j] = {fun[j].first};
271                weights[j] = fun[j].second;
272        }
273}
274
275double SimilMeasureDistribution::EMD(std::pair<double, float> *fun1, std::pair<double, float> *fun2)
276{
277        feature_t *points[2];
278        float *weights[2];
279       
280        for (int i = 0; i < 2; i++)
281        {
282                points[i] = new feature_t[bin_num];
283                weights[i] = new float[bin_num]();
284        }
285        SimilMeasureDistribution::fillPointsWeights(fun1, points[0], weights[0]);
286        SimilMeasureDistribution::fillPointsWeights(fun2, points[1], weights[1]);
287       
288        signature_t sig1 = {bin_num, points[0], weights[0]},
289        sig2 = {bin_num, points[1], weights[1]};
290
291        float e = emd(&sig1, &sig2, dist, 0, 0);
292       
293        for (int i = 0; i < 2; i++)
294        {
295                delete[] points[i];
296                delete[] weights[i];
297        }
298
299        return e;
300}
Note: See TracBrowser for help on using the repository browser.