change perceptionL interface

.gitignore

@@ -4,4 +4,5 @@
 *.o
 *.dSYM
 build/
+*.dat
 .vscode/

include/cma-es/cmaes.h

@@ -1524,7 +1524,7 @@ template <typename Scalar, typename Ordered = Scalar> class CMAES {
    * that contains the matched stop criteria via getStopMessage().
    * @return Does any stop criterion match?
    */
-  bool testForTermination(bool makeMessage = false) {
+  unsigned long testForTermination(bool makeMessage = false) {
     Ordered range;
     Scalar fac;
     int iAchse, iKoo;
@@ -1658,6 +1658,8 @@ template <typename Scalar, typename Ordered = Scalar> class CMAES {
                     << "] ";
           stopCriteria[6].second = message.str();
         }
+      } else if (makeMessage) {
+        stopCriteria[6].second = "NoEffectAxis: (empty)";
       }
     }
 
@@ -1677,9 +1679,7 @@ template <typename Scalar, typename Ordered = Scalar> class CMAES {
                 << iKoo << " without effect";
         stopCriteria[7].second = message.str();
       } else if (makeMessage) {
-        std::stringstream message;
-        message << "NoEffectCoordinate: (empty)";
-        stopCriteria[7].second = message.str();
+        stopCriteria[7].second = "NoEffectCoordinate: (empty)";
       }
     }
 
@@ -1705,11 +1705,12 @@ template <typename Scalar, typename Ordered = Scalar> class CMAES {
       }
     }
 
-    return std::accumulate(
-        stopCriteria.begin(), stopCriteria.end(), false,
-        [](bool accu, const std::pair<bool, std::string> crit) {
-          return accu || crit.first;
-        });
+    size_t crit_n = stopCriteria.size();
+    unsigned long flags = 0;
+    for (size_t i = 0; i < crit_n; ++i) {
+      flags |= (unsigned long)stopCriteria[i].first << i;
+    }
+    return flags;
   }
 
   /**

include/fitness.h

@@ -3,15 +3,36 @@
 
 #include "color.h"
 #include "lexi.h"
+#include <ostream>
 #include <vector>
 
 inline double offRange(double x, double a, double b);
 
 double offRGB(const color::LAB &lab);
 
-LexiProduct<double> fitnessFunc(double L, const std::vector<double> &x);
+LexiProduct<double> fitnessFunc(const std::vector<color::LAB> &lab);
 
-std::pair<std::vector<color::RGB>, LexiProduct<double>>
-perceptionL(double L, size_t M, bool quiet = false);
+class PerceptionResult {
+public:
+  friend std::ostream &operator<<(std::ostream &os,
+                                  const PerceptionResult &res) {
+    os << "{" << std::endl << "  flags: " << res.flags << std::endl;
+    os << "  L: " << res.L << std::endl;
+    os << "  rgb: [" << std::endl;
+    for (const auto &rgb : res.rgb)
+      os << "    0: " << rgb << std::endl;
+    os << "  ]" << std::endl;
+    os << "  fitness: " << res.fitness << std::endl << "}";
+    return os;
+  }
+
+  unsigned long flags;
+  double L;
+  std::vector<color::RGB> rgb;
+  LexiProduct<double> fitness;
+};
+
+PerceptionResult perceptionL(color::LAB foreground, color::LAB background,
+                             size_t M, bool quiet = false);
 
 #endif

include/lexi.h

@@ -64,9 +64,19 @@ template <typename Scalar> class LexiProduct {
 
   Scalar &operator[](size_t i) { return this->prd[i]; }
 
-private:
+  friend std::ostream &operator<<(std::ostream &os, const LexiProduct &lexi) {
+    const size_t n = lexi.prd.size();
+    os << "[";
+    for (size_t i = 0; i < n - 1; ++i) {
+      os << lexi.prd[i] << ",";
+    }
+    os << lexi.prd[n - 1] << "]";
+    return os;
+  }
+
   Product prd;
 
+private:
   inline Scalar lexicmp(const LexiProduct<Scalar> &lhs,
                         const LexiProduct<Scalar> &rhs) const {
     const Scalar lencmp = (Scalar)(lhs.prd.size() - rhs.prd.size());
@@ -79,16 +89,6 @@ template <typename Scalar> class LexiProduct {
     }
     return lencmp;
   }
-
-  friend std::ostream &operator<<(std::ostream &os, const LexiProduct &lexi) {
-    const size_t n = lexi.prd.size();
-    os << "[";
-    for (size_t i = 0; i < n - 1; ++i) {
-      os << lexi.prd[i] << ",";
-    }
-    os << lexi.prd[n - 1] << "]";
-    return os;
-  }
 };
 
 #endif

src/fitness.cpp

@@ -26,23 +26,22 @@ double offRGB(const color::LAB &lab) {
   return std::sqrt(dr * dr + dg * dg + db * db);
 }
 
-LexiProduct<double> fitnessFunc(double L, const std::vector<double> &x) {
+LexiProduct<double> fitnessFunc(const std::vector<color::LAB> &lab) {
   // number of colors
-  size_t M = x.size() / 2;
-  std::vector<color::LAB> lab(M);
-  // off-RGB penalty
+  size_t M = lab.size();
   std::vector<double> penalty(M);
   for (size_t i = 0; i < M; ++i) {
-    lab[i] = color::LAB{L, x[2 * i], x[2 * i + 1]};
     penalty[i] =
         offRGB(lab[i]) * 300 * M; /* 300 is the approx. diameter of ab-plane */
   }
 
-  size_t K = M * (M - 1) / 2; /* number of combinations between colors */
+  // number of combinations between colors, except the last two
+  size_t K = M * (M - 1) / 2 - 1; 
   static size_t cached_K = 0;
   static std::vector<size_t> cached_ref; /* for almost-sorted O(n) sorting */
   static std::vector<std::pair<size_t, size_t>> ref_combi;
 
+  // Cache combination order
   if (cached_K != K) {
     cached_K = K;
     cached_ref.resize(K);
@@ -51,7 +50,7 @@ LexiProduct<double> fitnessFunc(double L, const std::vector<double> &x) {
     }
     ref_combi.resize(K);
     size_t combi_i = 0;
-    for (size_t i = 0; i < M - 1; ++i) {
+    for (size_t i = 0; i < M - 2; ++i) {
       for (size_t j = i + 1; j < M; ++j) {
         ref_combi[combi_i] = std::pair<size_t, size_t>(i, j);
         ++combi_i;
@@ -81,51 +80,74 @@ LexiProduct<double> fitnessFunc(double L, const std::vector<double> &x) {
   return LexiProduct<double>(std::move(fitness));
 }
 
+void fill_lab(const std::vector<double> &x, std::vector<color::LAB> &lab,
+              size_t M) {
+  for (size_t j = 0; j < M; ++j) {
+    lab[j] = color::LAB{x[2 * M], x[2 * j], x[2 * j + 1]};
+  }
+}
+
 /*
- * @param L luminosity
+ * @param foreground
+ * @param background
  * @param M numbers of colors
  * @param quiet write info to stdout
- * @return pair(RGBs, bestFitness)
+ * @return PerceptionResult
  */
-std::pair<std::vector<color::RGB>, LexiProduct<double>>
-perceptionL(double L, size_t M, bool quiet) {
+PerceptionResult perceptionL(color::LAB foreground, color::LAB background,
+                             size_t M, bool quiet) {
   CMAES<double, LexiProduct<double>> evo;
   Individual<double, LexiProduct<double>> *pop;
 
-  const size_t N = M * 2;           //!< number of variables
-  const size_t K = M * (M - 1) / 2; //!< number of combinations
+  const size_t N = M * 2 + 1; //!< number of variables
+  std::vector<double> xstart(N);
+  const double initL = (foreground.l + background.l) / 2.;
+  std::vector<double> stddev(N, std::min(100. - initL, initL));
+  xstart[2 * M] = initL; // luminocity
+  stddev[2 * M] = 100.;  // luminocity
+
   Parameters<double, LexiProduct<double>> parameters;
-  parameters.lambda = (int)(100. * log(N));
-  // parameters.diagonalCov = 1;
-  // parameters.stopTolFun = 1e-12;
-  // parameters.stopTolFunHist = 1e-13;
-  const std::vector<double> xstart(N, 0.), stddev(N, std::min(100. - L, L));
+  parameters.lambda = (int)(300. * log(N)); /* 100x default */
+  // parameters.stopTolFun = {0, 0, 1e-13};
+  // parameters.stopTolFunHist = {0, 0, 1e-13};
   parameters.init((int)N, xstart.data(), stddev.data());
   evo.init(parameters);
 
   if (!quiet)
     std::cout << evo.sayHello() << std::endl;
 
+  std::vector<color::LAB> lab(M + 2);
+  lab[M] = foreground;
+  lab[M + 1] = background;
+
   std::vector<LexiProduct<double>> populationFitness((size_t)parameters.lambda);
-  while (!evo.testForTermination()) {
+  unsigned long flags = 0;
+  while (!(flags = evo.testForTermination())) {
     pop = evo.samplePopulation(); // Do not change content of pop
-    for (size_t i = 0; i < evo.get(CMAES<double, LexiProduct<double>>::Lambda);
-         ++i)
-      populationFitness[i] = fitnessFunc(50., pop[i].x);
+    for (size_t i = 0; i < parameters.lambda; ++i) {
+      fill_lab(pop[i].x, lab, M);
+      populationFitness[i] = fitnessFunc(lab);
+    }
     evo.updateDistribution(populationFitness);
   }
 
-  if (!quiet)
-    std::cout << "Stop:" << std::endl << evo.getStopMessage();
-
-  auto x = evo.getVector(CMAES<double, LexiProduct<double>>::XMean);
+  const auto xfinal = evo.getVector(CMAES<double, LexiProduct<double>>::XMean);
+  fill_lab(xfinal, lab, M);
+  std::sort(lab.begin(), lab.begin() + (long)M, [](const color::LAB &c1, const color::LAB &c2) {
+    const double h1 = std::atan2(-c1.b, -c1.a);
+    const double h2 = std::atan2(-c2.b, -c2.a);
+    return h1 < h2;
+  });
   auto fitness = evo.getValue(CMAES<double, LexiProduct<double>>::FBestEver);
 
+  if (!quiet)
+    std::cout << "Stop: flags " << flags << std::endl << evo.getStopMessage();
+
   std::vector<color::RGB> rgb(M);
   for (size_t i = 0; i < M; ++i) {
-    rgb[i] = color::LABtoRGB(color::LAB{50., x[2 * i], x[2 * i + 1]});
+    rgb[i] = color::LABtoRGB(lab[i]);
   }
 
-  return std::pair<std::vector<color::RGB>, LexiProduct<double>>(
-      std::move(rgb), std::move(fitness));
+  return PerceptionResult{flags, xfinal[2 * M], std::move(rgb),
+                          std::move(fitness)};
 }

test/testfitness.cpp

@@ -43,11 +43,16 @@ int testcolor() {
 
   check(offRGB(color::LAB{200, 0, 0}), 2.147648383, "offRGB");
 
-  check(fitnessFunc(50, std::vector<double>{0, 0, 10, 10, -10, 10, -10, -10, 10, -10})[0],
-  -12.8001, "fitnessFunc");
+  check(fitnessFunc(std::vector<color::LAB>{{50, 0, 0},
+                                            {50, 10, 10},
+                                            {50, -10, 10},
+                                            {50, -10, -10},
+                                            {50, 10, -10}})[0],
+        -12.8001, "fitnessFunc");
 
-  check(fitnessFunc(50, std::vector<double>{-200, -200, +200, +200})[0],
-  1626.445599, "fitnessFunc (off boundary)");
+  check(
+      fitnessFunc(std::vector<color::LAB>{{50, -200, -200}, {50, 200, 200}})[0],
+      1626.445599, "fitnessFunc (off boundary)");
 
   std::cout << std::endl;
   int ret = EXIT_SUCCESS;

test/testfitness2.cpp

@@ -6,14 +6,37 @@
 #include <cma-es/cmaes.h>
 #include <color.h>
 #include <fitness.h>
+#include <fstream>
 #include <lexi.h>
 
+unsigned compress(double x, unsigned max = 256) {
+  return std::min((unsigned)std::max(0., x * max), max - 1);
+}
+
 int main(int argc, char *argv[]) {
-  const auto pair = perceptionL(50., 7);
+  const size_t M = 7;
+  const color::LAB foreground{0, 0., 0.};
+  const color::LAB background{100, 0., 0.};
+  const auto result = perceptionL(foreground, background, M);
+  std::cout << result << std::endl;
 
-  for (size_t i = 0; i < pair.first.size(); ++i) {
-    std::cout << pair.first[i] << " ";
+  // Write a HTML demo
+  const color::RGB frgb(color::LABtoRGB(foreground));
+  const color::RGB brgb(color::LABtoRGB(background));
+  std::ofstream myfile;
+  myfile.open("index.html");
+  myfile << "<!DOCTYPE html><html><head><style>body {background-color: rgb("
+         << compress(brgb.r) << "," << compress(brgb.g) << ","
+         << compress(brgb.b)
+         << ");} div {margin: 10px; float: left; "
+            "width: 30px; height: 30px;}</style></head><body>";
+  myfile << "<div style=\"background-color: rgb(" << compress(frgb.r) << ","
+         << compress(frgb.g) << "," << compress(frgb.b) << ");\"></div>";
+  for (const auto &rgb : result.rgb) {
+    myfile << "<div style=\"background-color: rgb(" << compress(rgb.r) << ","
+           << compress(rgb.g) << "," << compress(rgb.b) << ");\"></div>";
   }
-  std::cout << std::endl << pair.second << std::endl;
+  myfile << "</body></html>";
+  myfile.close();
   return 0;
 }