add submodule enercast
[epclust.git] / epclust / src / computeMedoidsIndices.cpp
1 #include <Rcpp.h>
2
3 // [[Rcpp::depends(BH, bigmemory)]]
4 #include <bigmemory/MatrixAccessor.hpp>
5
6 #include <numeric>
7 #include <cfloat>
8
9 using namespace Rcpp;
10
11 //' computeMedoidsIndices
12 //'
13 //' Compute medoids indices
14 //'
15 //' @param pMedoids External pointer
16 //' @param ref_series reference series
17 //'
18 //' @return A map serie number -> medoid index
19 // [[Rcpp::export]]
20 IntegerVector computeMedoidsIndices(SEXP pMedoids, NumericMatrix ref_series)
21 {
22 XPtr<BigMatrix> pMed(pMedoids);
23 MatrixAccessor<double> medoids = MatrixAccessor<double>(*pMed);
24 int nb_series = ref_series.nrow(),
25 K = pMed->nrow(),
26 L = pMed->ncol();
27 IntegerVector mi(nb_series);
28
29 for (int i=0; i<nb_series ; i++)
30 {
31 // mi[i] <- which.min( rowSums( sweep(medoids, 2, ref_series[i,], '-')^2 ) )
32 mi[i] = 0;
33 double best_dist = DBL_MAX;
34 for (int j=0; j<K; j++)
35 {
36 double dist_ij = 0.;
37 for (int k=0; k<L; k++)
38 {
39 double delta = ref_series(i,k) - *(medoids[k]+j);
40 dist_ij += delta * delta;
41 }
42 if (dist_ij < best_dist)
43 {
44 mi[i] = j+1; //R indices start at 1
45 best_dist = dist_ij;
46 }
47 }
48 }
49
50 return mi;
51 }