3 * $Id: maurer.c,v 1.2 2000/08/11 21:34:59 mdw Exp $
5 * Maurer's universal statistical test
7 * (c) 2000 Straylight/Edgeware
10 /*----- Licensing notice --------------------------------------------------*
12 * This file is part of Catacomb.
14 * Catacomb is free software; you can redistribute it and/or modify
15 * it under the terms of the GNU Library General Public License as
16 * published by the Free Software Foundation; either version 2 of the
17 * License, or (at your option) any later version.
19 * Catacomb is distributed in the hope that it will be useful,
20 * but WITHOUT ANY WARRANTY; without even the implied warranty of
21 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
22 * GNU Library General Public License for more details.
24 * You should have received a copy of the GNU Library General Public
25 * License along with Catacomb; if not, write to the Free
26 * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
30 /*----- Revision history --------------------------------------------------*
33 * Revision 1.2 2000/08/11 21:34:59 mdw
34 * New restartable interface to Maurer testing.
36 * Revision 1.1 2000/06/17 11:29:49 mdw
37 * Maurer's universal statistical test.
41 /*----- Header files ------------------------------------------------------*/
48 #include <mLib/alloc.h>
49 #include <mLib/bits.h>
53 /*----- Data structures ---------------------------------------------------*/
55 typedef struct bitscan
{
61 /*----- Main code ---------------------------------------------------------*/
65 * Arguments: @bitscan *b@ = pointer to a bitscanner structure
66 * @unsigned l@ = number of bits wanted
68 * Returns: An integer representing the next @l@ bits of the stream.
70 * Use: Fetches the next @l@ bits of a data stream. The bitscanner
71 * has no way of knowing whether it's run out of data, so make
72 * sure that doesn't happen.
75 static unsigned more(bitscan
*b
, unsigned l
)
80 b
->a
|= U8(*b
->p
++) << b
->b
;
90 /* --- @maurer_init@ --- *
92 * Arguments: @maurer_ctx *m@ = pointer to a context to initialize
93 * @unsigned l@ = number of bits to take at a time
95 * Returns: Minimum recommended amount of data to test.
97 * Use: Initializes a Maurer testing context. Note that @l@ values
98 * larger than 16 are not supported, and values less than 6 can
99 * give bizarre results.
102 unsigned long maurer_init(maurer_ctx
*m
, unsigned l
)
106 assert(((void)"(maurer_init): chunks larger than 16 bits not supported",
115 m
->t
= xmalloc(sizeof(unsigned long) << l
);
116 for (i
= 0; i
< (1 << l
); i
++)
118 return (((1000ul << l
) * l
+ 7)/8);
123 * Arguments: @maurer_ctx *m@ = pointer to testing context
124 * @const octet **p@ = address of a buffer pointer
125 * @const octet *q@ = limit of the buffer pointer
126 * @unsigned *x@ = address to store next chunk
128 * Returns: Nonzero if some more bits arrived.
130 * Use: Reads some bits from a buffer.
133 static int bits(maurer_ctx
*m
, const octet
**p
, const octet
*q
, unsigned *x
)
135 while (m
->b
< m
->l
) {
138 m
->a
|= U8(*(*p
)++) << m
->b
;
141 *x
= m
->a
& ((1 << m
->l
) - 1);
148 /* --- @maurer_test@ --- *
150 * Arguments: @maurer_ctx *m@ = pointer to a Maurer testing context
151 * @const void *buf@ = pointer to a buffer of data
152 * @size_t sz@ = size of the buffer
156 * Use: Scans a buffer of data and updates the testing context.
159 void maurer_test(maurer_ctx
*m
, const void *buf
, size_t sz
)
161 const octet
*p
= buf
, *l
= p
+ sz
;
162 unsigned long q
= 10 << m
->l
;
165 /* --- Initialize the table --- */
168 if (!bits(m
, &p
, l
, &x
))
173 /* --- Update the statistic --- */
175 while (bits(m
, &p
, l
, &x
)) {
176 m
->x
+= log(m
->n
- m
->t
[x
]);
181 /* --- @maurer_done@ --- *
183 * Arguments: @maurer_ctx *m@ = pointer to a Maurer testing context
185 * Returns: The statistic %$Z_u = (X_u - \mu)/\sigma$%, which should be
186 * normally distributed with a mean of 0 and variance of 1.
188 * Use: Returns the result of Maurer's universal statistical test.
189 * Also frees the memory allocated during initialization.
191 * If insufficient data was received, the value @HUGE_VAL@ is
195 double maurer_done(maurer_ctx
*m
)
197 double x
, mu
, c
, sigma
;
198 unsigned long q
= 10 << m
->l
, k
;
200 /* --- Table for means and variances --- */
202 struct { double mu
, var_1
; } vtab
[] = {
203 { 0.7326495, 0.690 },
204 { 1.5374383, 1.338 },
205 { 2.4016068, 1.901 },
206 { 3.3112247, 2.358 },
207 { 4.2534266, 2.705 },
208 { 5.2177052, 2.945 },
209 { 6.1962507, 3.125 },
210 { 7.1836656, 3.238 },
211 { 8.1764248, 3.311 },
212 { 9.1723243, 3.356 },
213 { 10.170032 , 3.384 },
214 { 11.168765 , 3.401 },
215 { 12.168070 , 3.410 },
216 { 13.167693 , 3.416 },
217 { 14.167488 , 3.419 },
218 { 15.167379 , 3.421 }
221 /* --- Check for bogus requests --- */
228 /* --- Do the maths --- *
230 * This produces %$X_u$%, which should be normally distributed with a mean
231 * and variance we can compute.
235 x
= m
->x
/(k
* log(2));
237 /* --- Now translate this into %$Z_u$% distributed as %$N(0, 1)$% --- */
239 mu
= vtab
[m
->l
- 1].mu
;
240 c
= 0.7 - 0.8/m
->l
+ (1.6 + 12.8/m
->l
) * pow(k
, -4.0/m
->l
);
241 sigma
= sqrt(c
* c
* vtab
[m
->l
- 1].var_1
/ k
);
244 /* --- Done (phew!) --- */
251 /* --- @maurer@ --- *
253 * Arguments: @const octet *p@ = pointer to a buffer of data
254 * @size_t sz@ = size of the buffer
255 * @unsigned l@ = number of bits to take at a time
257 * Returns: The statistic %$Z_u$%.
259 * Use: Simple interface to Maurer's universal statistical test. For
260 * large %$l$% you should use the more complicated restartable
264 double maurer(const octet
*p
, size_t sz
, unsigned l
)
269 maurer_test(&m
, p
, sz
);
270 return (maurer_done(&m
));
273 /*----- That's all, folks -------------------------------------------------*/