freqtest.c

随机数算法

/*
          Frequency test within a block.

  This is slightly different than the NIST test, but it does the
  same thing using their suggestion of -1 for 0 bits and +1 for
  1 bits in the block.  This moves the average from 1/2 as they 
  have, to zero.  This simplifies the computation of the chi^2
  observed value to

              chi^2(observed) = (sum(P_i)^2)/M

  where M = number of bits in a block
        N = number of blocks to do
	P_i = sum( bits )  (+1 if bit set, -1 if bit clear)

  Enter with source pointer, endianness, size of each block and
    number of blocks to work on.
  Returns chi^2 value and source pointer updated.
*/

#include "randtest.h"

double freqtest( SRCPTR *src, int endian, int blocksize, int numblocks)
{
  int sumeps, i;
  unsigned char mask;
  double sum = 0.0;

  mask = 1 << src->bitpos;
  while( numblocks)
  {
    sumeps = 0;
    for( i=0; i<blocksize; i++)
    {
      if( *(src->byteptr) & mask) sumeps += 1;
      else sumeps -= 1;
      if( endian)
      {
	mask >>= 1;
	src->bitpos--;
	if (src->bitpos < 0)
	{
	  mask = 0x80;
	  src->bitpos = 7;
	  src->byteptr++;
	}
      }
      else
      {
	mask <<= 1;
	src->bitpos++;
	if (src->bitpos > 7)
	{
	  mask = 1;
	  src->bitpos = 0;
	  src->byteptr++;
	}
      }
    }
    sum += sumeps*sumeps;
    numblocks--;
  }
  sum /= blocksize;
  return (sum);
}

/*  check math - ensure we get the same answer that theory says
    we ought to.  It does work, exactly the same as expected :-)

double freqtestck( SRCPTR *src, int endian, int blocksize, int numblocks)
{
  int i;
  unsigned char mask;
  double sum = 0.0;
  double sumeps;

  mask = 1 << src->bitpos;
  while( numblocks)
  {
    sumeps = 0.0;
    for( i=0; i<blocksize; i++)
    {
      if( *(src->byteptr) & mask) sumeps += 1.0;
      if( endian)
      {
	mask >>= 1;
	src->bitpos--;
	if (src->bitpos < 0)
	{
	  mask = 0x80;
	  src->bitpos = 7;
	  src->byteptr++;
	}
      }
      else
      {
	mask <<= 1;
	src->bitpos++;
	if (src->bitpos > 7)
	{
	  mask = 1;
	  src->bitpos = 0;
	  src->byteptr++;
	}
      }
    }
    numblocks--;
    sumeps = sumeps/blocksize - 0.5;
    sum += sumeps*sumeps;
  }
  sum *= 4.0*blocksize;
  return sum;
}
*/