Index: tol/bmath/bprdist/prd.cpp
===================================================================
--- tol/bmath/bprdist/prd.cpp	(revision 902a6ae84beadb522e9d754c98bf15147db67896)
+++ tol/bmath/bprdist/prd.cpp	(revision fb55369f204a636e7dd8a7d4115f35de74f2849f)
@@ -2415,5 +2415,5 @@
 
 //--------------------------------------------------------------------
-BBinomialDist::BBinomialDist(BInt N, BDat p)
+BBinomialDist::BBinomialDist(BDat N, BDat p)
 
 /*! Creates an object aimed to work whith a binomial distribution. 
@@ -2430,5 +2430,5 @@
     logp_ = Log(p_);
     logq_ = Log(q_);
-    logGammaN_ = LogFact(N_);
+    logGammaN_ = LogGamma(1+N_);
     if((N_<=0)||(p_<0)||(p_>1))
     {
@@ -2451,5 +2451,5 @@
     BInt k = (BInt)Floor(x).Value();
     if((k<0)||(k>N_)) { return(0); }
-    return(Exp(logGammaN_-LogFact(N_-k)-LogFact(k)+k*logp_+(N_-k)*logq_));
+    return(Exp(logGammaN_-LogGamma(1+N_-k)-LogFact(k)+k*logp_+(N_-k)*logq_));
 }
 
@@ -2471,5 +2471,5 @@
 	
 	s = k;
-	xn = N_;
+	xn = N_.Value();
 	pr = p_.Value();
 	ompr = 1 - pr;
@@ -2494,5 +2494,5 @@
     p = x.Value();
     q = 1-p;
-    xn = (double) N_;
+    xn = N_.Value();
     pr = p_.Value();
     ompr = 1-pr;
@@ -2542,9 +2542,9 @@
 
 //--------------------------------------------------------------------
-BNegBinomialDist::BNegBinomialDist(BInt N, BDat P)
+BNegBinomialDist::BNegBinomialDist(BDat N, BDat P)
 
 /*! Creates an object aimed to work whith a negative binomial distribution. 
- *  A negative binomial distribution counts the number of failures
- *  before a fixed number of successes whith a fixed success probability.
+ *  A negative binomial distribution counts the number of successes
+ *  before a fixed number of failures whith a fixed success probability.
  * \param N Number of succeses.
  * \param P Success probability.
@@ -2576,5 +2576,5 @@
 
 /*! Returns NegBinomial(x,N_,P_) 
- *  = prob of getting x failures before the N_-th success
+ *  = prob of getting x succeses before the N_-th failures
  *
  *    ALGORITHM: NB(x,N_,P_) = P_ * DensBinomial(N_-1,k+N_-1,P_)
@@ -2583,12 +2583,7 @@
 //--------------------------------------------------------------------
 {
-    if(wrongParameter_ || x.IsUnknown()) { return(BDat::Unknown()); }
-    BInt k = (BInt)Floor(x).Value();
-    if((k<0)) { return(0); }
-    if((k==0)&&(N_==1)){return P_;}
-    else{
-	BBinomialDist b(k+N_-1,P_);
-	return(b.Dens(N_-1)*P_);    
-    } 
+  if(wrongParameter_ || x.IsUnknown()) { return(BDat::Unknown()); }
+  BInt k = (BInt)Floor(x).Value();
+  return(gsl_ran_negative_binomial_pdf(k,P_.Value(),N_.Value()));
 }
 
@@ -2601,8 +2596,13 @@
 //--------------------------------------------------------------------
 {
+  if(wrongParameter_ || x.IsUnknown()) { return(BDat::Unknown()); }
+  BInt k = (BInt)Floor(x).Value();
+//return(1-gsl_sf_beta_inc (k+1, N_.Value(), 1-P_.Value()));
+//return(gsl_cdf_negative_binomial_P(k,P_.Value(),N_.Value()));
+/* */
     double p, q, s, xn, pr, ompr, bound;
     int which, status; 
     s = x.Value();
-    xn = (double) N_;
+    xn = N_.Value();
     pr = P_.Value();
     ompr = 1-pr;
@@ -2613,4 +2613,5 @@
     else
 	return(BDat::Unknown());
+/* */
 }
 
@@ -2623,9 +2624,11 @@
 {
   if(wrongParameter_ || x.IsUnknown()) { return(BDat::Unknown()); }
+//return(BDiscreteDist::Inverse(x,tolerance));
+
     double p, q, s, xn, pr, ompr, bound;
     int which, status; 
     p = x.Value();
     q = 1-p;
-    xn = (double) N_;
+    xn = N_.Value();
     pr = P_.Value();
     ompr = 1-pr;
@@ -2641,9 +2644,9 @@
 	
 	cumnbn(&x0,&xn,&pr,&ompr,&fx0,&q);
-	/* may be, there is an incremental way */
 	cumnbn(&x1,&xn,&pr,&ompr,&fx1,&q);
 	int K1 = 1;
 	return fabs(fx1-p)<=spmpar(&K1) ? x1 : x0;
     }
+
 }
 
@@ -2778,4 +2781,5 @@
 //--------------------------------------------------------------------
 {
+  if(wrongParameter_ || p.IsUnknown()) { return(BDat::Unknown()); }
   return(BDiscreteDist::Inverse(p,tolerance));
 /*
@@ -2865,4 +2869,80 @@
 
 
+
+
+//--------------------------------------------------------------------
+BGenCountDist::BGenCountDist(BDat a, BDat v)
+
+/*! Creates an object aimed to deal with a Generic Counting distribution
+ *  depending just in average and variance parameters and that could be 
+ *    - Binomial:                 average>variance
+ *    - Poisson:                  average=variance
+ *    - Negative Binomial:        average<variance
+ */
+//--------------------------------------------------------------------
+:  BDiscreteDist(), a_(a), v_(v), 
+   cd_(NULL), m_(BDat::Unknown()), p_(BDat::Unknown())
+{
+  subType_ = Sign(a_-v_).Value();
+  m_ = a_*a_/Abs(a_-v_);
+  p_ = 1.0 - min(v_/a_,a_/v_);
+  if(subType_==1)
+  {
+    cd_ = new BBinomialDist(m_.Value(),p_);
+  }
+  else if(subType_==0)
+  {
+    cd_ = new BPoissonDist(a_);
+  }
+  if(subType_==-1)
+  {
+    cd_ = new BNegBinomialDist(m_.Value(),1-p_);
+  }
+}
+
+//--------------------------------------------------------------------
+ BGenCountDist::~BGenCountDist()
+//--------------------------------------------------------------------
+{
+  if(cd_) { delete cd_; }
+}
+
+//--------------------------------------------------------------------
+BDat BGenCountDist::Dens(BDat x)
+
+/*! Returns PDF
+ */
+//--------------------------------------------------------------------
+{
+  BDat d;
+  if(cd_) { d = cd_->Dens(x); }
+  return(d);
+}
+
+
+//--------------------------------------------------------------------
+BDat BGenCountDist::Dist(BDat x)
+
+/*! Returns CDF
+ */
+//--------------------------------------------------------------------
+{
+  BDat d;
+  if(cd_) { d = cd_->Dist(x); }
+  return(d);
+}
+
+//--------------------------------------------------------------------
+BDat BGenCountDist::Inverse(BDat p, BDat tolerance)
+
+/*! Returns the maximum x matching Poisson(x,fi) <= p
+ */
+//--------------------------------------------------------------------
+{
+  BDat i;
+  if(cd_) { i = cd_->Inverse(p,tolerance); }
+  return(i);
+}
+
 //--------------------------------------------------------------------
 BHypergeometricDist::BHypergeometricDist(BInt N, BInt P, BInt n)
Index: tol/bmath/bprdist/tol_bprdist.h
===================================================================
--- tol/bmath/bprdist/tol_bprdist.h	(revision 902a6ae84beadb522e9d754c98bf15147db67896)
+++ tol/bmath/bprdist/tol_bprdist.h	(revision fb55369f204a636e7dd8a7d4115f35de74f2849f)
@@ -873,5 +873,5 @@
 {
 private:
-  BInt N_;
+  BDat N_;
   BDat p_;
   BDat q_;
@@ -882,5 +882,5 @@
 public:
   // contructors and destructors:
-  BBinomialDist(BInt N, BDat p);
+  BBinomialDist(BDat N, BDat p);
 
     // Implementation: prd.cpp
@@ -901,5 +901,5 @@
 {
 private:
-  BInt N_;
+  BDat N_;
   BDat p_;
   BDat q_;
@@ -912,5 +912,5 @@
 public:
   // contructors and destructors:
-  BNegBinomialDist(BInt N, BDat P);
+  BNegBinomialDist(BDat N, BDat P);
 
   // Implementation: prd.cpp
@@ -959,4 +959,31 @@
 };
 
+//--------------------------------------------------------------------
+class TOL_API BGenCountDist : public BDiscreteDist
+
+/*! Generic counting probability distribution definition.
+ */
+//--------------------------------------------------------------------
+{
+private:
+  BDat a_; //average
+  BDat v_; //variance
+  int subType_;
+  BDiscreteDist* cd_;
+  BDat m_;
+  BDat p_;
+
+public:
+  // contructors and destructors:
+  BGenCountDist(BDat a, BDat v);
+ ~BGenCountDist();
+
+  // Implementation: prd.cpp
+  BDat Average	() { return(a_); }
+  BDat Varianze () { return(v_); }
+  BDat Dist	(BDat x);
+  BDat Dens	(BDat x);
+  BDat Inverse  (BDat prob, BDat tolerance = 0.1);
+};
 
 //--------------------------------------------------------------------