ACloudViewer  3.9.4
A Modern Library for 3D Data Processing
CutPursuit_KL.h
Go to the documentation of this file.
1 // ----------------------------------------------------------------------------
2 // - CloudViewer: www.cloudViewer.org -
3 // ----------------------------------------------------------------------------
4 // Copyright (c) 2018-2024 www.cloudViewer.org
5 // SPDX-License-Identifier: MIT
6 // ----------------------------------------------------------------------------
7 
8 #pragma once
9 #include "Common.h"
10 #include "CutPursuit.h"
11 
12 namespace CP {
13 template <typename T>
14 struct CutPursuit_KL : public CutPursuit<T> {
15  std::pair<T, T> compute_energy() override {
16  VertexAttributeMap<T> vertex_attribute_map =
17  boost::get(boost::vertex_bundle, this->main_graph);
18  EdgeAttributeMap<T> edge_attribute_map =
19  boost::get(boost::edge_bundle, this->main_graph);
20  std::pair<T, T> pair_energy;
21  T energy = 0, smoothedObservation, smoothedValue;
22  // #pragma omp parallel if (this->parameter.parallel)
23  for (VertexIterator<T> i_ver = boost::vertices(this->main_graph).first;
24  i_ver != this->lastIterator; ++i_ver) {
25  for (uint32_t i_dim = 0; i_dim < this->dim;
26  i_dim++) { // smoothing as a linear combination with the
27  // uniform probability
28  smoothedObservation =
29  this->parameter.smoothing / this->dim +
30  (1 - this->parameter.smoothing) *
31  vertex_attribute_map(*i_ver).observation[i_dim];
32  smoothedValue =
33  this->parameter.smoothing / this->dim +
34  (1 - this->parameter.smoothing) *
35  vertex_attribute_map(*i_ver).value[i_dim];
36  energy += smoothedObservation *
37  (log(smoothedObservation) - log(smoothedValue)) *
38  vertex_attribute_map(*i_ver).weight;
39  }
40  }
41  pair_energy.first = energy;
42  energy = 0;
43  EdgeIterator<T> i_edg_end = boost::edges(this->main_graph).second;
44  for (EdgeIterator<T> i_edg = boost::edges(this->main_graph).first;
45  i_edg != i_edg_end; ++i_edg) {
46  if (!edge_attribute_map(*i_edg).realEdge) {
47  continue;
48  }
49  energy += .5 * edge_attribute_map(*i_edg).isActive *
50  this->parameter.reg_strenth *
51  edge_attribute_map(*i_edg).weight;
52  }
53  pair_energy.second = energy;
54  return pair_energy;
55  }
56 
57  //=============================================================================================
58  //============================= SPLIT
59  //===========================================
60  //=============================================================================================
61  size_t split()
62  override { // split the graph by trying to find the best binary
63  // partition each components is split into B and notB
64  // for each components we associate the value h_1 and
65  // h_2 to vertices in B or notB the affectation as well
66  // as h_1 and h_2 are computed alternatively
67  // tic();
68  //--------loading
69  // structures---------------------------------------------------------------
70  TimeStack ts;
71  ts.tic();
72  uint32_t nb_comp = static_cast<uint32_t>(this->components.size());
73  VertexAttributeMap<T> vertex_attribute_map =
74  boost::get(boost::vertex_bundle, this->main_graph);
75  VertexIndexMap<T> vertex_index_map =
76  boost::get(boost::vertex_index, this->main_graph);
77  // initialize h_1 and h_2 with kmeans
78  // stores wether each vertex is B or notB
79  std::vector<bool> binary_label(this->nVertex);
80  this->init_labels(binary_label);
81  VectorOfCentroids<T> centers(nb_comp, this->dim);
82  //-----main
83  // loop----------------------------------------------------------------
84  // the optimal flow is iteratively approximated
85  for (uint32_t i_step = 1; i_step <= this->parameter.flow_steps;
86  i_step++) {
87  // compute h_1 and h_2
88  centers = VectorOfCentroids<T>(nb_comp, this->dim);
89  this->compute_centers(centers, binary_label);
90  // update the capacities of the flow graph
91  this->set_capacities(centers);
92  // compute flow
93  boost::boykov_kolmogorov_max_flow(
94  this->main_graph,
95  get(&EdgeAttribute<T>::capacity, this->main_graph),
96  get(&EdgeAttribute<T>::residualCapacity, this->main_graph),
97  get(&EdgeAttribute<T>::edge_reverse, this->main_graph),
98  get(&VertexAttribute<T>::color, this->main_graph),
99  get(boost::vertex_index, this->main_graph), this->source,
100  this->sink);
101 
102  for (uint32_t i_com = 0; i_com < nb_comp; i_com++) {
103  if (this->saturated_components[i_com]) {
104  continue;
105  }
106  for (uint32_t i_ver = 0; i_ver < this->components[i_com].size();
107  i_ver++) {
108  binary_label[vertex_index_map(
109  this->components[i_com][i_ver])] =
110  (vertex_attribute_map(
111  this->components[i_com][i_ver])
112  .color ==
113  vertex_attribute_map(this->sink).color);
114  }
115  }
116  }
117  size_t saturation = this->activate_edges();
118  return saturation;
119  }
120 
121  //=============================================================================================
122  //============================= INIT_KL
123  //===================================================
124  //=============================================================================================
125  inline void init_labels(
126  std::vector<bool>&
127  binary_label) { //-----initialize the labelling for each
128  // components with
129  // kmeans------------------------------
130  VertexAttributeMap<T> vertex_attribute_map =
131  boost::get(boost::vertex_bundle, this->main_graph);
132  VertexIndexMap<T> vertex_index_map =
133  boost::get(boost::vertex_index, this->main_graph);
134  std::vector<std::vector<T>> kernels(2, std::vector<T>(this->dim));
135  std::vector<std::vector<T>> smooth_kernels(2,
136  std::vector<T>(this->dim));
137  T total_weight[2];
138  uint32_t nb_comp = static_cast<uint32_t>(this->components.size());
139  T best_energy, current_energy;
140  // #pragma omp parallel for private(kernels, total_weight, best_energy,
141  // current_energy) if (this->parameter.parallel && nb_comp>8)
142  // schedule(dynamic)
143  for (uint32_t i_com = 0; i_com < nb_comp; i_com++) {
144  uint32_t comp_size =
145  static_cast<uint32_t>(this->components[i_com].size());
146  std::vector<bool> potential_label(comp_size);
147  std::vector<T> energy_array(comp_size);
148  std::vector<T> constant_part(comp_size);
149  std::vector<std::vector<T>> smooth_obs(comp_size,
150  std::vector<T>(2, 0));
151  if (this->saturated_components[i_com] || comp_size <= 1) {
152  continue;
153  }
154  // KL fidelity has a part that depends
155  // purely on the observation that can be precomputed
156  // #pragma omp parallel for if (this->parameter.parallel &&
157  // nb_comp<=8) schedule(dynamic)
158  for (uint32_t i_ver = 0; i_ver < comp_size; i_ver++) {
159  constant_part[i_ver] = 0;
160  for (uint32_t i_dim = 0; i_dim < this->dim; i_dim++) {
161  smooth_obs[i_ver][i_dim] = 0;
162  }
163  for (uint32_t i_dim = 0; i_dim < this->dim; i_dim++) {
164  smooth_obs[i_ver][i_dim] =
165  this->parameter.smoothing / this->dim +
166  (1 - this->parameter.smoothing) *
167  vertex_attribute_map(
168  this->components[i_com][i_ver])
169  .observation[i_dim];
170  constant_part[i_ver] +=
171  smooth_obs[i_ver][i_dim] *
172  log(smooth_obs[i_ver][i_dim]) *
173  vertex_attribute_map(this->components[i_com][i_ver])
174  .weight;
175  }
176  }
177  for (uint32_t init_kmeans = 0;
178  init_kmeans < this->parameter.kmeans_resampling;
179  init_kmeans++) {
180  //----- initialization with KM++ ------------------
181  // first kernel chosen randomly
182  uint32_t first_kernel = std::rand() % comp_size,
183  second_kernel = 0;
184  for (uint32_t i_dim = 0; i_dim < this->dim;
185  i_dim++) { // fill the first kernel
186  kernels[0][i_dim] =
187  vertex_attribute_map(
188  this->components[i_com][first_kernel])
189  .observation[i_dim];
190  smooth_kernels[0][i_dim] =
191  this->parameter.smoothing / this->dim +
192  (1 - this->parameter.smoothing) * kernels[0][i_dim];
193  }
194  // now compute the square distance of each pouint32_t to this
195  // kernel
196  best_energy = 0; // energy total
197  // #pragma omp parallel for if (this->parameter.parallel &&
198  // nb_comp<=8) schedule(dynamic)
199  for (uint32_t i_ver = 0; i_ver < comp_size; i_ver++) {
200  energy_array[i_ver] = constant_part[i_ver];
201  for (uint32_t i_dim = 0; i_dim < this->dim; i_dim++) {
202  energy_array[i_ver] -=
203  smooth_obs[i_ver][i_dim] *
204  log(smooth_kernels[0][i_dim]) *
205  vertex_attribute_map(
206  this->components[i_com][i_ver])
207  .weight;
208  }
209  energy_array[i_ver] = pow(energy_array[i_ver], 2);
210  best_energy += energy_array[i_ver];
211  } // we now generate a random number to determinate which node
212  // will be the second kernel
213  if (best_energy ==
214  0) { // all the points in this components are identical
215  for (uint32_t i_ver = 0; i_ver < comp_size; i_ver++) {
216  binary_label[vertex_index_map(
217  this->components[i_com][i_ver])] = false;
218  }
219  break;
220  }
221  // we now choose the second kernel with a probability
222  // proportional to the square distance
223  T random_sample = ((T)(rand())) / ((T)(RAND_MAX));
224  current_energy = best_energy * random_sample;
225  for (uint32_t i_ver = 0; i_ver < comp_size; i_ver++) {
226  current_energy -= energy_array[i_ver];
227  if (current_energy <
228  0) { // we have selected the second kernel
229  second_kernel = i_ver;
230  break;
231  }
232  }
233  for (uint32_t i_dim = 0; i_dim < this->dim;
234  i_dim++) { // now fill the second kernel
235  kernels[1][i_dim] =
236  vertex_attribute_map(
237  this->components[i_com][second_kernel])
238  .observation[i_dim];
239  smooth_kernels[1][i_dim] =
240  this->parameter.smoothing / this->dim +
241  (1 - this->parameter.smoothing) * kernels[1][i_dim];
242  }
243  //----main kmeans loop-----
244  for (uint32_t ite_kmeans = 0;
245  ite_kmeans < this->parameter.kmeans_ite; ite_kmeans++) {
246  //--affectation step: associate each node with its closest
247  // kernel------------------- #pragma omp parallel for if
248  //(this->parameter.parallel && nb_comp<=8) schedule(dynamic)
249  for (uint32_t i_ver = 0; i_ver < comp_size; i_ver++) {
250  // the distance to each kernel
251  std::vector<T> distance_kernels(2,
252  constant_part[i_ver]);
253  for (uint32_t i_dim = 0; i_dim < this->dim; i_dim++) {
254  distance_kernels[0] -=
255  smooth_obs[i_ver][i_dim] *
256  log(smooth_kernels[0][i_dim]) *
257  vertex_attribute_map(
258  this->components[i_com][i_ver])
259  .weight;
260  distance_kernels[1] -=
261  smooth_obs[i_ver][i_dim] *
262  log(smooth_kernels[1][i_dim]) *
263  vertex_attribute_map(
264  this->components[i_com][i_ver])
265  .weight;
266  }
267  potential_label[i_ver] =
268  distance_kernels[0] > distance_kernels[1];
269  }
270  //-----computation of the new
271  // kernels----------------------------
272  total_weight[0] = 0.;
273  total_weight[1] = 0.;
274  for (uint32_t i_dim = 0; i_dim < this->dim; i_dim++) {
275  kernels[0][i_dim] = 0;
276  kernels[1][i_dim] = 0;
277  }
278  for (uint32_t i_ver = 0; i_ver < comp_size; i_ver++) {
279  if (vertex_attribute_map(this->components[i_com][i_ver])
280  .weight == 0) {
281  continue;
282  }
283  if (potential_label[i_ver]) {
284  total_weight[0] +=
285  vertex_attribute_map(
286  this->components[i_com][i_ver])
287  .weight;
288  for (uint32_t i_dim = 0; i_dim < this->dim;
289  i_dim++) {
290  kernels[0][i_dim] +=
291  vertex_attribute_map(
292  this->components[i_com][i_ver])
293  .observation[i_dim] *
294  vertex_attribute_map(
295  this->components[i_com][i_ver])
296  .weight;
297  }
298  } else {
299  total_weight[1] +=
300  vertex_attribute_map(
301  this->components[i_com][i_ver])
302  .weight;
303  for (uint32_t i_dim = 0; i_dim < this->dim;
304  i_dim++) {
305  kernels[1][i_dim] +=
306  vertex_attribute_map(
307  this->components[i_com][i_ver])
308  .observation[i_dim] *
309  vertex_attribute_map(
310  this->components[i_com][i_ver])
311  .weight;
312  }
313  }
314  }
315  if ((total_weight[0] == 0) || (total_weight[1] == 0)) {
316  std::cout << "kmeans error" << std::endl;
317  }
318  for (uint32_t i_dim = 0; i_dim < this->dim; i_dim++) {
319  kernels[0][i_dim] = kernels[0][i_dim] / total_weight[0];
320  kernels[1][i_dim] = kernels[1][i_dim] / total_weight[1];
321  smooth_kernels[0][i_dim] =
322  this->parameter.smoothing / this->dim +
323  (1 - this->parameter.smoothing) *
324  kernels[0][i_dim];
325  smooth_kernels[1][i_dim] =
326  this->parameter.smoothing / this->dim +
327  (1 - this->parameter.smoothing) *
328  kernels[1][i_dim];
329  }
330  }
331  //----compute the associated energy ------
332  current_energy = 0;
333  for (uint32_t i_ver = 0; i_ver < comp_size; i_ver++) {
334  current_energy += constant_part[i_ver];
335  if (potential_label[i_ver]) {
336  for (uint32_t i_dim = 0; i_dim < this->dim; i_dim++) {
337  current_energy -=
338  smooth_obs[i_ver][i_dim] *
339  log(smooth_kernels[0][i_dim]) *
340  vertex_attribute_map(
341  this->components[i_com][i_ver])
342  .weight;
343  }
344  } else {
345  for (uint32_t i_dim = 0; i_dim < this->dim; i_dim++) {
346  current_energy -=
347  smooth_obs[i_ver][i_dim] *
348  log(smooth_kernels[1][i_dim]) *
349  vertex_attribute_map(
350  this->components[i_com][i_ver])
351  .weight;
352  }
353  }
354  }
355  if (current_energy < best_energy) {
356  best_energy = current_energy;
357  for (uint32_t i_ver = 0; i_ver < comp_size; i_ver++) {
358  binary_label[vertex_index_map(
359  this->components[i_com][i_ver])] =
360  potential_label[i_ver];
361  }
362  }
363  }
364  }
365  }
366 
367  //=============================================================================================
368  //============================= COMPUTE_CENTERS_KL
369  //==========================================
370  //=============================================================================================
371  inline void compute_centers(VectorOfCentroids<T>& centers,
372  const std::vector<bool>& binary_label) {
373  // compute for each component the values of h_1 and h_2
374  VertexAttributeMap<T> vertex_attribute_map =
375  boost::get(boost::vertex_bundle, this->main_graph);
376  VertexIndexMap<T> vertex_index_map =
377  boost::get(boost::vertex_index, this->main_graph);
378  uint32_t nb_comp = static_cast<uint32_t>(this->components.size());
379  // #pragma omp parallel for if (this->parameter.parallel)
380  // schedule(dynamic)
381  for (uint32_t i_com = 0; i_com < nb_comp; i_com++) {
382  if (this->saturated_components[i_com]) {
383  continue;
384  }
385  T total_weight[2];
386  total_weight[0] = 0.;
387  total_weight[1] = 0.;
388  for (uint32_t i_dim = 0; i_dim < this->dim; i_dim++) {
389  centers.centroids[i_com][0][i_dim] = 0.;
390  centers.centroids[i_com][1][i_dim] = 0.;
391  }
392  for (uint32_t i_ver = 0; i_ver < this->components[i_com].size();
393  i_ver++) {
394  if (vertex_attribute_map(this->components[i_com][i_ver])
395  .weight == 0) {
396  continue;
397  }
398  if (binary_label[vertex_index_map(
399  this->components[i_com][i_ver])]) {
400  total_weight[0] +=
401  vertex_attribute_map(this->components[i_com][i_ver])
402  .weight;
403  for (uint32_t i_dim = 0; i_dim < this->dim; i_dim++) {
404  centers.centroids[i_com][0][i_dim] +=
405  vertex_attribute_map(
406  this->components[i_com][i_ver])
407  .observation[i_dim] *
408  vertex_attribute_map(
409  this->components[i_com][i_ver])
410  .weight;
411  }
412  } else {
413  total_weight[1] +=
414  vertex_attribute_map(this->components[i_com][i_ver])
415  .weight;
416  for (uint32_t i_dim = 0; i_dim < this->dim; i_dim++) {
417  centers.centroids[i_com][1][i_dim] +=
418  vertex_attribute_map(
419  this->components[i_com][i_ver])
420  .observation[i_dim] *
421  vertex_attribute_map(
422  this->components[i_com][i_ver])
423  .weight;
424  }
425  }
426  }
427  if ((total_weight[0] == 0) || (total_weight[1] == 0)) {
428  // the component is saturated
429  this->saturateComponent(i_com);
430  for (uint32_t i_dim = 0; i_dim < this->dim; i_dim++) {
431  centers.centroids[i_com][0][i_dim] =
432  vertex_attribute_map(this->components[i_com].back())
433  .value[i_dim];
434  centers.centroids[i_com][1][i_dim] =
435  vertex_attribute_map(this->components[i_com].back())
436  .value[i_dim];
437  }
438  } else {
439  for (uint32_t i_dim = 0; i_dim < this->dim; i_dim++) {
440  centers.centroids[i_com][0][i_dim] =
441  centers.centroids[i_com][0][i_dim] /
442  total_weight[0];
443  centers.centroids[i_com][1][i_dim] =
444  centers.centroids[i_com][1][i_dim] /
445  total_weight[1];
446  }
447  }
448  }
449  }
450 
451  //=============================================================================================
452  //============================= SET_CAPACITIES
453  //==========================================
454  //=============================================================================================
455  inline void set_capacities(const VectorOfCentroids<T>& centers) {
456  VertexAttributeMap<T> vertex_attribute_map =
457  boost::get(boost::vertex_bundle, this->main_graph);
458  EdgeAttributeMap<T> edge_attribute_map =
459  boost::get(boost::edge_bundle, this->main_graph);
460  VertexDescriptor<T> desc_v;
461  EdgeDescriptor desc_source2v, desc_v2sink, desc_v2source;
462  uint32_t nb_comp = static_cast<uint32_t>(this->components.size());
463  T cost_B, cost_notB, smoothedValueB, smoothedValueNotB,
464  smoothedObservation; // the cost of being in B or not B, local
465  // for each component
466  //----first compute the capacity in sink/node
467  // edges------------------------------------ #pragma omp parallel for if
468  //(this->parameter.parallel) schedule(dynamic)
469  for (uint32_t i_com = 0; i_com < nb_comp; i_com++) {
470  if (this->saturated_components[i_com]) {
471  continue;
472  }
473  for (uint32_t i_ver = 0; i_ver < this->components[i_com].size();
474  i_ver++) {
475  desc_v = this->components[i_com][i_ver];
476  // because of the adjacency structure NEVER access edge
477  // (source,v) directly!
478  desc_v2source =
479  boost::edge(desc_v, this->source, this->main_graph)
480  .first;
481  desc_source2v =
482  edge_attribute_map(desc_v2source)
483  .edge_reverse; // use edge_reverse instead
484  desc_v2sink =
485  boost::edge(desc_v, this->sink, this->main_graph).first;
486  cost_B = 0;
487  cost_notB = 0;
488  if (vertex_attribute_map(desc_v).weight == 0) {
489  edge_attribute_map(desc_source2v).capacity = 0;
490  edge_attribute_map(desc_v2sink).capacity = 0;
491  continue;
492  }
493  for (uint32_t i_dim = 0; i_dim < this->dim; i_dim++) {
494  smoothedObservation =
495  this->parameter.smoothing / this->dim +
496  (1 - this->parameter.smoothing) *
497  vertex_attribute_map(desc_v)
498  .observation[i_dim];
499  smoothedValueB = this->parameter.smoothing / this->dim +
500  (1 - this->parameter.smoothing) *
501  centers.centroids[i_com][0][i_dim];
502  smoothedValueNotB =
503  this->parameter.smoothing / this->dim +
504  (1 - this->parameter.smoothing) *
505  centers.centroids[i_com][1][i_dim];
506  cost_B += smoothedObservation *
507  (log(smoothedObservation) - log(smoothedValueB));
508  cost_notB +=
509  smoothedObservation *
510  (log(smoothedObservation) - log(smoothedValueNotB));
511  }
512  if (cost_B > cost_notB) {
513  edge_attribute_map(desc_source2v).capacity =
514  cost_B - cost_notB;
515  edge_attribute_map(desc_v2sink).capacity = 0.;
516  } else {
517  edge_attribute_map(desc_source2v).capacity = 0.;
518  edge_attribute_map(desc_v2sink).capacity =
519  cost_notB - cost_B;
520  }
521  }
522  }
523  //----then set the vertex to vertex edges
524  //---------------------------------------------
525  EdgeIterator<T> i_edg, i_edg_end;
526  for (boost::tie(i_edg, i_edg_end) = boost::edges(this->main_graph);
527  i_edg != i_edg_end; ++i_edg) {
528  if (!edge_attribute_map(*i_edg).realEdge) {
529  continue;
530  }
531  if (!edge_attribute_map(*i_edg).isActive) {
532  edge_attribute_map(*i_edg).capacity =
533  edge_attribute_map(*i_edg).weight *
534  this->parameter.reg_strenth;
535  } else {
536  edge_attribute_map(*i_edg).capacity = 0;
537  }
538  }
539  }
540 
541  //=============================================================================================
542  //================================= COMPUTE_VALUE
543  //=========================================
544  //=============================================================================================
545  std::pair<std::vector<T>, T> compute_value(const uint32_t& i_com) override {
546  VertexAttributeMap<T> vertex_attribute_map =
547  boost::get(boost::vertex_bundle, this->main_graph);
548  T total_weight = 0;
549  std::vector<T> compValue(this->dim);
550  std::fill((compValue.begin()), (compValue.end()), 0);
551  for (uint32_t ind_ver = 0; ind_ver < this->components[i_com].size();
552  ++ind_ver) {
553  total_weight +=
554  vertex_attribute_map(this->components[i_com][ind_ver])
555  .weight;
556  for (uint32_t i_dim = 0; i_dim < this->dim; i_dim++) {
557  compValue[i_dim] +=
558  vertex_attribute_map(this->components[i_com][ind_ver])
559  .observation[i_dim] *
560  vertex_attribute_map(this->components[i_com][ind_ver])
561  .weight;
562  }
563  vertex_attribute_map(this->components[i_com][ind_ver])
564  .in_component = i_com;
565  }
566  for (uint32_t i_dim = 0; i_dim < this->dim; i_dim++) {
567  compValue[i_dim] = compValue[i_dim] / total_weight;
568  }
569  for (uint32_t ind_ver = 0; ind_ver < this->components[i_com].size();
570  ++ind_ver) {
571  for (uint32_t i_dim = 0; i_dim < this->dim; i_dim++) {
572  vertex_attribute_map(this->components[i_com][ind_ver])
573  .value[i_dim] = compValue[i_dim];
574  }
575  }
576  return std::pair<std::vector<T>, T>(compValue, total_weight);
577  }
578 
579  //=============================================================================================
580  //================================= COMPUTE_MERGE_GAIN
581  //=========================================
582  //=============================================================================================
583  std::pair<std::vector<T>, T> compute_merge_gain(
584  const VertexDescriptor<T>& comp1,
585  const VertexDescriptor<T>& comp2) override {
586  VertexAttributeMap<T> reduced_vertex_attribute_map =
587  boost::get(boost::vertex_bundle, this->reduced_graph);
588  std::vector<T> merge_value(this->dim);
589  T gain = 0, smoothedValue1, smoothedValue2, smoothedValueMerged;
590  // compute the value obtained by mergeing the two connected components
591  for (uint32_t i_dim = 0; i_dim < this->dim; i_dim++) {
592  merge_value[i_dim] =
593  (reduced_vertex_attribute_map(comp1).weight *
594  reduced_vertex_attribute_map(comp1).value[i_dim] +
595  reduced_vertex_attribute_map(comp2).weight *
596  reduced_vertex_attribute_map(comp2).value[i_dim]) /
597  (reduced_vertex_attribute_map(comp1).weight +
598  reduced_vertex_attribute_map(comp2).weight);
599  smoothedValue1 =
600  this->parameter.smoothing / this->dim +
601  (1 - this->parameter.smoothing) *
602  reduced_vertex_attribute_map(comp1).value[i_dim];
603  smoothedValue2 =
604  this->parameter.smoothing / this->dim +
605  (1 - this->parameter.smoothing) *
606  reduced_vertex_attribute_map(comp2).value[i_dim];
607  smoothedValueMerged =
608  this->parameter.smoothing / this->dim +
609  (1 - this->parameter.smoothing) * merge_value[i_dim];
610  gain -= reduced_vertex_attribute_map(comp1).weight *
611  smoothedValue1 *
612  (log(smoothedValue1) - log(smoothedValueMerged)) +
613  reduced_vertex_attribute_map(comp2).weight *
614  smoothedValue2 *
615  (log(smoothedValue2) - log(smoothedValueMerged));
616  }
617  return std::pair<std::vector<T>, T>(merge_value, gain);
618  }
619 };
620 } // namespace CP
vertex_index
long vertex_index
Definition: FileIOFactory.cpp:47
TimeStack::tic
void tic()
Definition: Common.h:57
VectorOfCentroids::centroids
std::vector< std::vector< std::vector< T > > > centroids
Definition: Common.h:109
CP
Definition: API.h:123
CP::VertexAttributeMap
boost::vec_adj_list_vertex_property_map< Graph< T >, Graph< T > *, VertexAttribute< T >, VertexAttribute< T > &, boost::vertex_bundle_t > VertexAttributeMap
Definition: Graph.h:92
CP::EdgeIterator
typename boost::graph_traits< Graph< T > >::edge_iterator EdgeIterator
Definition: Graph.h:84
CP::VertexIterator
typename boost::graph_traits< Graph< T > >::vertex_iterator VertexIterator
Definition: Graph.h:82
CP::EdgeDescriptor
boost::graph_traits< boost::adjacency_list< boost::vecS, boost::vecS, boost::directedS > >::edge_descriptor EdgeDescriptor
Definition: Graph.h:22
CP::VertexDescriptor
typename boost::graph_traits< CP::Graph< T > >::vertex_descriptor VertexDescriptor
Definition: Graph.h:75
CP::VertexIndexMap
typename boost::property_map< Graph< T >, boost::vertex_index_t >::type VertexIndexMap
Definition: Graph.h:103
CP::EdgeAttributeMap
boost::adj_list_edge_property_map< boost::directed_tag, EdgeAttribute< T >, EdgeAttribute< T > &, uint64_t, CP::EdgeAttribute< T >, boost::edge_bundle_t > EdgeAttributeMap
Definition: Graph.h:100
QtCompat::endl
QTextStream & endl(QTextStream &stream)
Definition: QtCompat.h:718
flann::random_sample
Matrix< T > random_sample(Matrix< T > &srcMatrix, size_t size, bool remove=false)
Definition: sampling.h:40
CP::CutPursuit_KL
Definition: CutPursuit_KL.h:14
CP::CutPursuit_KL::init_labels
void init_labels(std::vector< bool > &binary_label)
Definition: CutPursuit_KL.h:125
CP::CutPursuit_KL::set_capacities
void set_capacities(const VectorOfCentroids< T > &centers)
Definition: CutPursuit_KL.h:455
CP::CutPursuit_KL::compute_value
std::pair< std::vector< T >, T > compute_value(const uint32_t &i_com) override
Definition: CutPursuit_KL.h:545
CP::CutPursuit_KL::compute_energy
std::pair< T, T > compute_energy() override
Definition: CutPursuit_KL.h:15
CP::CutPursuit_KL::compute_merge_gain
std::pair< std::vector< T >, T > compute_merge_gain(const VertexDescriptor< T > &comp1, const VertexDescriptor< T > &comp2) override
Definition: CutPursuit_KL.h:583
CP::CutPursuit_KL::compute_centers
void compute_centers(VectorOfCentroids< T > &centers, const std::vector< bool > &binary_label)
Definition: CutPursuit_KL.h:371
CP::CutPursuit_KL::split
size_t split() override
Definition: CutPursuit_KL.h:61
CP::CutPursuit::parameter
CPparameter< T > parameter
Definition: CutPursuit.h:65
CP::CutPursuit::main_graph
Graph< T > main_graph
Definition: CutPursuit.h:45
CP::CutPursuit::lastIterator
CP::VertexIterator< T > lastIterator
Definition: CutPursuit.h:63
CP::CutPursuit::components
std::vector< std::vector< VertexDescriptor< T > > > components
Definition: CutPursuit.h:49
CP::CutPursuit::sink
VertexDescriptor< T > sink
Definition: CutPursuit.h:58
CP::CutPursuit::source
VertexDescriptor< T > source
Definition: CutPursuit.h:57
CP::CutPursuit::activate_edges
size_t activate_edges(bool allows_saturation=true)
Definition: CutPursuit.h:298
CP::CutPursuit::dim
uint32_t dim
Definition: CutPursuit.h:59
CP::CutPursuit::nVertex
uint32_t nVertex
Definition: CutPursuit.h:60
CP::CutPursuit::saturated_components
std::vector< bool > saturated_components
Definition: CutPursuit.h:53
CP::CutPursuit::reduced_graph
Graph< T > reduced_graph
Definition: CutPursuit.h:46
CP::CutPursuit::saturateComponent
void saturateComponent(const uint32_t &ind_com)
Definition: CutPursuit.h:835