added jensen shannon, fixe kullback leibler when P=0

math/divergence/JensenShannon.h

@@ -0,0 +1,28 @@
+#ifndef JENSENSHANNON_H
+#define JENSENSHANNON_H
+
+#include "KullbackLeibler.h"
+
+#include "../../Assertions.h"
+#include <string>
+
+
+namespace Divergence {
+
+template <typename Scalar> class JensenShannon {
+
+    /** Calculate the Jensen Shannon Divergece from a set of sample densities
+     *  Info: https://en.wikipedia.org/wiki/Jensen–Shannon_divergence
+     * @param P is the vector containing the densities of a set of samples
+     * @param Q is a vector containg the densities of the same samples set then P
+    */
+    static inline Scalar getGeneralFromSamples(Eigen::VectorXd P, Eigen::VectorXd Q, LOGMODE mode){
+        Eigen::VectorXd M = 0.5 * (P + Q);
+
+        return (0.5 * KullbackLeibler::getGeneralFromSamples(P, M, mode)) + (0.5 * KullbackLeibler::getGeneralFromSamples(Q, M, mode));
+    }
+};
+
+}
+
+#endif // JENSENSHANNON_H

math/divergence/KullbackLeibler.h

@@ -45,13 +45,18 @@ namespace Divergence {
                 //sum up the logarithmic difference between P and Q, also called kullback leibler...
                 for(int i = 0; i < P.size(); ++i){
 
+                    // if both prob are near zero we assume a 0 distance since lim->0 = 0
                     if((P[i] == 0.0) && (Q[i] == 0.0)){
-                            dist += 0.0;
+                        dist += 0.0;
-                    } else {
+                    }
+                    //if prob of P is 0 we also assume a 0 distance since lim->0 0 * log(0) = 0
+                    else if ((P[i] == 0.0) && (Q[i] > 0.0)){
+                        dist += 0.0;
+                    } else{
 
                         // calc PQratio
                         if(Q[i] == 0.0){
-                            Assert::doThrow("Division by zero is not allowed ;). TODO: What if the densities are to small?");
+                            Assert::doThrow("Division by zero is not allowed ;).");
                             //PQratio = P[i] / 0.00001;
                         } else {
                             PQratio = P[i] / Q[i];