add comments, fix some things. TODO: comment tests, finish computeWerDists, test it

[epclust.git] / epclust / src / computeMedoidsIndices.cpp

diff --git a/epclust/src/computeMedoidsIndices.cpp b/epclust/src/computeMedoidsIndices.cpp
index f247584..895031a 100644 (file)
--- a/epclust/src/computeMedoidsIndices.cpp
+++ b/epclust/src/computeMedoidsIndices.cpp
@@ -10,33 +10,39 @@ using namespace Rcpp;
 
 //' computeMedoidsIndices
 //'
-//' Compute medoids indices
+//' For each column of the 'series' matrix input, search for the closest medoid
+//' (euclidian distance) and store its index
 //'
-//' @param pMedoids External pointer
-//' @param ref_series reference series
+//' @param pMedoids External pointer (a big.matrix 'address' slot in R)
+//' @param series (reference) series, a matrix of size Lxn
 //'
-//' @return A map serie number -> medoid index
+//' @return An integer vector of the closest medoids indices, for each (column) serie
 // [[Rcpp::export]]
-IntegerVector computeMedoidsIndices(SEXP pMedoids, NumericMatrix ref_series)
+IntegerVector computeMedoidsIndices(SEXP pMedoids, NumericMatrix series)
 {
+       // Turn SEXP external pointer into BigMatrix (description) object
        XPtr<BigMatrix> pMed(pMedoids);
+       // medoids: access to the content of the BigMatrix object
        MatrixAccessor<double> medoids = MatrixAccessor<double>(*pMed);
-       int nb_series = ref_series.ncol(),
+
+       int nb_series = series.ncol(),
                K = pMed->ncol(),
                L = pMed->nrow();
        IntegerVector mi(nb_series);
 
-       for (int i=0; i<nb_series ; i++)
+       for (int i=0; i<nb_series ; i++) //column index in series
        {
-               // mi[i] <- which.min( rowSums( sweep(medoids, 2, ref_series[i,], '-')^2 ) )
+               // In R: mi[i] <- which.min( rowSums( sweep(medoids, 2, series[i,], '-')^2 ) )
+               // In C(++), computations must be unrolled
                mi[i] = 0;
                double best_dist = DBL_MAX;
-               for (int j=0; j<K; j++)
+               for (int j=0; j<K; j++) //column index in medoids
                {
                        double dist_ij = 0.;
-                       for (int k=0; k<L; k++)
+                       for (int k=0; k<L; k++) //common row index (medoids, series)
                        {
-                               double delta = ref_series(k,i) - *(medoids[j]+k);
+                               // Accessing values for a big matrix is a bit weird; see Rcpp gallery/bigmemory
+                               double delta = series(k,i) - *(medoids[j]+k);
                                dist_ij += delta * delta;
                        }
                        if (dist_ij < best_dist)