Small cleanups

* j'ai essayé d'utiliser les notations de l'article
* il faut faire attention que c'est N-I et non pas juste N
* tu as l'air de beaucoup aimer push_back. C'est très pratique
  quand on ne connait pas la taille du tableau à l'avance, mais
  si on la connait, je trouve que la notation a[i]=x est au moins
  aussi lisible.

1 files changed, 16 insertions, 23 deletions

diff --git a/stratified_sampling.cpp b/stratified_sampling.cpp
index 4703aba..af716bc 100644
--- a/stratified_sampling.cpp
+++ b/stratified_sampling.cpp
@@ -4,13 +4,9 @@
 #include "rtnorm.hpp"
 #include <gsl/gsl_cdf.h>
 #include <gsl/gsl_math.h>
-#define LOW_DISCREPANCY_HPP
+
 #include <cmath>
-#include <climits>
-#include <list>
 #include <algorithm>
-#include <numeric>
-
 
 //--génération quantiles--
 std::vector<double> quantile_norm(int n, double sigma){
@@ -45,32 +41,29 @@ std::pair<double, double> mean_var( std::vector<double> r){
 //actualisation du nombre de tirages à faire par strates
 std::vector<int> update_sampling (std::vector<double> p, 
 std::vector<double> sigma, int n) {
-	std::vector<int> r; // notre vecteur final à retourner
-	std::vector<double> m; //le vecteur des Mi idéals 
-	int s = p.size();
-	if (sigma.empty()) {
-		for (int i = 0; i<s ; i++) {
-			m.push_back(double(n)*p[i]);
+    int I = p.size();
+	std::vector<int> M(I, 1); // notre vecteur final à retourner
+	std::vector<double> m(I, 0); //le vecteur des m_i idéals
+
+    if (sigma.empty()) {
+		for (int i=0; i<I; i++) {
+			m[i] = (Nk-I)*p[i];
 		}
 	}
 	else {
 		double scal = std::inner_product(p.begin(), p.end(), sigma.begin(), (double) 0);
-		for (int i=0; i<s; i++) {
-			double v = p[i]*sigma[i]/scal;
-			m.push_back(v*n);
+		for (int i=0; i < I; i++) {
+			m[i] = (Nk-I)*p[i]*sigma[i]/scal;
 			std::cout<<m[i]<<std::endl;
 		}
 	}
-	r.push_back((int)(m[0]));
-	double v1 = 0;
-	for (int i=1; i<s-1; i++){
-		for (int j=0; j<i; j++) {
-			v1=v1+m[j];
-		}
-		r.push_back((int)v1 - (int)(v1-m[i]));
+    M[0]+=floor(m[0]);
+    double current = m[0];
+	for (int i=1; i<I; i++){
+        M[i] += floor(current+m[i]) - floor(current);
+        current += m[i];
 	}
-	r.push_back(n - std::accumulate( r.begin(), r.end(), 0 ));
-	return r;
+	return M;
 }