Move ComputeP from gmmreg_utils to utils/em_utils

bing-jian · Jun 14, 2019 · 0f22c7e · 0f22c7e
1 parent c50b51d
commit 0f22c7e
Show file tree

Hide file tree

Showing 4 changed files with 65 additions and 51 deletions.
diff --git a/C++/gmmreg_cpd.cpp b/C++/gmmreg_cpd.cpp
@@ -1,15 +1,11 @@
 #include "gmmreg_cpd.h"
 
-#include <cstdlib>
-#include <cstring>
-#include <fstream>
-#include <iostream>
-#include <vnl/algo/vnl_determinant.h>
 #include <vnl/algo/vnl_qr.h>
 #include <vnl/algo/vnl_svd.h>
 #include <vnl/vnl_trace.h>
 
 #include "gmmreg_utils.h"
+#include "utils/em_utils.h"
 #include "utils/io_utils.h"
 
 namespace gmmreg {

diff --git a/C++/gmmreg_utils.h b/C++/gmmreg_utils.h
@@ -33,9 +33,6 @@ template <typename T>
 void Denormalize(vnl_matrix<T>& x, const vnl_vector<T>& centroid,
                  const T scale);
 
-template <typename T>
-void ComputeP(const vnl_matrix<T>& x, const vnl_matrix<T>& y, vnl_matrix<T>& P,
-              T& E, T sigma, int outliers);
 
 #define SQR(X) ((X) * (X))
 
@@ -273,49 +270,6 @@ void Denormalize(vnl_matrix<T>& x, const vnl_vector<T>& centroid,
   }
 }
 
-template <typename T>
-void ComputeP(const vnl_matrix<T>& x, const vnl_matrix<T>& y, vnl_matrix<T>& P,
-              T& E, T sigma, int outliers) {
-  T k;
-  k = -2 * sigma * sigma;
-
-  vnl_vector<T> column_sum;
-  int m = x.rows();
-  int s = y.rows();
-  int d = x.cols();
-  column_sum.set_size(s);
-  column_sum.fill(0);
-  T outlier_term = outliers * pow((2 * sigma * sigma * 3.1415926), 0.5 * d);
-#pragma omp for
-  for (int i = 0; i < m; ++i) {
-    for (int j = 0; j < s; ++j) {
-      T r = 0;
-      for (int t = 0; t < d; ++t) {
-        r += (x(i, t) - y(j, t)) * (x(i, t) - y(j, t));
-      }
-      P(i, j) = exp(r / k);
-      column_sum[j] += P(i, j);
-    }
-  }
-
-  if (outliers != 0) {
-#pragma omp for
-    for (int i = 0; i < s; ++i) column_sum[i] += outlier_term;
-  }
-  if (column_sum.min_value() > (1e-12)) {
-    E = 0;
-#pragma omp for
-    for (int i = 0; i < s; ++i) {
-      for (int j = 0; j < m; ++j) {
-        P(j, i) = P(j, i) / column_sum[i];
-      }
-      E -= log(column_sum[i]);
-    }
-  } else {
-    P.empty();
-  }
-}
-
 }  // namespace gmmreg
 
 #endif  // GMMREG_UTILS_H_
diff --git a/C++/utils/em_utils.cc b/C++/utils/em_utils.cc
@@ -0,0 +1,48 @@
+#include "em_utils.h"
+
+namespace gmmreg {
+
+template <typename T>
+void ComputeP(const vnl_matrix<T>& x, const vnl_matrix<T>& y, vnl_matrix<T>& P,
+              T& E, T sigma, int outliers) {
+  T k;
+  k = -2 * sigma * sigma;
+
+  vnl_vector<T> column_sum;
+  int m = x.rows();
+  int s = y.rows();
+  int d = x.cols();
+  column_sum.set_size(s);
+  column_sum.fill(0);
+  T outlier_term = outliers * pow((2 * sigma * sigma * 3.1415926), 0.5 * d);
+#pragma omp for
+  for (int i = 0; i < m; ++i) {
+    for (int j = 0; j < s; ++j) {
+      T r = 0;
+      for (int t = 0; t < d; ++t) {
+        r += (x(i, t) - y(j, t)) * (x(i, t) - y(j, t));
+      }
+      P(i, j) = exp(r / k);
+      column_sum[j] += P(i, j);
+    }
+  }
+
+  if (outliers != 0) {
+#pragma omp for
+    for (int i = 0; i < s; ++i) column_sum[i] += outlier_term;
+  }
+  if (column_sum.min_value() > (1e-12)) {
+    E = 0;
+#pragma omp for
+    for (int i = 0; i < s; ++i) {
+      for (int j = 0; j < m; ++j) {
+        P(j, i) = P(j, i) / column_sum[i];
+      }
+      E -= log(column_sum[i]);
+    }
+  } else {
+    P.empty();
+  }
+}
+
+}  // namespace gmmreg
diff --git a/C++/utils/em_utils.h b/C++/utils/em_utils.h
@@ -0,0 +1,16 @@
+#ifndef GMMREG_UTILS_EM_UTILS_H_
+#define GMMREG_UTILS_EM_UTILS_H_
+
+#include <vnl/vnl_matrix.h>
+
+namespace gmmreg {
+
+template <typename T>
+void ComputeP(const vnl_matrix<T>& x, const vnl_matrix<T>& y, vnl_matrix<T>& P,
+              T& E, T sigma, int outliers);
+
+}  // namespace gmmreg
+
+#include "em_utils.cc"
+
+#endif // GMMREG_UTILS_EM_UTILS_H_