projects
/
pspp-builds.git
/ blobdiff
commit
grep
author
committer
pickaxe
?
search:
re
summary
|
shortlog
|
log
|
commit
|
commitdiff
|
tree
raw
|
inline
| side by side
Reduce number of multiplications for higher moments.
[pspp-builds.git]
/
src
/
math
/
moments.c
diff --git
a/src/math/moments.c
b/src/math/moments.c
index ce2f825c8e0156cc0e3d27dca54415da74512906..17040e9f699335fe91021cb5331ce1f257574ef5 100644
(file)
--- a/
src/math/moments.c
+++ b/
src/math/moments.c
@@
-1,6
+1,5
@@
/* PSPP - computes sample statistics.
Copyright (C) 1997-9, 2000 Free Software Foundation, Inc.
/* PSPP - computes sample statistics.
Copyright (C) 1997-9, 2000 Free Software Foundation, Inc.
- Written by Ben Pfaff <blp@gnu.org>.
This program is free software; you can redistribute it and/or
modify it under the terms of the GNU General Public License as
This program is free software; you can redistribute it and/or
modify it under the terms of the GNU General Public License as
@@
-22,9
+21,9
@@
#include <assert.h>
#include <math.h>
#include <stdlib.h>
#include <assert.h>
#include <math.h>
#include <stdlib.h>
-#include
"alloc.h"
-#include
"misc.h"
-#include
"value.h"
+#include
<libpspp/alloc.h>
+#include
<libpspp/misc.h>
+#include
<data/value.h>
#include "gettext.h"
#define _(msgid) gettext (msgid)
#include "gettext.h"
#define _(msgid) gettext (msgid)
@@
-44,11
+43,7
@@
calc_moments (enum moment max_moment,
if (max_moment >= MOMENT_VARIANCE && w > 1.)
{
if (max_moment >= MOMENT_VARIANCE && w > 1.)
{
- double s2;
-
- /* From _Numerical Recipes in C_, 2nd ed., 0-521-43108-5,
- section 14.1. */
- s2 = (d2 - pow2 (d1) / w) / (w - 1.);
+ double s2 = (d2 - pow2 (d1) / w) / (w - 1.);
if (variance != NULL)
*variance = s2;
if (variance != NULL)
*variance = s2;
@@
-156,8
+151,6
@@
moments_pass_one (struct moments *m, double value, double weight)
void
moments_pass_two (struct moments *m, double value, double weight)
{
void
moments_pass_two (struct moments *m, double value, double weight)
{
- double d, d_power;
-
assert (m != NULL);
if (m->pass == 1)
assert (m != NULL);
if (m->pass == 1)
@@
-169,28
+162,25
@@
moments_pass_two (struct moments *m, double value, double weight)
if (value != SYSMIS && weight >= 0.)
{
if (value != SYSMIS && weight >= 0.)
{
- m->w2 += weight;
-
- d = d_power = value - m->mean;
- m->d1 += d_power * weight;
-
+ double d = value - m->mean;
+ double d1_delta = d * weight;
+ m->d1 += d1_delta;
if (m->max_moment >= MOMENT_VARIANCE)
{
if (m->max_moment >= MOMENT_VARIANCE)
{
- d_power *= d;
- m->d2 += d_power * weight;
-
+ double d2_delta = d1_delta * d;
+ m->d2 += d2_delta;
if (m->max_moment >= MOMENT_SKEWNESS)
{
if (m->max_moment >= MOMENT_SKEWNESS)
{
- d_power *= d;
- m->d3 += d_power * weight;
-
+ double d3_delta = d2_delta * d;
+ m->d3 += d3_delta;
if (m->max_moment >= MOMENT_KURTOSIS)
{
if (m->max_moment >= MOMENT_KURTOSIS)
{
- d
_power *=
d;
- m->d4 += d
_power * weight
;
+ d
ouble d4_delta = d3_delta *
d;
+ m->d4 += d
4_delta
;
}
}
}
}
}
}
+ m->w2 += weight;
}
}
}
}