[Libreoffice-commits] core.git: sc/source

[Winfried Donkers]

 sc/source/core/tool/interpr3.cxx |   73 +++++++++++++++++++++++++++++++--------
 1 file changed, 59 insertions(+), 14 deletions(-)

New commits:
commit fa7085bd182e1968237b742c544de3ff0a5b9742
Author: Winfried Donkers <winfrieddonkers at libreoffice.org>
Date:   Mon Nov 4 10:32:40 2013 +0100

    improve calc function CRITBINOM
    
    CRITBINOM could not handle large sample quantities, e.g.
    CRITBINOM(1080,0.5,0.1).
    This patch fixes that and also improves the efficiency of
    the function by always iterating from the end nearest to
    the cumulative binomial distribution.
    
    Change-Id: I35f965acc83e26141fb2cf5a7e3303ec0791f04a
    Reviewed-on: https://gerrit.libreoffice.org/6559
    Reviewed-by: Eike Rathke <erack at redhat.com>
    Tested-by: Eike Rathke <erack at redhat.com>

diff --git a/sc/source/core/tool/interpr3.cxx b/sc/source/core/tool/interpr3.cxx
index e7b071f..3f665f3 100644
--- a/sc/source/core/tool/interpr3.cxx
+++ b/sc/source/core/tool/interpr3.cxx
@@ -1383,15 +1383,50 @@ void ScInterpreter::ScCritBinom()
             PushIllegalArgument();
         else
         {
+            double fFactor;
             double q = (0.5 - p) + 0.5;           // get one bit more for p near 1.0
-            double fFactor = pow(q,n);
-            if (fFactor <= ::std::numeric_limits<double>::min())
+            if ( q > p )                          // work from the side where the cumulative curve is
             {
-                fFactor = pow(p, n);
-                if (fFactor <= ::std::numeric_limits<double>::min())
-                    PushNoValue();
+                // work from 0 upwards
+                fFactor = pow(q,n);
+                if (fFactor > ::std::numeric_limits<double>::min())
+                {
+                    double fSum = fFactor;
+                    sal_uInt32 max = static_cast<sal_uInt32> (n), i;
+                    for (i = 0; i < max && fSum < alpha; i++)
+                    {
+                        fFactor *= (n-i)/(i+1)*p/q;
+                        fSum += fFactor;
+                    }
+                    PushDouble(i);
+                }
                 else
                 {
+                    // accumulate BinomDist until accumulated BinomDist reaches alpha
+                    double fSum = 0.0, x;
+                    sal_uInt32 max = static_cast<sal_uInt32> (n), i;
+                    for (i = 0; i < max && fSum < alpha; i++)
+                    {
+                        x = GetBetaDistPDF( p, ( i + 1 ), ( n - i + 1 ) )/( n + 1 );
+                        if ( !nGlobalError )
+                        {
+                            fSum += x;
+                        }
+                        else
+                        {
+                            PushNoValue();
+                            return;
+                        }
+                    }
+                    PushDouble( i - 1 );
+                }
+            }
+            else
+            {
+                // work from n backwards
+                fFactor = pow(p, n);
+                if (fFactor > ::std::numeric_limits<double>::min())
+                {
                     double fSum = 1.0 - fFactor;
                     sal_uInt32 max = static_cast<sal_uInt32> (n), i;
                     for (i = 0; i < max && fSum >= alpha; i++)
@@ -1401,17 +1436,27 @@ void ScInterpreter::ScCritBinom()
                     }
                     PushDouble(n-i);
                 }
-            }
-            else
-            {
-                double fSum = fFactor;
-                sal_uInt32 max = static_cast<sal_uInt32> (n), i;
-                for (i = 0; i < max && fSum < alpha; i++)
+                else
                 {
-                    fFactor *= (n-i)/(i+1)*p/q;
-                    fSum += fFactor;
+                    // accumulate BinomDist until accumulated BinomDist reaches alpha
+                    double fSum = 0.0, x;
+                    sal_uInt32 max = static_cast<sal_uInt32> (n), i;
+                    alpha = 1 - alpha;
+                    for (i = 0; i < max && fSum < alpha; i++)
+                    {
+                        x = GetBetaDistPDF( q, ( i + 1 ), ( n - i + 1 ) )/( n + 1 );
+                        if ( !nGlobalError )
+                        {
+                            fSum += x;
+                        }
+                        else
+                        {
+                            PushNoValue();
+                            return;
+                        }
+                    }
+                    PushDouble( n - i + 1 );
                 }
-                PushDouble(i);
             }
         }
     }