1 /* PSPP - computes sample statistics.
2 Copyright (C) 1997-9, 2000 Free Software Foundation, Inc.
3 Written by Ben Pfaff <blp@gnu.org>.
5 This program is free software; you can redistribute it and/or
6 modify it under the terms of the GNU General Public License as
7 published by the Free Software Foundation; either version 2 of the
8 License, or (at your option) any later version.
10 This program is distributed in the hope that it will be useful, but
11 WITHOUT ANY WARRANTY; without even the implied warranty of
12 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
13 General Public License for more details.
15 You should have received a copy of the GNU General Public License
16 along with this program; if not, write to the Free Software
17 Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
22 #include <gsl/gsl_rng.h>
27 #include <data/procedure.h>
28 #include <data/variable.h>
29 #include <language/command.h>
30 #include <language/lexer/lexer.h>
31 #include <libpspp/alloc.h>
32 #include <libpspp/compiler.h>
33 #include <libpspp/message.h>
34 #include <libpspp/str.h>
35 #include <math/random.h>
38 #define _(msgid) gettext (msgid)
40 /* The two different types of samples. */
43 TYPE_A_FROM_B, /* 5 FROM 10 */
44 TYPE_FRACTION /* 0.5 */
47 /* SAMPLE transformation. */
50 int type; /* One of TYPE_*. */
51 int n, N; /* TYPE_A_FROM_B: n from N. */
52 int m, t; /* TYPE_A_FROM_B: # picked so far; # so far. */
53 unsigned frac; /* TYPE_FRACTION: a fraction of UINT_MAX. */
56 static trns_proc_func sample_trns_proc;
57 static trns_free_func sample_trns_free;
60 cmd_sample (struct dataset *ds)
62 struct sample_trns *trns;
68 if (!lex_force_num ())
70 if (!lex_is_integer ())
72 unsigned long min = gsl_rng_min (get_rng ());
73 unsigned long max = gsl_rng_max (get_rng ());
76 if (tokval <= 0 || tokval >= 1)
78 msg (SE, _("The sampling factor must be between 0 and 1 "
83 frac = tokval * (max - min) + min;
91 if (!lex_force_match_id ("FROM"))
93 if (!lex_force_int ())
98 msg (SE, _("Cannot sample %d observations from a population of "
108 trns = xmalloc (sizeof *trns);
112 trns->m = trns->t = 0;
114 add_transformation (ds, sample_trns_proc, sample_trns_free, trns);
116 return lex_end_of_command ();
119 /* Executes a SAMPLE transformation. */
121 sample_trns_proc (void *t_, struct ccase *c UNUSED,
122 casenumber case_num UNUSED)
124 struct sample_trns *t = t_;
127 if (t->type == TYPE_FRACTION)
129 if (gsl_rng_get (get_rng ()) <= t->frac)
130 return TRNS_CONTINUE;
132 return TRNS_DROP_CASE;
136 return TRNS_DROP_CASE;
138 U = gsl_rng_uniform (get_rng ());
139 if ((t->N - t->t) * U >= t->n - t->m)
142 return TRNS_DROP_CASE;
148 return TRNS_CONTINUE;
153 sample_trns_free (void *t_)
155 struct sample_trns *t = t_;