X-Git-Url: https://git.auder.net/?a=blobdiff_plain;f=epclust%2Fsrc%2FcomputeMedoidsIndices.cpp;h=895031a90f71376f10be93d4efc16a3a1faab94e;hb=a52836b23adb4bfa6722642ec6426fb7b5f39650;hp=693418118f9f3d38d1007c3923ea534c44703ec9;hpb=363ae13430cdee6ba76b42b7316aa4b292b04d93;p=epclust.git

diff --git a/epclust/src/computeMedoidsIndices.cpp b/epclust/src/computeMedoidsIndices.cpp
index 6934181..895031a 100644
--- 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.nrow(),
-		K = pMed->nrow(),
-		L = pMed->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(i,k) - *(medoids[k]+j);
+				// 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)