/* KonStat2.c. Copyright (C) 2002, 2003, 2004, 2005 by J. David Sexton. Some rights NOT reserved! This source code will be more readable if tabs are set equal to 4 spaces. KonStat2.c (and Konton2.c) should compile with any C compiler that meets the ISO 9899:1990 standard. KonStat2.c is not public domain or Open Source but I hereby license you to use it and give it to others so long as you don't profit by doing so. In addition, you may include the functions that actually implement the tests (the ones whose names end with "Stat"), HexTimesX, LDFract, ATimesLog2OfX, Log2OfX, TwoToTheX, Log2OfGamma, IncompleteGammaFactor, IncompleteGammaBySeries, IncompleteGammaByFraction, IncompleteGamma, ChiSquaredPValue, and the stopwatch functions in your own software; so long as that software is not allowed to be sold or used for profit. If I sold this, I would guarantee it; but as I don't, I won't. KonStat2.c will probably always be a work-in-progress. Keep looking for new revisions of it. This version was revised on Wednesday, February 9, 2005. More information is available at: http://www.geocities.com/da5id65536/ KonStat2.c, compiled and linked with with Konton2.c, performs statistical tests of randomness on random or pseudo-random bit streams. These bit streams can be files or they can be generated by stream ciphers: Konton2, RC4, or (optionally) Snow 2.0. To find out how to specify a file to be evaluated, invoke the compiled program with "-h" as the command-line argument. In the reports created by KonStat2.c: "P" is the tail probability [i.e., a P-value less than 0.1 indicates failure at the 0.1 level of significance] "MAD" is the mean absolute deviation of the P-values from their theoretical median of 0.5 "tidbit" is a 2-bit unsigned integer [4 possible values] "nibble" is a 4-bit unsigned integer [16 possible values] "byte" is an 8-bit unsigned integer [256 possible values] The P-value is a percentile; it is the probability of a statistic equal to or greater than the one reported. P-values for a given test should be more-or-less evenly distributed between 0 and 1. One asterisk (*) appears to the right of a P-value less than 0.1, two to the right of a P-value less than 0.01, three to the right of a P-value less than 0.001, and four to the right of a P-value less than 0.0001. One caret (^) appears to the right of a P-value greater than 0.9, two to the right of a P-value greater than 0.99, three to the right of a P-value greater than 0.999, and four to the right of a P-value greater than 0.9999. All the statistics calculated are chi-squared statistics. An effort was made to keep all the expected counts in each chi-squared calculation high. Expected counts in all the tests (with a minimum test sequence length of one megabyte) are never lower than 2048. The "cumulative sum" chi-squared statistic for each test is calculated by summing the statistics for all the keys (or sequences) tested. The degrees-of-freedom for the resulting chi-squared statistic is the degrees-of-freedom for the test multiplied by the number of statistics summed. When more than one key or sequence is tested, the Mean Absolute Deviation (MAD) of the P-values for each test is calculated. This is the mean absolute deviation of the P-values from their theoretical median of 0.5. If a large number of keys or sequences is tested, the MAD should be very close to 0.25. A value much higher than 0.25 indicates that too many P-values are far from 0.5. A value much lower than 0.25 indicates that too many P-values are close to 0.5. This MAD has two weaknesses. First, it will look very good if all the P-values are very close to 0.25 and 0.75 (I'm not sure if this is even possible, but it is a weakness); and, second, it is difficult to calculate just how probable a given MAD is. I can give some general guidelines. If the number of sequences or keys tested is at least 256, any MADs greater than 0.30 or less than 0.20 should be viewed with suspicion. If at least 1024 keys or sequences are tested, the MADs should all be 0.24, 0.25, or 0.26. Compile and link this file with Konton2.c to create KonStat2. Optionally you may include Snow 2.0 (I do). See the macro definition and comment after the #includes below. I once made a rather odd discovery about the output of Snow 2.0's PRNG, so I keep trying new tests on it as well as my own cipher and RC4 (which has tiny but real biases). Compile and run the source code in http://www.geocities.com/da5id65536/snowtest.zip to see the oddness for yourself. The source code for Snow 2.0 (snow2.c, snow2.h, and snow2tab.h) is in http://www.it.lth.se/cryptology/snow/snow2.0.zip */ #include #include #include #include #include #include #include #include #include #include #include "Konton2.h" /* If KONSTAT_STD_C is not defined by the following lines, the program may still compile, but warnings that nonstandard C was used will be printed to user. */ #ifdef __STDC__ #if __STDC__ #define KONSTAT_STD_C #endif #endif /* To test Snow2 as well as Konton2 and RC4, define TEST_SNOW2 as 1; else, define it as 0. */ #ifndef TEST_SNOW2 #define TEST_SNOW2 0 #endif #if (TEST_SNOW2 != 0 && TEST_SNOW2 != 1) #error "TEST_SNOW2 must equal 0 or 1." #endif #if TEST_SNOW2 #include "snow2.h" #define kMySnowKeyBits 256 /* must be 256 or 128 */ #endif /* In many implementations of C, the getchar and scanf functions are unsuitable for interactive programs because they wait for the enter (or return) key to be pressed before returning but leave the enter (or return) character in the input queue. This should be considered a bug, and it should be fixed. In any case, if your implementation has this problem, define EAT_THE_ENTER as 1; else, define it as 0. I'd suggest initially defining it as 1; then, if you find that you have to press enter (or return) twice every time you enter something, define it as 0. */ #ifndef EAT_THE_ENTER #define EAT_THE_ENTER 1 #endif #if (EAT_THE_ENTER != 0 && EAT_THE_ENTER != 1) #error "EAT_THE_ENTER must equal 0 or 1." #endif #if (UCHAR_MAX == 0xFF) typedef unsigned char OneByteInt; #else #error "OneByteInt must be unsigned and 1 byte (8 bits) long." #endif #if (USHRT_MAX == 0xFFFFFFFF) typedef unsigned short FourByteInt; #elif (UINT_MAX == 0xFFFFFFFF) typedef unsigned int FourByteInt; #elif (ULONG_MAX >= 0xFFFFFFFF) typedef unsigned long FourByteInt; #else #error "FourByteInt must be unsigned and at least 4 bytes (32 bits) long." #endif #define TO_FILE 1 #define NOT_TO_FILE 0 typedef struct {size_t indexA; size_t indexB; size_t permArray [256];} RC4EnvStruct; typedef struct {clock_t startTime; double accumTime; unsigned char watchMoving;} WatchStruct; int main (int argc, char **argv); void PrepareForExit (void); void TheSignalHandler (int theSignal); char BetterGetChar (void); unsigned char BetterScanULong (unsigned long *theNumber); void PrintAHorzLine (char writeToFile); void PrintToScrnAndFile (char writeToFile, const char *theString, ...); void RunTests (const void *theData, FourByteInt theDataLength); void WriteChiSquaredStatistic (const char *theString, double theDF, double theStat, unsigned long theIndex); void WriteSumStatistic (unsigned long theIndex); void ClearSums (void); void WriteSums (char *generatorName, unsigned long keyLen, double encryptionSpeed); double BitRunsStat (const void *theData, FourByteInt theDataLength, double *theDF); double AltBitRunsStat (const void *theData, FourByteInt theDataLength, double *theDF, FourByteInt theLength); double BlockBitFreqStat (const void *theData, FourByteInt theDataLength, double *theDF); double BitFreqStat (const void *theData, FourByteInt theDataLength, double *theDF); double TidbitFreqStat (const void *theData, FourByteInt theDataLength, double *theDF); double NibbleFreqStat (const void *theData, FourByteInt theDataLength, double *theDF); double ByteFreqStat (const void *theData, FourByteInt theDataLength, double *theDF); double BitStat (const void *theData, FourByteInt theDataLength, double *theDF, OneByteInt theBit); double EightBitSumStat (const void *theData, FourByteInt theDataLength, double *theDF); double SixteenBitSumStat (const void *theData, FourByteInt theDataLength, double *theDF); double ThirtyTwoBitSumStat (const void *theData, FourByteInt theDataLength, double *theDF); double OddEvenBitSumStat (const void *theData, FourByteInt theDataLength, double *theDF); double Matrix64Stat (const void *theData, FourByteInt theDataLength, double *theDF); double Matrix128Stat (const void *theData, FourByteInt theDataLength, double *theDF); double Matrix256Stat (const void *theData, FourByteInt theDataLength, double *theDF); double Matrix512Stat (const void *theData, FourByteInt theDataLength, double *theDF); double Matrix1KStat (const void *theData, FourByteInt theDataLength, double *theDF); double Matrix2KStat (const void *theData, FourByteInt theDataLength, double *theDF); double BitPredictAStat (const void *theData, FourByteInt theDataLength, double *theDF); double BitPredictStat (const void *theData, FourByteInt theDataLength, double *theDF, FourByteInt theLength); double BytePredictAStat (const void *theData, FourByteInt theDataLength, double *theDF); double BytePredictBStat (const void *theData, FourByteInt theDataLength, double *theDF); double BytePredictCStat (const void *theData, FourByteInt theDataLength, double *theDF); double BytePredictDStat (const void *theData, FourByteInt theDataLength, double *theDF); double ByteRepStat (const void *theData, FourByteInt theDataLength, double *theDF); double TwoByteANDStat (const void *theData, FourByteInt theDataLength, double *theDF); double FourByteANDStat (const void *theData, FourByteInt theDataLength, double *theDF); double FourByteAND_ORStat (const void *theData, FourByteInt theDataLength, double *theDF); double TwoByteUpDownStat (const void *theData, FourByteInt theDataLength, double *theDF); double FourByteUpDwnStat (const void *theData, FourByteInt theDataLength, double *theDF); double RectDistStat (const void *theData, FourByteInt theDataLength, double *theDF); double ByteCollisionStat (const void *theData, FourByteInt theDataLength, double *theDF, FourByteInt segLength); double OffsetByteXORStat (const void *theData, FourByteInt theDataLength, double *theDF, FourByteInt theOffset); double ModNStat (const void *theData, FourByteInt theDataLength, double *theDF, OneByteInt theDivisor); FourByteInt BitSum (OneByteInt theByte); FourByteInt GCD (FourByteInt numA, FourByteInt numB); void KeyGenerator (void *theKey, size_t theKeyLength); void ZeroDataBytes (void *theData, size_t theDataLength); long double HexTimesX (char *theString, long double theMultiplier); long double LDFract (long double theX); long double ATimesLog2OfX (long double theA, long double theX, long double *theInt); long double Log2OfX (long double theX, long double *theInt); long double TwoToTheX (long double theX); long double Log2OfGamma (long double theX, long double *theInt); long double IncompleteGammaFactor (long double theA, long double theX); double IncompleteGammaBySeries (double theA, double theX); double IncompleteGammaByFraction (double theA, double theX); double IncompleteGamma (double theA, double theX); double ChiSquaredPValue (double degOfFreedom, double theChiSquared); void ResetStopwatch (WatchStruct *theWatch); void StartStopwatch (WatchStruct *theWatch); void StopStopwatch (WatchStruct *theWatch); double ReadStopwatch (WatchStruct *theWatch); void RC4SetKey (RC4EnvStruct *theEnv, const void *theKey, size_t theKeyLength); void RC4 (RC4EnvStruct *theEnv, void *theData, size_t theDataLength); #if TEST_SNOW2 void Snow2 (void *theData, size_t theDataLength); #endif #define kMyMaxStats 70 RC4EnvStruct gRC4Env; KonEnvStruct gKonEnv; unsigned long gLCGSeed; unsigned long gTheN; unsigned long gNumOfStats; double gDevSums [kMyMaxStats]; double gStatSums [kMyMaxStats]; double gDFSums [kMyMaxStats]; char gStatString [kMyMaxStats][64]; FourByteInt gCountArray [256]; OneByteInt gByteArray [65]; void *gBlock; FILE *gTextFile; FILE *gBinFile; char gAllToFile; int main (int argc, char **argv) {char theOpt; unsigned char evalKonton2, evalKonton2Idx, evalRC4; #if TEST_SNOW2 unsigned char evalSnow2; #endif unsigned long theIndex, theLength, numOfKeys, lcgSeed; double bytesEncrypted; OneByteInt theKey [256]; char dftlTextFileName [] = "KonStat.txt"; char dfltBinFileName [] = "KonStat.bin"; char *textFileName; char *binFileName; WatchStruct theWatch; errno = 0; gBlock = (void *) 0; gTextFile = (FILE *) 0; gBinFile = (FILE *) 0; signal (SIGINT, TheSignalHandler); signal (SIGTERM, TheSignalHandler); textFileName = dftlTextFileName; binFileName = dfltBinFileName; printf ("\n"); /* These are all the ways that the command line can be bad. */ theOpt = 0; if (argc == 4 || argc > 5) theOpt = 1; else if (argc < 3) {if (argc == 2 && strcmp (argv [1], "-h") && strcmp (argv [1], "--help")) theOpt = 1;} else if (strcmp (argv [1], "-o") && strcmp (argv [1], "-i")) theOpt = 1; else if (argc == 5) {if ((strcmp (argv [3], "-o") && strcmp (argv [3], "-i")) || !strcmp (argv [1], argv [3]) || !strcmp (argv [2], argv [4])) theOpt = 1;} if (theOpt || argc == 2) {if (theOpt) printf ("The command-line is inappropriate.\a\n\n"); printf ("KonStat2 command-line arguments:\n\n"); printf ("-h or --help displays this help and exits;\n"); printf (" must be the only argument\n"); printf ("-i specifies the file to be evaluated\n"); printf ("-o specifies the file to which the\n"); printf (" results should be output\n\n"); printf ("* KonStat2 can evaluate a file, or it can evaluate\n"); #if TEST_SNOW2 printf (" stream ciphers: either Konton2, RC4, Snow2, or\n"); printf (" all three. Of course, if ciphers are evaluated,\n"); printf (" it's not necessary to specify a filename with -i.\n"); #else printf (" stream ciphers: either Konton2, RC4, or both.\n"); printf (" Of course, if ciphers are evaluated, it's not\n"); printf (" necessary to specify a filename with -i.\n"); #endif printf ("* The results are always output to the screen;\n"); printf (" specifying a filename with -o is optional.\n"); printf ("* When you don't specify filenames, the default names\n"); printf (" are \"%s\" for -i and \"%s\" for -o.\n\n", binFileName, textFileName); printf ("KonStat2.\n"); printf ("Copyright (C) 2002, 2003, 2004, 2005 by J. David Sexton.\n"); printf ("Revised on Wednesday, February 9, 2005 at 14:30 UTC.\n"); #ifdef __GNUC__ printf ("Compiled with the GNU C compiler %d.%d", __GNUC__, __GNUC_MINOR__); #ifdef __GNUC_PATCHLEVEL__ printf (".%d", __GNUC_PATCHLEVEL__); #endif #ifdef __GNU_LIBRARY__ printf (", library %d", __GNU_LIBRARY__); #endif printf (".\n"); #endif #ifndef KONSTAT_STD_C printf ("\nWARNING: Compiled with nonstandard C!\n"); printf ("(__STDC__ is not defined as required by ISO 9899:1990.)\n\n"); #endif #if defined (__linux__) || defined (__linux) || defined (linux) || defined (__gnu_linux__) printf ("Thank you for using GNU/Linux!\n"); #endif printf ("KonStat2 is not public domain or Open Source,\n"); printf ("but I hereby license you to use it and give it to\n"); printf ("others so long as you don't profit by doing so.\n"); printf ("The tests are documented in the source code file,\n"); printf ("KonStat2.c. More information is available at:\n"); printf ("http://www.geocities.com/da5id65536/\n"); printf ("If I sold this, I would guarantee it; but as I don't,\n"); printf ("I won't.\n\n"); if (theOpt) return (EXIT_FAILURE); return (EXIT_SUCCESS);} if (argc > 2) {if (strcmp (argv [1], "-o")) binFileName = argv [2]; else textFileName = argv [2]; if (argc == 5) {if (strcmp (argv [3], "-o")) binFileName = argv [4]; else textFileName = argv [4];}} printf ("KonStat2.\n"); printf ("Copyright (C) 2002, 2003, 2004, 2005 by J. David Sexton.\n"); #ifndef KONSTAT_STD_C printf ("\nWARNING: Compiled with nonstandard C!\n"); printf ("(__STDC__ is not defined as required by ISO 9899:1990.)\n\n"); #endif printf (" For help, type \"%s -h\" at the command line prompt.\n\n", argv [0]); printf ("1. Evaluate Konton2.\n"); printf ("2. Evaluate Konton2 indexed encryption.\n"); printf ("3. Evaluate RC4.\n"); #if TEST_SNOW2 printf ("4. Evaluate Snow2.\n"); printf ("5. Evaluate Konton2, Konton2 (indexed),\n"); printf (" RC4, and Snow2.\n"); printf ("6. Evaluate the file \"%s.\"\n\n", binFileName); printf ("Enter your choice (1, 2, 3, 4, 5, or 6): "); theOpt = BetterGetChar (); evalKonton2 = (theOpt == '1' || theOpt == '5'); evalKonton2Idx = (theOpt == '2' || theOpt == '5'); evalRC4 = (theOpt == '3' || theOpt == '5'); evalSnow2 = (theOpt == '4' || theOpt == '5'); if (evalKonton2 || evalKonton2Idx || evalRC4 || evalSnow2) theOpt = '1'; else if (theOpt == '6') theOpt = '2'; #else printf ("4. Evaluate Konton2, Konton2 (indexed),\n"); printf (" and RC4.\n"); printf ("5. Evaluate the file \"%s.\"\n\n", binFileName); printf ("Enter your choice (1, 2, 3, 4, or 5): "); theOpt = BetterGetChar (); evalKonton2 = (theOpt == '1' || theOpt == '4'); evalKonton2Idx = (theOpt == '2' || theOpt == '4'); evalRC4 = (theOpt == '3' || theOpt == '4'); if (evalKonton2 || evalKonton2Idx || evalRC4) theOpt = '1'; else if (theOpt == '5') theOpt = '2'; #endif numOfKeys = 0UL; if (theOpt == '1') {printf ("The keys are generated by a truncated\n"); printf ("multiplicative linear congruent generator.\n"); printf ("The value entered below is the seed value\n"); printf ("for this generator.\n\n"); printf ("Enter the key-generating-LCG seed: "); if (BetterScanULong (&lcgSeed)) {printf ("Enter the number of keys to test: "); if (!BetterScanULong (&numOfKeys)) numOfKeys = 0UL;}} else if (theOpt == '2') {if (!(gBinFile = fopen (binFileName, "rb"))) {printf ("\"%s\" could not be opened.\n", binFileName); printf ("System: \"%s\"\n\n", strerror (errno));}} if (!numOfKeys && !gBinFile) {printf ("\a"); return (EXIT_FAILURE);} printf ("Enter the test sequence length\n"); printf ("in megabytes (1 through 4095): "); if (!BetterScanULong (&theLength)) theLength = 0UL; else if (theLength > 4095UL) theLength = 0UL; if (theLength) {theLength <<= 20; if (!(gBlock = malloc (theLength))) {printf ("Memory to run the tests couldn't be allocated.\n"); printf ("System: \"%s\"\n\n", strerror (errno)); theLength = 0UL;}} if (!theLength) {printf ("\a"); PrepareForExit (); return (EXIT_FAILURE);} printf ("1. Output all the results only to the screen.\n"); printf ("2. Output all the results to the screen and\n"); printf (" to the file \"%s.\"\n", textFileName); printf ("3. Output all the results to the screen and\n"); printf (" only the cumulative sums to the file\n"); printf (" \"%s.\"\n\n", textFileName); printf ("Enter your choice (1, 2, or 3): "); theOpt = BetterGetChar (); if (theOpt == '1') gAllToFile = NOT_TO_FILE; else if (theOpt == '2' || theOpt == '3') {if ((gTextFile = fopen (textFileName, "w"))) {if (theOpt == '2') gAllToFile = TO_FILE; else gAllToFile = NOT_TO_FILE;} else {printf ("\"%s\" could not be opened.\n", textFileName); printf ("System: \"%s\"\n\n", strerror (errno)); theOpt = 0;}} else theOpt = 0; if (!theOpt) {printf ("\a"); PrepareForExit (); return (EXIT_FAILURE);} gNumOfStats = 0UL; if (!gBinFile) PrintToScrnAndFile (TO_FILE, "the key-generating-LCG seed = %lu\n", lcgSeed); PrintToScrnAndFile (TO_FILE, "the test sequence length = %lu MB\n", theLength >> 20); PrintAHorzLine (TO_FILE); if (gBinFile) {ClearSums (); theIndex = 0; while (fread (gBlock, 1, theLength, gBinFile) == theLength) {PrintToScrnAndFile (gAllToFile, "\n %s, sequence %lu\n", binFileName, ++theIndex); RunTests (gBlock, theLength);} fclose (gBinFile); gBinFile = (FILE *) 0; if (theIndex) WriteSums (binFileName, 0UL, 0.0); else {PrintToScrnAndFile (TO_FILE, "\n** \"%s\" appears to be less than %lu MB long.\n", binFileName, theLength >> 20); printf ("\a"); PrintAHorzLine (TO_FILE);}} else {bytesEncrypted = theLength; bytesEncrypted *= numOfKeys; if (evalKonton2) {gLCGSeed = lcgSeed; ResetStopwatch (&theWatch); ClearSums (); theIndex = 0UL; do {KeyGenerator (theKey, (size_t) kKonCryptKeyLength); KonSetCryptKey (&gKonEnv, theKey, (size_t) kKonCryptKeyLength); ZeroDataBytes (gBlock, theLength); PrintToScrnAndFile (gAllToFile, "\n"); StartStopwatch (&theWatch); KonCryptSynch (&gKonEnv, gBlock, theLength); StopStopwatch (&theWatch); PrintToScrnAndFile (gAllToFile, " Konton2, (%lu-byte) key %lu\n", (unsigned long) kKonCryptKeyLength, ++theIndex); RunTests (gBlock, theLength);} while (theIndex < numOfKeys); WriteSums ("Konton2", (unsigned long) kKonCryptKeyLength, bytesEncrypted / ReadStopwatch (&theWatch));} if (evalKonton2Idx) {gLCGSeed = lcgSeed; ResetStopwatch (&theWatch); ClearSums (); theIndex = 0UL; do {KeyGenerator (theKey, (size_t) kKonCryptKeyLength); KonSetCryptKey (&gKonEnv, theKey, (size_t) kKonCryptKeyLength); ZeroDataBytes (gBlock, theLength); PrintToScrnAndFile (gAllToFile, "\n"); StartStopwatch (&theWatch); KonCryptSynch (&gKonEnv, gBlock, theLength); StopStopwatch (&theWatch); KonIndexState (&gKonEnv, ++theIndex); StartStopwatch (&theWatch); KonCryptSynch (&gKonEnv, gBlock, theLength); StopStopwatch (&theWatch); PrintToScrnAndFile (gAllToFile, " Konton2 (indexed), (%lu-byte) key %lu\n", (unsigned long) kKonCryptKeyLength, theIndex); RunTests (gBlock, theLength);} while (theIndex < numOfKeys); WriteSums ("Konton2 (indexed)", (unsigned long) kKonCryptKeyLength, (bytesEncrypted * 2.0) / ReadStopwatch (&theWatch));} if (evalRC4) {gLCGSeed = lcgSeed; ResetStopwatch (&theWatch); ClearSums (); theIndex = 0UL; do {KeyGenerator (theKey, sizeof (theKey)); RC4SetKey (&gRC4Env, theKey, sizeof (theKey)); ZeroDataBytes (gBlock, theLength); PrintToScrnAndFile (gAllToFile, "\n"); StartStopwatch (&theWatch); RC4 (&gRC4Env, gBlock, theLength); StopStopwatch (&theWatch); PrintToScrnAndFile (gAllToFile, " RC4, (%lu-byte) key %lu\n", (unsigned long) sizeof (theKey), ++theIndex); RunTests (gBlock, theLength);} while (theIndex < numOfKeys); WriteSums ("RC4", sizeof (theKey), bytesEncrypted / ReadStopwatch (&theWatch));} #if TEST_SNOW2 if (evalSnow2) {gLCGSeed = lcgSeed; ResetStopwatch (&theWatch); ClearSums (); theIndex = 0UL; do {KeyGenerator (theKey, (size_t) kMySnowKeyBits / 8); snow_loadkey (theKey, kMySnowKeyBits, 0UL, 0UL, 0UL, 0UL); ZeroDataBytes (gBlock, theLength); PrintToScrnAndFile (gAllToFile, "\n"); StartStopwatch (&theWatch); Snow2 (gBlock, theLength); StopStopwatch (&theWatch); PrintToScrnAndFile (gAllToFile, " Snow2, (%lu-byte) key %lu\n", kMySnowKeyBits / 8UL, ++theIndex); RunTests (gBlock, theLength);} while (theIndex < numOfKeys); WriteSums ("Snow2", kMySnowKeyBits / 8UL, bytesEncrypted / ReadStopwatch (&theWatch));} #endif } PrintToScrnAndFile (TO_FILE, "\n"); printf ("\n"); PrepareForExit (); return (EXIT_SUCCESS);} void PrepareForExit (void) {if (gBlock) free (gBlock); if (gBinFile) fclose (gBinFile); if (gTextFile) {fflush (gTextFile); fclose (gTextFile);} return;} void TheSignalHandler (int theSignal) {PrepareForExit (); exit (EXIT_SUCCESS);} char BetterGetChar (void) {char outputChar; outputChar = getchar (); #if EAT_THE_ENTER getchar (); #endif printf ("\n"); return (outputChar);} unsigned char BetterScanULong (unsigned long *theNumber) {int argsAssigned; argsAssigned = scanf ("%lu", theNumber); #if EAT_THE_ENTER getchar (); #endif printf ("\n"); return (argsAssigned == 1);} void PrintAHorzLine (char writeToFile) {PrintToScrnAndFile (writeToFile, "---------------------------------------"); PrintToScrnAndFile (writeToFile, "---------------------------------------"); if (writeToFile == TO_FILE && gTextFile) fflush (gTextFile); return;} void PrintToScrnAndFile (char writeToFile, const char *theString, ...) {va_list theVAList; va_start (theVAList, theString); vprintf (theString, theVAList); if (writeToFile == TO_FILE && gTextFile) {if (vfprintf (gTextFile, theString, theVAList) < 0) {fclose (gTextFile); gTextFile = (FILE *) 0L;}} va_end (theVAList); return;} void RunTests (const void *theData, FourByteInt theDataLength) {unsigned long theIndex; double theStat, theDF; double cumStat, cumDF; ++gTheN; theIndex = 0UL; theStat = BitRunsStat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the bit runs statistic ", theDF, theStat, theIndex++); theStat = AltBitRunsStat (theData, theDataLength, &theDF, 1); WriteChiSquaredStatistic ("the alt. 1 bit runs statistic ", theDF, theStat, theIndex++); theStat = AltBitRunsStat (theData, theDataLength, &theDF, 2); WriteChiSquaredStatistic ("the alt. 2 bit runs statistic ", theDF, theStat, theIndex++); theStat = AltBitRunsStat (theData, theDataLength, &theDF, 3); WriteChiSquaredStatistic ("the alt. 3 bit runs statistic ", theDF, theStat, theIndex++); theStat = AltBitRunsStat (theData, theDataLength, &theDF, 4); WriteChiSquaredStatistic ("the alt. 4 bit runs statistic ", theDF, theStat, theIndex++); theStat = AltBitRunsStat (theData, theDataLength, &theDF, 5); WriteChiSquaredStatistic ("the alt. 5 bit runs statistic ", theDF, theStat, theIndex++); theStat = BlockBitFreqStat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the block bit freq. statistic ", theDF, theStat, theIndex++); theStat = BitFreqStat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the bit frequency statistic ", theDF, theStat, theIndex++); theStat = TidbitFreqStat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the tidbit frequency statistic ", theDF, theStat, theIndex++); theStat = NibbleFreqStat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the nibble frequency statistic ", theDF, theStat, theIndex++); theStat = ByteFreqStat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the byte frequency statistic ", theDF, theStat, theIndex++); theStat = BitStat (theData, theDataLength, &theDF, 0); WriteChiSquaredStatistic ("the bit 0 statistic ", theDF, theStat, theIndex++); cumStat = theStat; cumDF = theDF; theStat = BitStat (theData, theDataLength, &theDF, 1); WriteChiSquaredStatistic ("the bit 1 statistic ", theDF, theStat, theIndex++); cumStat += theStat; cumDF += theDF; theStat = BitStat (theData, theDataLength, &theDF, 2); WriteChiSquaredStatistic ("the bit 2 statistic ", theDF, theStat, theIndex++); cumStat += theStat; cumDF += theDF; theStat = BitStat (theData, theDataLength, &theDF, 3); WriteChiSquaredStatistic ("the bit 3 statistic ", theDF, theStat, theIndex++); cumStat += theStat; cumDF += theDF; theStat = BitStat (theData, theDataLength, &theDF, 4); WriteChiSquaredStatistic ("the bit 4 statistic ", theDF, theStat, theIndex++); cumStat += theStat; cumDF += theDF; theStat = BitStat (theData, theDataLength, &theDF, 5); WriteChiSquaredStatistic ("the bit 5 statistic ", theDF, theStat, theIndex++); cumStat += theStat; cumDF += theDF; theStat = BitStat (theData, theDataLength, &theDF, 6); WriteChiSquaredStatistic ("the bit 6 statistic ", theDF, theStat, theIndex++); cumStat += theStat; cumDF += theDF; theStat = BitStat (theData, theDataLength, &theDF, 7); WriteChiSquaredStatistic ("the bit 7 statistic ", theDF, theStat, theIndex++); cumStat += theStat; cumDF += theDF; WriteChiSquaredStatistic ("the overall bit statistic ", cumDF, cumStat, theIndex++); theStat = EightBitSumStat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the 8-bit-sum statistic ", theDF, theStat, theIndex++); theStat = SixteenBitSumStat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the 16-bit-sum statistic ", theDF, theStat, theIndex++); theStat = ThirtyTwoBitSumStat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the 32-bit-sum statistic ", theDF, theStat, theIndex++); theStat = OddEvenBitSumStat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the odd/even-bit-sum statistic ", theDF, theStat, theIndex++); theStat = Matrix64Stat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the 64-bit matrix statistic ", theDF, theStat, theIndex++); theStat = Matrix128Stat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the 128-bit matrix statistic ", theDF, theStat, theIndex++); theStat = Matrix256Stat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the 256-bit matrix statistic ", theDF, theStat, theIndex++); theStat = Matrix512Stat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the 512-bit matrix statistic ", theDF, theStat, theIndex++); theStat = Matrix1KStat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the 1024-bit matrix statistic ", theDF, theStat, theIndex++); theStat = Matrix2KStat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the 2048-bit matrix statistic ", theDF, theStat, theIndex++); theStat = BitPredictAStat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the bit prediction A statistic ", theDF, theStat, theIndex++); theStat = BitPredictStat (theData, theDataLength, &theDF, 4); WriteChiSquaredStatistic ("the bit prediction B statistic ", theDF, theStat, theIndex++); theStat = BitPredictStat (theData, theDataLength, &theDF, 8); WriteChiSquaredStatistic ("the bit prediction C statistic ", theDF, theStat, theIndex++); theStat = BitPredictStat (theData, theDataLength, &theDF, 16); WriteChiSquaredStatistic ("the bit prediction D statistic ", theDF, theStat, theIndex++); theStat = BitPredictStat (theData, theDataLength, &theDF, 32); WriteChiSquaredStatistic ("the bit prediction E statistic ", theDF, theStat, theIndex++); theStat = BytePredictAStat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the byte prediction A statistic", theDF, theStat, theIndex++); theStat = BytePredictBStat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the byte prediction B statistic", theDF, theStat, theIndex++); theStat = BytePredictCStat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the byte prediction C statistic", theDF, theStat, theIndex++); theStat = BytePredictDStat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the byte prediction D statistic", theDF, theStat, theIndex++); theStat = ByteRepStat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the byte repetition statistic ", theDF, theStat, theIndex++); theStat = TwoByteANDStat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the 2-byte AND statistic ", theDF, theStat, theIndex++); theStat = FourByteANDStat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the 4-byte AND statistic ", theDF, theStat, theIndex++); theStat = FourByteAND_ORStat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the 4-byte AND-OR statistic ", theDF, theStat, theIndex++); theStat = TwoByteUpDownStat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the 2-byte up/down statistic ", theDF, theStat, theIndex++); theStat = FourByteUpDwnStat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the 4-byte up/down statistic ", theDF, theStat, theIndex++); theStat = RectDistStat (theData, theDataLength, &theDF); WriteChiSquaredStatistic ("the rect. distance statistic ", theDF, theStat, theIndex++); theStat = ByteCollisionStat (theData, theDataLength, &theDF, 4); WriteChiSquaredStatistic ("the 4-byte collision statistic ", theDF, theStat, theIndex++); theStat = ByteCollisionStat (theData, theDataLength, &theDF, 6); WriteChiSquaredStatistic ("the 6-byte collision statistic ", theDF, theStat, theIndex++); theStat = ByteCollisionStat (theData, theDataLength, &theDF, 8); WriteChiSquaredStatistic ("the 8-byte collision statistic ", theDF, theStat, theIndex++); theStat = ByteCollisionStat (theData, theDataLength, &theDF, 12); WriteChiSquaredStatistic ("the 12-byte collision statistic", theDF, theStat, theIndex++); theStat = ByteCollisionStat (theData, theDataLength, &theDF, 16); WriteChiSquaredStatistic ("the 16-byte collision statistic", theDF, theStat, theIndex++); theStat = ByteCollisionStat (theData, theDataLength, &theDF, 24); WriteChiSquaredStatistic ("the 24-byte collision statistic", theDF, theStat, theIndex++); theStat = ByteCollisionStat (theData, theDataLength, &theDF, 32); WriteChiSquaredStatistic ("the 32-byte collision statistic", theDF, theStat, theIndex++); theStat = OffsetByteXORStat (theData, theDataLength, &theDF, 1); WriteChiSquaredStatistic ("the 1-byte offset XOR statistic", theDF, theStat, theIndex++); theStat = OffsetByteXORStat (theData, theDataLength, &theDF, 2); WriteChiSquaredStatistic ("the 2-byte offset XOR statistic", theDF, theStat, theIndex++); theStat = OffsetByteXORStat (theData, theDataLength, &theDF, 6); WriteChiSquaredStatistic ("the 6-byte offset XOR statistic", theDF, theStat, theIndex++); theStat = OffsetByteXORStat (theData, theDataLength, &theDF, 24); WriteChiSquaredStatistic ("the 24-byte offset XOR stat. ", theDF, theStat, theIndex++); theStat = OffsetByteXORStat (theData, theDataLength, &theDF, 120); WriteChiSquaredStatistic ("the 120-byte offset XOR stat. ", theDF, theStat, theIndex++); theStat = OffsetByteXORStat (theData, theDataLength, &theDF, 720); WriteChiSquaredStatistic ("the 720-byte offset XOR stat. ", theDF, theStat, theIndex++); theStat = OffsetByteXORStat (theData, theDataLength, &theDF, 5040); WriteChiSquaredStatistic ("the 5040-byte offset XOR stat. ", theDF, theStat, theIndex++); theStat = OffsetByteXORStat (theData, theDataLength, &theDF, 40320); WriteChiSquaredStatistic ("the 40320-byte offset XOR stat.", theDF, theStat, theIndex++); theStat = OffsetByteXORStat (theData, theDataLength, &theDF, 362880); WriteChiSquaredStatistic ("the 362880-byte off. XOR stat. ", theDF, theStat, theIndex++); theStat = ModNStat (theData, theDataLength, &theDF, 3); WriteChiSquaredStatistic ("the mod 3 statistic ", theDF, theStat, theIndex++); theStat = ModNStat (theData, theDataLength, &theDF, 5); WriteChiSquaredStatistic ("the mod 5 statistic ", theDF, theStat, theIndex++); theStat = ModNStat (theData, theDataLength, &theDF, 7); WriteChiSquaredStatistic ("the mod 7 statistic ", theDF, theStat, theIndex++); theStat = ModNStat (theData, theDataLength, &theDF, 11); WriteChiSquaredStatistic ("the mod 11 statistic ", theDF, theStat, theIndex++); theStat = ModNStat (theData, theDataLength, &theDF, 13); WriteChiSquaredStatistic ("the mod 13 statistic ", theDF, theStat, theIndex++); theStat = ModNStat (theData, theDataLength, &theDF, 17); WriteChiSquaredStatistic ("the mod 17 statistic ", theDF, theStat, theIndex++); theStat = ModNStat (theData, theDataLength, &theDF, 19); WriteChiSquaredStatistic ("the mod 19 statistic ", theDF, theStat, theIndex++); theStat = ModNStat (theData, theDataLength, &theDF, 23); WriteChiSquaredStatistic ("the mod 23 statistic ", theDF, theStat, theIndex++); PrintAHorzLine (gAllToFile); return;} void WriteChiSquaredStatistic (const char *theString, double theDF, double theStat, unsigned long theIndex) {double thePValue; if (theIndex == gNumOfStats) {if (gNumOfStats == kMyMaxStats) return; ++gNumOfStats; strncpy (gStatString [theIndex], theString, sizeof (gStatString [0]) - 1); gStatString [theIndex][sizeof (gStatString [0]) - 1] = '\0';} thePValue = ChiSquaredPValue (theDF, theStat); gDevSums [theIndex] += (thePValue < 0.5 ? 0.5 - thePValue : thePValue - 0.5); gStatSums [theIndex] += theStat; gDFSums [theIndex] += theDF; PrintToScrnAndFile (gAllToFile, "%s\t= %10.5f [P = %6.4f]", gStatString [theIndex], theStat, thePValue); if (thePValue < 0.1) {PrintToScrnAndFile (gAllToFile, "*"); if (thePValue < 0.01) {PrintToScrnAndFile (gAllToFile, "*"); if (thePValue < 0.001) {PrintToScrnAndFile (gAllToFile, "*"); if (thePValue < 0.0001) PrintToScrnAndFile (gAllToFile, "*");}}} else if (thePValue > 0.9) {PrintToScrnAndFile (gAllToFile, "^"); if (thePValue > 0.99) {PrintToScrnAndFile (gAllToFile, "^"); if (thePValue > 0.999) {PrintToScrnAndFile (gAllToFile, "^"); if (thePValue > 0.9999) PrintToScrnAndFile (gAllToFile, "^");}}} PrintToScrnAndFile (gAllToFile, "\n"); return;} void WriteSumStatistic (unsigned long theIndex) {double thePValue; thePValue = ChiSquaredPValue (gDFSums [theIndex], gStatSums [theIndex]); PrintToScrnAndFile (TO_FILE, "%s\t= %14.5f ", gStatString [theIndex], gStatSums [theIndex]); if (gTheN > 1UL) PrintToScrnAndFile (TO_FILE, "[MAD = %4.2f]", gDevSums [theIndex] / gTheN); PrintToScrnAndFile (TO_FILE, "[P = %6.4f]", thePValue); if (thePValue < 0.1) {PrintToScrnAndFile (TO_FILE, "*"); if (thePValue < 0.01) {PrintToScrnAndFile (TO_FILE, "*"); if (thePValue < 0.001) {PrintToScrnAndFile (TO_FILE, "*"); if (thePValue < 0.0001) PrintToScrnAndFile (TO_FILE, "*");}}} else if (thePValue > 0.9) {PrintToScrnAndFile (TO_FILE, "^"); if (thePValue > 0.99) {PrintToScrnAndFile (TO_FILE, "^"); if (thePValue > 0.999) {PrintToScrnAndFile (TO_FILE, "^"); if (thePValue > 0.9999) PrintToScrnAndFile (TO_FILE, "^");}}} PrintToScrnAndFile (TO_FILE, "\n"); return;} void ClearSums (void) {unsigned long theIndex; theIndex = 0UL; do {gDevSums [theIndex] = 0.0; gStatSums [theIndex] = 0.0; gDFSums [theIndex] = 0.0;} while (++theIndex < kMyMaxStats); gTheN = 0UL; return;} void WriteSums (char *generatorName, unsigned long keyLen, double encryptionSpeed) {unsigned long theIndex; PrintToScrnAndFile (TO_FILE, "\n %s cumulative sums, %lu ", generatorName, gTheN); if (keyLen) {PrintToScrnAndFile (TO_FILE, "(%lu-byte) ", keyLen); if (gTheN == 1UL) PrintToScrnAndFile (TO_FILE, "key"); else PrintToScrnAndFile (TO_FILE, "keys");} else if (gTheN == 1UL) PrintToScrnAndFile (TO_FILE, "sequence"); else PrintToScrnAndFile (TO_FILE, "sequences"); if (encryptionSpeed > 0.0) PrintToScrnAndFile (TO_FILE, "\n speed = %.2f KB/second", encryptionSpeed / 1024.0); PrintToScrnAndFile (TO_FILE, "\n"); theIndex = 0UL; do WriteSumStatistic (theIndex); while (++theIndex < gNumOfStats); PrintAHorzLine (TO_FILE); return;} /* A run is one or more consecutive bits with the same value (either 0 or 1). Whenever a bit has a value different from the preceding bit, a new run is begun. For example, 01001011010100 consists of eleven runs: 0, 1, 00, 1, 0, 11, 0, 1, 0, 1, and 00. In a random sequence, approximately half the runs should be one bit long, a fourth should be two bits long, an eighth should be three bits long, etc. Runs of 1 to 10 bits are counted. Runs of 11 bits and longer are lumped together. A chi-squared statistic is calculated. The degrees-of-freedom is 10. */ double BitRunsStat (const void *theData, FourByteInt theDataLength, double *theDF) {OneByteInt theByte, currentBit; const OneByteInt *theBytePtr; FourByteInt byteIndex, bitIndex, runLength, theAdjLength; double theStat, theExpected, theDiff; /* gCountArray [0] stores the total number of runs. gCountArray [1]...gCountArray [10] store the numbers of runs of length 1...10, respectively. */ *theDF = 10.0; if (!theDataLength) return (0.0); theAdjLength = theDataLength; if (theAdjLength > 0x1FFFFFFF) theAdjLength = 0x1FFFFFFF; byteIndex = 11; do gCountArray [--byteIndex] = 0; while (byteIndex); theBytePtr = theData; runLength = 0; byteIndex = theAdjLength; currentBit = (*theBytePtr & 0x01); do {theByte = *(theBytePtr++); bitIndex = 8; do {if ((theByte & 0x01) == currentBit) ++runLength; else {++gCountArray [0]; if (runLength < 11) ++gCountArray [runLength]; runLength = 1; currentBit ^= 0x01;} theByte >>= 1;} while (--bitIndex);} while (--byteIndex); ++gCountArray [0]; if (runLength < 11) ++gCountArray [runLength]; theStat = 0.0; runLength = gCountArray [0]; theExpected = runLength; bitIndex = 1; do {theExpected /= 2.0; byteIndex = gCountArray [bitIndex]; theDiff = byteIndex; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; runLength -= byteIndex;} while (++bitIndex < 11); theDiff = runLength; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; return (theStat);} /* A run is one or more consecutive bits that alternate in value between 0 and 1 in a regular pattern. The pattern is specified by theLength. If theLength is 1, the pattern is 010101... or its complement 101010...; if theLength is 2, the pattern is 001100110011... or its complement 110011001100...; if theLength is 3, the pattern is 000111000111... or its complement 111000111000...; and so on. Whenever the pattern is not followed, a new run is begun. For example, if theLength is 1, the sequence 01001011010100 consists of the following four runs: 010, 0101, 101010, and 0. If theLength is 2, the same sequence consists of the following nine runs: 0, 1, 001, 0, 110, 1, 0, 1, and 00. In a random sequence, approximately half the runs should be one bit long, a fourth should be two bits long, an eighth should be three bits long, etc. Runs of 1 to 10 bits are counted. Runs of 11 bits and longer are lumped together. A chi-squared statistic is calculated. The degrees-of-freedom is 10. */ double AltBitRunsStat (const void *theData, FourByteInt theDataLength, double *theDF, FourByteInt theLength) {OneByteInt theByte, currentBit; const OneByteInt *theBytePtr; FourByteInt byteIndex, bitIndex, runLength, theAdjLength, xIndex; double theStat, theExpected, theDiff; /* gCountArray [0] stores the total number of runs. gCountArray [1]...gCountArray [10] store the numbers of runs of length 1...10, respectively. */ *theDF = 10.0; if (!theLength || !theDataLength) return (0.0); theAdjLength = theDataLength; if (theAdjLength > 0x1FFFFFFF) theAdjLength = 0x1FFFFFFF; byteIndex = 11; do gCountArray [--byteIndex] = 0; while (byteIndex); theBytePtr = theData; runLength = 0; byteIndex = theAdjLength; currentBit = (*theBytePtr & 0x01) ^ 0x01; xIndex = 1; do {theByte = *(theBytePtr++); bitIndex = 8; do {if (!(--xIndex)) {currentBit ^= 0x01; xIndex = theLength;} if ((theByte & 0x01) == currentBit) ++runLength; else {++gCountArray [0]; if (runLength < 11) ++gCountArray [runLength]; runLength = 1; currentBit ^= 0x01; xIndex = theLength;} theByte >>= 1;} while (--bitIndex);} while (--byteIndex); ++gCountArray [0]; if (runLength < 11) ++gCountArray [runLength]; theStat = 0.0; runLength = gCountArray [0]; theExpected = runLength; bitIndex = 1; do {theExpected /= 2.0; byteIndex = gCountArray [bitIndex]; theDiff = byteIndex; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; runLength -= byteIndex;} while (++bitIndex < 11); theDiff = runLength; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; return (theStat);} /* The sequence is divided into 256 blocks of equal length. The total number of bits of each possible value in each block is counted. In a random sequence the probability of each possible bit value is 1/2. A chi-squared statistic is calculated. The degrees-of-freedom is 256. */ double BlockBitFreqStat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt theAdjLength; FourByteInt blockLength, blockIndex, theBlock, theTotal; const OneByteInt *theBytePtr; double theStat, theExpected, theDiff; *theDF = 256.0; if (!(theAdjLength = theDataLength & 0xFFFFFF00)) return (0.0); theBytePtr = theData; blockLength = theAdjLength / 256; theBlock = 256; do {theTotal = 0; blockIndex = blockLength; do theTotal += BitSum (*(theBytePtr++)); while (--blockIndex); gCountArray [--theBlock] = theTotal;} while (theBlock); theStat = 0.0; theExpected = blockLength * 4; theBlock = 256; do {theDiff = gCountArray [--theBlock]; theDiff -= theExpected; theDiff *= theDiff; theStat += theDiff;} while (theBlock); theStat /= theExpected; theStat *= 2.0; return (theStat);} /* The total number of bits of each possible value is counted. In a random sequence the probability of each possible bit value is 1/2. A chi-squared statistic is calculated. The degrees-of-freedom is 1. */ double BitFreqStat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt theIndex, theTotal, theAdjLength; const OneByteInt *theBytePtr; double theStat, theExpected; *theDF = 1.0; if (!(theAdjLength = theDataLength)) return (0.0); if (theAdjLength > 0x1FFFFFFF) theAdjLength = 0x1FFFFFFF; theTotal = 0; theIndex = theAdjLength; theBytePtr = theData; do theTotal += BitSum (*(theBytePtr++)); while (--theIndex); theExpected = theAdjLength * 4; theStat = theTotal; theStat -= theExpected; theStat *= theStat; theStat /= theExpected; theStat *= 2.0; return (theStat);} /* The total number of tidbits of each possible value is counted. In a random sequence the probability of each possible tidbit value is 1/4. A chi-squared statistic is calculated. The degrees-of-freedom is 3. */ double TidbitFreqStat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt theIndex, theAdjLength; OneByteInt theByte; const OneByteInt *theBytePtr; double theStat, theExpected, theDiff; *theDF = 3.0; if (!(theAdjLength = theDataLength)) return (0.0); if (theAdjLength > 0x3FFFFFFF) theAdjLength = 0x3FFFFFFF; theIndex = 4; do gCountArray [--theIndex] = 0; while (theIndex); theBytePtr = theData; theIndex = theAdjLength; do {theByte = *(theBytePtr++); ++gCountArray [theByte & 0x03]; theByte >>= 2; ++gCountArray [theByte & 0x03]; theByte >>= 2; ++gCountArray [theByte & 0x03]; theByte >>= 2; ++gCountArray [theByte];} while (--theIndex); theStat = 0.0; theExpected = theAdjLength; theIndex = 4; do {theDiff = gCountArray [--theIndex]; theDiff -= theExpected; theDiff *= theDiff; theStat += theDiff;} while (theIndex); theStat /= theExpected; return (theStat);} /* The total number of nibbles of each possible value is counted. In a random sequence the probability of each possible nibble value is 1/16. A chi-squared statistic is calculated. The degrees-of-freedom is 15. */ double NibbleFreqStat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt theIndex, theAdjLength; OneByteInt theByte; const OneByteInt *theBytePtr; double theStat, theExpected, theDiff; *theDF = 15.0; if (!(theAdjLength = theDataLength & 0xFFFFFFF8)) return (0.0); if (theAdjLength > 0x7FFFFFF8) theAdjLength = 0x7FFFFFF8; theIndex = 16; do gCountArray [--theIndex] = 0; while (theIndex); theBytePtr = theData; theIndex = theAdjLength; do {theByte = *(theBytePtr++); ++gCountArray [theByte & 0x0F]; ++gCountArray [theByte >> 4];} while (--theIndex); theStat = 0.0; theExpected = theAdjLength / 8; theIndex = 16; do {theDiff = gCountArray [--theIndex]; theDiff -= theExpected; theDiff *= theDiff; theStat += theDiff;} while (theIndex); theStat /= theExpected; return (theStat);} /* The total number of bytes of each possible value is counted. In a random sequence the probability of each possible byte value is 1/256. A chi-squared statistic is calculated. The degrees-of-freedom is 255. */ double ByteFreqStat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt theIndex, theAdjLength; const OneByteInt *theBytePtr; double theStat, theExpected, theDiff; *theDF = 255.0; if (!(theAdjLength = theDataLength & 0xFFFFFF00)) return (0.0); theIndex = 256; do gCountArray [--theIndex] = 0; while (theIndex); theBytePtr = theData; theIndex = theAdjLength; do ++gCountArray [*(theBytePtr++)]; while (--theIndex); theStat = 0.0; theExpected = theAdjLength / 256; theIndex = 256; do {theDiff = gCountArray [--theIndex]; theDiff -= theExpected; theDiff *= theDiff; theStat += theDiff;} while (theIndex); theStat /= theExpected; return (theStat);} /* One bit of each byte of the sequence is evaluated in this test. The bit is specified by theBit, which will be a number from 0 to 7: where 0 specifies the lowest-order bit, and 7 specifies the highest-order bit. The value of the specified bit of each byte in the sequence will be either 0 or 1. In a random sequence the probability of each possible bit value is 1/2. A chi-squared statistic is calculated. The degrees-of-freedom is 1. The bit test described above is performed for each of the 8 bit positions. The resulting chi-squared statistics for each of the 8 tests are summed; the sum is an overall chi-squared statistic. The degrees-of-freedom is 8. */ double BitStat (const void *theData, FourByteInt theDataLength, double *theDF, OneByteInt theBit) {FourByteInt theIndex, theTotal, theAdjLength; const OneByteInt *theBytePtr; OneByteInt bitMask; double theStat, theExpected; *theDF = 1.0; if (theBit > 7) return (0.0); if (!(theAdjLength = theDataLength & 0xFFFFFFFE)) return (0.0); bitMask = (0x01 << theBit); theTotal = 0; theIndex = theAdjLength; theBytePtr = theData; do {if (*(theBytePtr++) & bitMask) ++theTotal;} while (--theIndex); theExpected = theAdjLength / 2; theStat = theTotal; theStat -= theExpected; theStat *= theStat; theStat /= theExpected; theStat *= 2.0; return (theStat);} /* The bits of each byte of the sequence are summed. The bit sum will be either 0, 1, 2, 3, 4, 5, 6, 7, or 8. In a random sequence, the probabilities of 0, 1, 2, 3, 4, 5, 6, 7, and 8 are 1/256, 8/256, 28/256, 56/256, 70/256, 56/256, 28/256, 8/256, and 1/256 respectively. The total number of bytes with each possible bit sum is counted. A chi-squared statistic is calculated. The degrees-of-freedom is 8. */ double EightBitSumStat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt theIndex, theAdjLength; const OneByteInt *theBytePtr; double theStat, theExpected, theDiff; *theDF = 8.0; if (!(theAdjLength = theDataLength & 0xFFFFFF00)) return (0.0); theIndex = 9; do gCountArray [--theIndex] = 0; while (theIndex); theBytePtr = theData; theIndex = theAdjLength; do ++gCountArray [BitSum (*(theBytePtr++))]; while (--theIndex); theExpected = theAdjLength / 256; theDiff = gCountArray [0]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat = theDiff; theDiff = gCountArray [8]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theExpected = (theAdjLength / 256) * 8; theDiff = gCountArray [1]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theDiff = gCountArray [7]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theExpected = (theAdjLength / 256) * 28; theDiff = gCountArray [2]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theDiff = gCountArray [6]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theExpected = (theAdjLength / 256) * 56; theDiff = gCountArray [3]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theDiff = gCountArray [5]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theExpected = (theAdjLength / 256) * 70; theDiff = gCountArray [4]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; return (theStat);} /* The sequence is divided into segments of 16 bits. The bits of each segment are summed. The bit sum (s) will fall within one of three ranges: [s<7], [s>9], or [7<=s<=9]. In a random sequence, the probabilities of [s<7] and [s>9] are both 14893/65536, and the probability of [7<=s<=9] is 35750/65536. The total number of segments with bit sums in each of the 3 ranges is counted. A chi-squared statistic is calculated. The degrees-of-freedom is 2. */ double SixteenBitSumStat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt theIndex, theAdjLength, theBitSum; const OneByteInt *theBytePtr; double theStat, theExpected, theDiff; *theDF = 2.0; if (!(theAdjLength = theDataLength & 0xFFFFFF00)) return (0.0); theIndex = 17; do gCountArray [--theIndex] = 0; while (theIndex); theBytePtr = theData; theIndex = theAdjLength / 2; do {theBitSum = BitSum (*(theBytePtr++)); theBitSum += BitSum (*(theBytePtr++)); ++gCountArray [theBitSum];} while (--theIndex); theExpected = theAdjLength / 2; theExpected *= 35750.0 / 65536.0; theDiff = gCountArray [7] + gCountArray [8] + gCountArray [9]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat = theDiff; theDiff = theAdjLength / 2; theExpected = (theDiff - theExpected) / 2.0; theDiff = gCountArray [0] + gCountArray [1] + gCountArray [2] + gCountArray [3] + gCountArray [4] + gCountArray [5] + gCountArray [6]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theDiff = gCountArray [10] + gCountArray [11] + gCountArray [12] + gCountArray [13] + gCountArray [14] + gCountArray [15] + gCountArray [16]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; return (theStat);} /* The sequence is divided into segments of 32 bits. The bits of each segment are summed. The bit sum (s) will fall within one of three ranges: [s<15], [s>17], or [15<=s<=17]. In a random sequence, the probabilities of [s<15] and [s>17] are both 1281220733/4294967296, and the probability of [15<=s<=17] is 1732525830/4294967296. The total number of segments with bit sums in each of the 3 ranges is counted. A chi-squared statistic is calculated. The degrees-of-freedom is 2. */ double ThirtyTwoBitSumStat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt theIndex, theAdjLength, theBitSum; const OneByteInt *theBytePtr; double theStat, theExpected, theDiff; *theDF = 2.0; if (!(theAdjLength = theDataLength & 0xFFFFFF00)) return (0.0); theIndex = 33; do gCountArray [--theIndex] = 0; while (theIndex); theBytePtr = theData; theIndex = theAdjLength / 4; do {theBitSum = BitSum (*(theBytePtr++)); theBitSum += BitSum (*(theBytePtr++)); theBitSum += BitSum (*(theBytePtr++)); theBitSum += BitSum (*(theBytePtr++)); ++gCountArray [theBitSum];} while (--theIndex); theExpected = theAdjLength / 4; theExpected *= 1732525830.0 / 4294967296.0; theDiff = gCountArray [15] + gCountArray [16] + gCountArray [17]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat = theDiff; theDiff = theAdjLength / 4; theExpected = (theDiff - theExpected) / 2.0; theDiff = gCountArray [0] + gCountArray [1] + gCountArray [2] + gCountArray [3] + gCountArray [4] + gCountArray [5] + gCountArray [6] + gCountArray [7] + gCountArray [8] + gCountArray [9] + gCountArray [10] + gCountArray [11] + gCountArray [12] + gCountArray [13] + gCountArray [14]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theDiff = gCountArray [18] + gCountArray [19] + gCountArray [20] + gCountArray [21] + gCountArray [22] + gCountArray [23] + gCountArray [24] + gCountArray [25] + gCountArray [26] + gCountArray [27] + gCountArray [28] + gCountArray [29] + gCountArray [30] + gCountArray [31] + gCountArray [32]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; return (theStat);} /* The sequence is divided into 6-byte segments. For each segment, the bits in each byte are summed. The bit sums will be either even or odd. In a random sequence the probability of both even and odd bit sums is 1/2. In each segment, there are 64 possible permutations of even and odd bit sums, each with a probability of 1/64. The number of segments with each possible permutation is counted. A chi-squared statistic is calculated. The degrees-of-freedom is 63. */ double OddEvenBitSumStat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt theIndex, theAdjLength, theValue; const OneByteInt *theBytePtr; double theStat, theExpected, theDiff; *theDF = 63.0; if (!(theAdjLength = (theDataLength / 384) * 384)) return (0.0); theIndex = 64; do gCountArray [--theIndex] = 0; while (theIndex); theIndex = theAdjLength / 6; theBytePtr = theData; do {theValue = BitSum (*(theBytePtr++)) & 0x01; theValue = (theValue << 1) | (BitSum (*(theBytePtr++)) & 0x01); theValue = (theValue << 1) | (BitSum (*(theBytePtr++)) & 0x01); theValue = (theValue << 1) | (BitSum (*(theBytePtr++)) & 0x01); theValue = (theValue << 1) | (BitSum (*(theBytePtr++)) & 0x01); theValue = (theValue << 1) | (BitSum (*(theBytePtr++)) & 0x01); ++gCountArray [theValue];} while (--theIndex); theStat = 0.0; theExpected = theAdjLength / 384; theIndex = 64; do {theDiff = gCountArray [--theIndex]; theDiff -= theExpected; theDiff *= theDiff; theStat += theDiff;} while (theIndex); theStat /= theExpected; return (theStat);} /* The sequence is divided into 64-bit segments. Eight bytes are built from these 64 bits. The first byte is built from the 1st, 9th, 17th, 25th, etc. bits in the segment; the second is built from the 2nd, 10th, 18th, 26th, etc. bits, and so on. The total number of resulting bytes of each possible value is totaled over the whole sequence. In a random sequence the probability of each possible byte value is 1/256. A chi-squared statistic is calculated. The degrees-of-freedom is 255. */ double Matrix64Stat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt seqIndex, theAdjLength; FourByteInt indexA, indexB; const OneByteInt *theSeqPtr; OneByteInt seqByte, theBytes [8], *theBytePtr; double theStat, theExpected, theDiff; *theDF = 255.0; if (!(theAdjLength = theDataLength & 0xFFFFFF00)) return (0.0); seqIndex = 256; do gCountArray [--seqIndex] = 0; while (seqIndex); theSeqPtr = theData; seqIndex = theAdjLength; do {indexA = 8; do {seqByte = *(theSeqPtr++); --seqIndex; theBytePtr = theBytes; indexB = 8; do {*theBytePtr <<= 1; *(theBytePtr++) |= seqByte & 0x01; seqByte >>= 1;} while (--indexB);} while (--indexA); theBytePtr = theBytes; indexB = 8; do ++gCountArray [*(theBytePtr++)]; while (--indexB);} while (seqIndex); theStat = 0.0; theExpected = theAdjLength / 256; seqIndex = 256; do {theDiff = gCountArray [--seqIndex]; theDiff -= theExpected; theDiff *= theDiff; theStat += theDiff;} while (seqIndex); theStat /= theExpected; return (theStat);} /* The sequence is divided into 128-bit segments. Sixteen bytes are built from these 128 bits. The first byte is built from the 1st, 17th, 33th, 49th, etc. bits in the segment; the second is built from the 2nd, 18th, 34th, 50th, etc. bits, and so on. The total number of resulting bytes of each possible value is totaled over the whole sequence. In a random sequence the probability of each possible byte value is 1/256. A chi-squared statistic is calculated. The degrees-of-freedom is 255. */ double Matrix128Stat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt seqIndex, theAdjLength; FourByteInt indexA, indexB; const OneByteInt *theSeqPtr; OneByteInt seqByteA, seqByteB, theBytes [16], *theBytePtr; double theStat, theExpected, theDiff; *theDF = 255.0; if (!(theAdjLength = theDataLength & 0xFFFFFF00)) return (0.0); seqIndex = 256; do gCountArray [--seqIndex] = 0; while (seqIndex); theSeqPtr = theData; seqIndex = theAdjLength / 2; do {indexA = 8; do {seqByteA = *(theSeqPtr++); seqByteB = *(theSeqPtr++); --seqIndex; theBytePtr = theBytes; indexB = 8; do {*theBytePtr <<= 1; *(theBytePtr++) |= seqByteA & 0x01; seqByteA >>= 1;} while (--indexB); indexB = 8; do {*theBytePtr <<= 1; *(theBytePtr++) |= seqByteB & 0x01; seqByteB >>= 1;} while (--indexB);} while (--indexA); theBytePtr = theBytes; indexB = 16; do ++gCountArray [*(theBytePtr++)]; while (--indexB);} while (seqIndex); theStat = 0.0; theExpected = theAdjLength / 256; seqIndex = 256; do {theDiff = gCountArray [--seqIndex]; theDiff -= theExpected; theDiff *= theDiff; theStat += theDiff;} while (seqIndex); theStat /= theExpected; return (theStat);} /* The sequence is divided into 256-bit segments. Thirty-two bytes are built from these 256 bits. The first byte is built from the 1st, 33th, 65th, 97th, etc. bits in the segment; the second is built from the 2nd, 34th, 66th, 98th, etc. bits, and so on. The total number of resulting bytes of each possible value is totaled over the whole sequence. In a random sequence the probability of each possible byte value is 1/256. A chi-squared statistic is calculated. The degrees-of-freedom is 255. */ double Matrix256Stat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt seqIndex, theAdjLength, seqLong; FourByteInt indexA, indexB; const OneByteInt *theSeqPtr; OneByteInt theBytes [32], *theBytePtr; double theStat, theExpected, theDiff; *theDF = 255.0; if (!(theAdjLength = theDataLength & 0xFFFFFF00)) return (0.0); seqIndex = 256; do gCountArray [--seqIndex] = 0; while (seqIndex); theSeqPtr = theData; seqIndex = theAdjLength / 4; do {indexA = 8; do {seqLong = *((FourByteInt *) theSeqPtr) & 0xFFFFFFFF; theSeqPtr += 4; --seqIndex; theBytePtr = theBytes; indexB = 32; do {*theBytePtr <<= 1; *(theBytePtr++) |= seqLong & 0x01; seqLong >>= 1;} while (--indexB);} while (--indexA); theBytePtr = theBytes; indexB = 32; do ++gCountArray [*(theBytePtr++)]; while (--indexB);} while (seqIndex); theStat = 0.0; theExpected = theAdjLength / 256; seqIndex = 256; do {theDiff = gCountArray [--seqIndex]; theDiff -= theExpected; theDiff *= theDiff; theStat += theDiff;} while (seqIndex); theStat /= theExpected; return (theStat);} /* The sequence is divided into 512-bit segments. Sixty-four bytes are built from these 512 bits. The first byte is built from the 1st, 65th, 129th, 193th, etc. bits in the segment; the second is built from the 2nd, 66th, 130th, 194th, etc. bits, and so on. The total number of resulting bytes of each possible value is totaled over the whole sequence. In a random sequence the probability of each possible byte value is 1/256. A chi-squared statistic is calculated. The degrees-of-freedom is 255. */ double Matrix512Stat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt seqIndex, theAdjLength, seqLongA, seqLongB; FourByteInt indexA, indexB; const OneByteInt *theSeqPtr; OneByteInt theBytes [64], *theBytePtr; double theStat, theExpected, theDiff; *theDF = 255.0; if (!(theAdjLength = theDataLength & 0xFFFFFF00)) return (0.0); seqIndex = 256; do gCountArray [--seqIndex] = 0; while (seqIndex); theSeqPtr = theData; seqIndex = theAdjLength / 8; do {indexA = 8; do {seqLongA = *((FourByteInt *) theSeqPtr) & 0xFFFFFFFF; theSeqPtr += 4; seqLongB = *((FourByteInt *) theSeqPtr) & 0xFFFFFFFF; theSeqPtr += 4; --seqIndex; theBytePtr = theBytes; indexB = 32; do {*theBytePtr <<= 1; *(theBytePtr++) |= seqLongA & 0x01; seqLongA >>= 1;} while (--indexB); indexB = 32; do {*theBytePtr <<= 1; *(theBytePtr++) |= seqLongB & 0x01; seqLongB >>= 1;} while (--indexB);} while (--indexA); theBytePtr = theBytes; indexB = 64; do ++gCountArray [*(theBytePtr++)]; while (--indexB);} while (seqIndex); theStat = 0.0; theExpected = theAdjLength / 256; seqIndex = 256; do {theDiff = gCountArray [--seqIndex]; theDiff -= theExpected; theDiff *= theDiff; theStat += theDiff;} while (seqIndex); theStat /= theExpected; return (theStat);} /* The sequence is divided into 1024-bit segments. One hundred and twenty-eight bytes are built from these 1024 bits. The first byte is built from the 1st, 129th, 257th, 385th, etc. bits in the segment; the second is built from the 2nd, 130th, 258th, 386th, etc. bits, and so on. The total number of resulting bytes of each possible value is totaled over the whole sequence. In a random sequence the probability of each possible byte value is 1/256. A chi-squared statistic is calculated. The degrees-of-freedom is 255. */ double Matrix1KStat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt seqIndex, theAdjLength; FourByteInt indexA, indexB; FourByteInt seqLongA, seqLongB, seqLongC, seqLongD; const OneByteInt *theSeqPtr; OneByteInt theBytes [128], *theBytePtr; double theStat, theExpected, theDiff; *theDF = 255.0; if (!(theAdjLength = theDataLength & 0xFFFFFF00)) return (0.0); seqIndex = 256; do gCountArray [--seqIndex] = 0; while (seqIndex); theSeqPtr = theData; seqIndex = theAdjLength / 16; do {indexA = 8; do {seqLongA = *((FourByteInt *) theSeqPtr) & 0xFFFFFFFF; theSeqPtr += 4; seqLongB = *((FourByteInt *) theSeqPtr) & 0xFFFFFFFF; theSeqPtr += 4; seqLongC = *((FourByteInt *) theSeqPtr) & 0xFFFFFFFF; theSeqPtr += 4; seqLongD = *((FourByteInt *) theSeqPtr) & 0xFFFFFFFF; theSeqPtr += 4; --seqIndex; theBytePtr = theBytes; indexB = 32; do {*theBytePtr <<= 1; *(theBytePtr++) |= seqLongA & 0x01; seqLongA >>= 1;} while (--indexB); indexB = 32; do {*theBytePtr <<= 1; *(theBytePtr++) |= seqLongB & 0x01; seqLongB >>= 1;} while (--indexB); indexB = 32; do {*theBytePtr <<= 1; *(theBytePtr++) |= seqLongC & 0x01; seqLongC >>= 1;} while (--indexB); indexB = 32; do {*theBytePtr <<= 1; *(theBytePtr++) |= seqLongD & 0x01; seqLongD >>= 1;} while (--indexB);} while (--indexA); theBytePtr = theBytes; indexB = 128; do ++gCountArray [*(theBytePtr++)]; while (--indexB);} while (seqIndex); theStat = 0.0; theExpected = theAdjLength / 256; seqIndex = 256; do {theDiff = gCountArray [--seqIndex]; theDiff -= theExpected; theDiff *= theDiff; theStat += theDiff;} while (seqIndex); theStat /= theExpected; return (theStat);} /* The sequence is divided into 2048-bit segments. Two hundred and fifty-six bytes are built from these 2048 bits. The first byte is built from the 1st, 257th, 513th, 769th, etc. bits in the segment; the second is built from the 2nd, 258th, 514th, 770th, etc. bits, and so on. The total number of resulting bytes of each possible value is totaled over the whole sequence. In a random sequence the probability of each possible byte value is 1/256. A chi-squared statistic is calculated. The degrees-of-freedom is 255. */ double Matrix2KStat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt seqIndex, theAdjLength; FourByteInt indexA, indexB; FourByteInt seqLongA, seqLongB, seqLongC, seqLongD; FourByteInt seqLongE, seqLongF, seqLongG, seqLongH; const OneByteInt *theSeqPtr; OneByteInt theBytes [256], *theBytePtr; double theStat, theExpected, theDiff; *theDF = 255.0; if (!(theAdjLength = theDataLength & 0xFFFFFF00)) return (0.0); seqIndex = 256; do gCountArray [--seqIndex] = 0; while (seqIndex); theSeqPtr = theData; seqIndex = theAdjLength / 32; do {indexA = 8; do {seqLongA = *((FourByteInt *) theSeqPtr) & 0xFFFFFFFF; theSeqPtr += 4; seqLongB = *((FourByteInt *) theSeqPtr) & 0xFFFFFFFF; theSeqPtr += 4; seqLongC = *((FourByteInt *) theSeqPtr) & 0xFFFFFFFF; theSeqPtr += 4; seqLongD = *((FourByteInt *) theSeqPtr) & 0xFFFFFFFF; theSeqPtr += 4; seqLongE = *((FourByteInt *) theSeqPtr) & 0xFFFFFFFF; theSeqPtr += 4; seqLongF = *((FourByteInt *) theSeqPtr) & 0xFFFFFFFF; theSeqPtr += 4; seqLongG = *((FourByteInt *) theSeqPtr) & 0xFFFFFFFF; theSeqPtr += 4; seqLongH = *((FourByteInt *) theSeqPtr) & 0xFFFFFFFF; theSeqPtr += 4; --seqIndex; theBytePtr = theBytes; indexB = 32; do {*theBytePtr <<= 1; *(theBytePtr++) |= seqLongA & 0x01; seqLongA >>= 1;} while (--indexB); indexB = 32; do {*theBytePtr <<= 1; *(theBytePtr++) |= seqLongB & 0x01; seqLongB >>= 1;} while (--indexB); indexB = 32; do {*theBytePtr <<= 1; *(theBytePtr++) |= seqLongC & 0x01; seqLongC >>= 1;} while (--indexB); indexB = 32; do {*theBytePtr <<= 1; *(theBytePtr++) |= seqLongD & 0x01; seqLongD >>= 1;} while (--indexB); indexB = 32; do {*theBytePtr <<= 1; *(theBytePtr++) |= seqLongE & 0x01; seqLongE >>= 1;} while (--indexB); indexB = 32; do {*theBytePtr <<= 1; *(theBytePtr++) |= seqLongF & 0x01; seqLongF >>= 1;} while (--indexB); indexB = 32; do {*theBytePtr <<= 1; *(theBytePtr++) |= seqLongG & 0x01; seqLongG >>= 1;} while (--indexB); indexB = 32; do {*theBytePtr <<= 1; *(theBytePtr++) |= seqLongH & 0x01; seqLongH >>= 1;} while (--indexB);} while (--indexA); theBytePtr = theBytes; indexB = 256; do ++gCountArray [*(theBytePtr++)]; while (--indexB);} while (seqIndex); theStat = 0.0; theExpected = theAdjLength / 256; seqIndex = 256; do {theDiff = gCountArray [--seqIndex]; theDiff -= theExpected; theDiff *= theDiff; theStat += theDiff;} while (seqIndex); theStat /= theExpected; return (theStat);} /* An algorithm is used to predict the value of each bit of the sequence from the beginning of the sequence to the end. In a random sequence the probability of success of any such algorithm is 1/2. The number of successes is counted. A chi-squared statistic is calculated. The degrees-of-freedom is 1. The algorithm is as follows: the numbers of zeros and ones in all the previous bits are counted. If the ones outnumber the zeros, a zero is predicted; if the zeros outnumber the ones, a one is predicted. Otherwise the prediction is the same as for the previous bit. The interesting part of this test is the prediction made when the numbers of zeros and ones are equal. */ double BitPredictAStat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt theIndex, theAdjLength; FourByteInt theTotalCorr, totalZeros, totalOnes; OneByteInt theByte, theBit, predictedBit, bitIndex; const OneByteInt *theBytePtr; double theStat, theExpected; *theDF = 1.0; if (!(theAdjLength = theDataLength)) return (0.0); if (theAdjLength > 0x1FFFFFFF) theAdjLength = 0x1FFFFFFF; theTotalCorr = totalZeros = totalOnes = 0; predictedBit = 0x00; theIndex = theAdjLength; theBytePtr = theData; do {theByte = *(theBytePtr++); bitIndex = 8; do {theBit = (theByte & 0x01); if (theBit == predictedBit) ++theTotalCorr; if (theBit) ++totalOnes; else ++totalZeros; if (totalZeros > totalOnes) predictedBit = 0x01; else if (totalOnes > totalZeros) predictedBit = 0x00; theByte >>= 1;} while (--bitIndex);} while (--theIndex); theExpected = theAdjLength * 4; theStat = theTotalCorr; theStat -= theExpected; theStat *= theStat; theStat /= theExpected; theStat *= 2.0; return (theStat);} /* An algorithm is used to predict the value of each bit of the sequence from the beginning of the sequence to the end. In a random sequence the probability of success of any such algorithm is 1/2. The number of successes is counted. A chi-squared statistic is calculated. The degrees-of-freedom is 1. The algorithm is as follows: the numbers of zeros and ones in the previous (theLength * 2) + 1 bits are counted. If the ones outnumber the zeros, a zero is predicted; if the zeros outnumber the ones, a one is predicted. */ double BitPredictStat (const void *theData, FourByteInt theDataLength, double *theDF, FourByteInt theLength) {FourByteInt theIndex, theAdjLength, theTotalCorr; FourByteInt arrayIndex, wLen, hLen, bitSum; OneByteInt theByte, theBit, predictedBit, bitIndex; OneByteInt *arrayPtr; const OneByteInt *theBytePtr; double theStat, theExpected; *theDF = 1.0; theAdjLength = theDataLength; hLen = theLength; if (!theAdjLength || hLen > (sizeof (gByteArray) - 1) / 2) return (0.0); if (theAdjLength > 0x1FFFFFFF) theAdjLength = 0x1FFFFFFF; theTotalCorr = 0; wLen = (hLen * 2) + 1; arrayIndex = wLen; do gByteArray [--arrayIndex] = 0x00; while (arrayIndex); arrayPtr = gByteArray; arrayIndex = wLen; predictedBit = 0x01; bitSum = 0; theIndex = theAdjLength; theBytePtr = theData; do {theByte = *(theBytePtr++); bitIndex = 8; do {theBit = (theByte & 0x01); if (theBit == predictedBit) ++theTotalCorr; bitSum -= *arrayPtr; *arrayPtr = theBit; bitSum += theBit; if (--arrayIndex) ++arrayPtr; else {arrayPtr = gByteArray; arrayIndex = wLen;} if (bitSum > hLen) predictedBit = 0x00; else predictedBit = 0x01; theByte >>= 1;} while (--bitIndex);} while (--theIndex); theExpected = theAdjLength * 4; theStat = theTotalCorr; theStat -= theExpected; theStat *= theStat; theStat /= theExpected; theStat *= 2.0; return (theStat);} /* An algorithm is used to predict the value of each byte of the sequence from the beginning of the sequence to the end. In a random sequence the probability of success of any such algorithm is 1/256. The number of successes is counted. A chi-squared statistic is calculated. The degrees-of-freedom is 1. The algorithm is as follows: the next byte is predicted to be equal to all the previous bytes bitwise XORed together. The first byte of the sequence is predicted to equal zero. */ double BytePredictAStat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt theIndex, theAdjLength, theTotalCorr; OneByteInt theByte, expByte; const OneByteInt *theBytePtr; double theStat, theExpected, theDiff; *theDF = 1.0; if (!(theAdjLength = theDataLength & 0xFFFFFF00)) return (0.0); theTotalCorr = 0; expByte = 0x00; theIndex = theAdjLength; theBytePtr = theData; do {theByte = *(theBytePtr++); if (theByte == expByte) ++theTotalCorr; expByte ^= theByte;} while (--theIndex); theExpected = theAdjLength / 256; theDiff = theTotalCorr; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat = theDiff; theExpected *= 255.0; theDiff = theAdjLength - theTotalCorr; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; return (theStat);} /* An algorithm is used to predict the value of each byte of the sequence from the beginning of the sequence to the end. In a random sequence the probability of success of any such algorithm is 1/256. The number of successes is counted. A chi-squared statistic is calculated. The degrees-of-freedom is 1. The algorithm is as follows: the next byte is predicted to be equal to the sum of all the previous bytes, modulo 256. The first byte of the sequence is predicted to equal zero. */ double BytePredictBStat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt theIndex, theAdjLength, theTotalCorr; OneByteInt theByte, expByte; const OneByteInt *theBytePtr; double theStat, theExpected, theDiff; *theDF = 1.0; if (!(theAdjLength = theDataLength & 0xFFFFFF00)) return (0.0); theTotalCorr = 0; expByte = 0x00; theIndex = theAdjLength; theBytePtr = theData; do {theByte = *(theBytePtr++); if (theByte == expByte) ++theTotalCorr; expByte += theByte;} while (--theIndex); theExpected = theAdjLength / 256; theDiff = theTotalCorr; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat = theDiff; theExpected *= 255.0; theDiff = theAdjLength - theTotalCorr; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; return (theStat);} /* An algorithm is used to predict the value of each byte of the sequence from the beginning of the sequence to the end. In a random sequence the probability of success of any such algorithm is 1/256. The number of successes is counted. A chi-squared statistic is calculated. The degrees-of-freedom is 1. The algorithm is as follows: the next byte value is predicted to be zero until the first zero is found. From that point on, the next byte value is predicted to be the byte value whose last appearance was furthest back in the sequence. */ double BytePredictCStat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt theIndex, theAdjLength, theTotal, theMax, arrayIndex; OneByteInt theByte, expByte; const OneByteInt *theBytePtr; double theStat, theExpected, theDiff; *theDF = 1.0; if (!(theAdjLength = theDataLength & 0xFFFFFF00)) return (0.0); arrayIndex = 256; do gCountArray [--arrayIndex] = 0; while (arrayIndex); theTotal = 0; expByte = 0x00; theIndex = theAdjLength; theBytePtr = theData; do {theByte = *(theBytePtr++); gCountArray [theByte] = --theIndex; if (theByte == expByte) {++theTotal; theMax = 0; arrayIndex = 256; do {if (theMax < gCountArray [--arrayIndex]) {theMax = gCountArray [arrayIndex]; expByte = arrayIndex;}} while (arrayIndex);}} while (theIndex); theExpected = theAdjLength / 256; theDiff = theTotal; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat = theDiff; theExpected *= 255.0; theDiff = theAdjLength - theTotal; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; return (theStat);} /* An algorithm is used to predict the value of each byte of the sequence from the beginning of the sequence to the end. In a random sequence the probability of success of any such algorithm is 1/256. The number of successes is counted. A chi-squared statistic is calculated. The degrees-of-freedom is 1. The algorithm is as follows: a given byte value is predicted to be followed by the same byte value it was followed by the last time it appeared in the sequence. A byte value that has not previously appeared in the sequence is predicted to be followed by the byte value of the first byte of the sequence. The first byte of the sequence is predicted to equal zero. */ double BytePredictDStat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt theIndex, theAdjLength, theTotal; OneByteInt theByte, expByte; const OneByteInt *theBytePtr; double theStat, theExpected, theDiff; *theDF = 1.0; if (!(theAdjLength = theDataLength & 0xFFFFFF00)) return (0.0); theIndex = 256; do gCountArray [--theIndex] = 0; while (theIndex); theTotal = 0; expByte = 0x00; theIndex = 0; theBytePtr = theData; do {theByte = theBytePtr [theIndex]; if (theByte == expByte) ++theTotal; expByte = theBytePtr [gCountArray [theByte]]; gCountArray [theByte] = ++theIndex;} while (theIndex < theAdjLength); theExpected = theAdjLength / 256; theDiff = theTotal; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat = theDiff; theExpected *= 255.0; theDiff = theAdjLength - theTotal; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; return (theStat);} /* Each byte in the sequence either is or is not identical to its preceding byte. In a random sequence the probability that a byte will be identical to its preceding byte is 1/256. The total number of bytes identical to their preceding bytes is counted. A chi-squared statistic is calculated. The degrees-of-freedom is 1. This test is equivalent to a byte prediction test where each byte is predicted to be equal to its preceding byte. The first byte of the sequence is predicted to equal the last byte of the sequence. */ double ByteRepStat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt theIndex, theAdjLength, theTotal; OneByteInt thisByte, lastByte; const OneByteInt *theBytePtr; double theStat, theExpected, theDiff; *theDF = 1.0; if (!(theAdjLength = theDataLength & 0xFFFFFF00)) return (0.0); theTotal = 0; theIndex = theAdjLength; theBytePtr = theData; lastByte = theBytePtr [theIndex - 1]; do {thisByte = *(theBytePtr++); if (thisByte == lastByte) ++theTotal; lastByte = thisByte;} while (--theIndex); theExpected = theAdjLength / 256; theDiff = theTotal; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat = theDiff; theExpected *= 255.0; theDiff = theAdjLength - theTotal; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; return (theStat);} /* The sequence is divided into two-byte segments. The two bytes in each segment bitwise ANDed together. The result (d) will fall into one of 7 ranges: [d<37], [37<=d<74], [74<=d<111], [111<=d<148], [148<=d<185], [185<=d<222], or [d>=222]. In a random sequence the probabilities of these ranges are 32265/65536, 10755/65536, 5355/65536, 8985/65536, 3969/65536, 3171/65536, and 1036/65536, respectively. The total number of segments in each of these 7 ranges is counted. A chi-squared statistic is calculated. The degrees-of-freedom is 6. */ double TwoByteANDStat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt theIndex, theAdjLength; OneByteInt theByte; const OneByteInt *theBytePtr; double theStat, theExpected, theDiff; *theDF = 6.0; if (!(theAdjLength = theDataLength & 0xFFFE0000)) return (0.0); theIndex = 7; do gCountArray [--theIndex] = 0; while (theIndex); theIndex = theAdjLength / 2; theBytePtr = theData; do {theByte = *(theBytePtr++); theByte &= *(theBytePtr++); ++gCountArray [theByte / 37];} while (--theIndex); theExpected = 32265 * (theAdjLength / 131072); theDiff = gCountArray [0]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat = theDiff; theExpected = 10755 * (theAdjLength / 131072); theDiff = gCountArray [1]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theExpected = 5355 * (theAdjLength / 131072); theDiff = gCountArray [2]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theExpected = 8985 * (theAdjLength / 131072); theDiff = gCountArray [3]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theExpected = 3969 * (theAdjLength / 131072); theDiff = gCountArray [4]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theExpected = 3171 * (theAdjLength / 131072); theDiff = gCountArray [5]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theExpected = 1036 * (theAdjLength / 131072); theDiff = gCountArray [6]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; return (theStat);} /* The sequence is divided into four-byte segments. The four bytes in each segment bitwise ANDed together. The result (d) will fall into one of 3 ranges: [d<73], [73<=d<146], or [d>=146]. In a random sequence the probabilities of these ranges are 3993624225/4294967296, 266241615/4294967296, and 35101456/4294967296, respectively. The total number of segments in each of these 3 ranges is counted. A chi-squared statistic is calculated. The degrees-of-freedom is 2. */ double FourByteANDStat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt theIndex, theAdjLength; OneByteInt theByte; const OneByteInt *theBytePtr; double theStat, theExpected, theDiff; *theDF = 2.0; if (!(theAdjLength = theDataLength & 0xFFFFFFFC)) return (0.0); theIndex = 4; do gCountArray [--theIndex] = 0; while (theIndex); theIndex = theAdjLength / 4; theBytePtr = theData; do {theByte = *(theBytePtr++); theByte &= *(theBytePtr++); theByte &= *(theBytePtr++); theByte &= *(theBytePtr++); ++gCountArray [theByte / 73];} while (--theIndex); theExpected = theAdjLength / 4; theExpected *= 3993624225.0 / 4294967296.0; theDiff = gCountArray [0]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat = theDiff; theExpected = theAdjLength / 4; theExpected *= 266241615.0 / 4294967296.0; theDiff = gCountArray [1]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theExpected = theAdjLength / 4; theExpected *= 35101456.0 / 4294967296.0; theDiff = gCountArray [2] + gCountArray [3]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; return (theStat);} /* The sequence is divided into four-byte segments. The first two bytes and the last two bytes in each segment are bitwise ANDed together. The results of these two AND operations are then bitwise ORed together. The result (d) will fall into one of 3 ranges: [d<73], [73<=d<146], or [d>=146]. In a random sequence the probabilities of these ranges are 1573112097/4294967296, 1223531631/4294967296, and 1498323568/4294967296, respectively. The total number of segments in each of these 3 ranges is counted. A chi-squared statistic is calculated. The degrees-of-freedom is 2. */ double FourByteAND_ORStat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt theIndex, theAdjLength; OneByteInt firstByte, secondByte; const OneByteInt *theBytePtr; double theStat, theExpected, theDiff; *theDF = 2.0; if (!(theAdjLength = theDataLength & 0xFFFFFFFC)) return (0.0); theIndex = 4; do gCountArray [--theIndex] = 0; while (theIndex); theIndex = theAdjLength / 4; theBytePtr = theData; do {firstByte = *(theBytePtr++); firstByte &= *(theBytePtr++); secondByte = *(theBytePtr++); secondByte &= *(theBytePtr++); ++gCountArray [(firstByte | secondByte) / 73];} while (--theIndex); theExpected = theAdjLength / 4; theExpected *= 1573112097.0 / 4294967296.0; theDiff = gCountArray [0]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat = theDiff; theExpected = theAdjLength / 4; theExpected *= 1223531631.0 / 4294967296.0; theDiff = gCountArray [1]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theExpected = theAdjLength / 4; theExpected *= 1498323568.0 / 4294967296.0; theDiff = gCountArray [2] + gCountArray [3]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; return (theStat);} /* The sequence is divided into two-byte segments. The second byte in each segment is either equal to, greater than, or less than the first. The total number of segments in each of these three categories is counted. In a random sequence the probability that the second byte will be equal to the first is 1/256. The probabilities that second byte will be greater than or less than the first are both 255/512. A chi-squared statistic is calculated. The degrees-of-freedom is 2. */ double TwoByteUpDownStat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt theIndex, theAdjLength; FourByteInt totalSame, totalUp, totalDown; OneByteInt firstByte, secondByte; const OneByteInt *theBytePtr; double theStat, theExpected, theDiff; *theDF = 2.0; if (!(theAdjLength = theDataLength & 0xFFFFFE00)) return (0.0); totalSame = totalUp = totalDown = 0; theIndex = theAdjLength / 2; theBytePtr = theData; do {firstByte = *(theBytePtr++); secondByte = *(theBytePtr++); if (firstByte > secondByte) ++totalDown; else if (firstByte < secondByte) ++totalUp; else ++totalSame;} while (--theIndex); theExpected = theAdjLength / 512; theDiff = totalSame; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat = theDiff; theExpected *= 127.5; theDiff = totalUp; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theDiff = totalDown; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; return (theStat);} /* The sequence is divided into four-byte segments. If the bytes of each segment are called a, b, c, and d; then comparing these four bytes will classify the segments into 8 categories: [a<=b<=c<=d], [a>b<=c<=d], [a<=b>c<=d], [a>b>c<=d], [a<=b<=c>d], [a>b<=c>d], [a<=b>c>d], and [a>b>c>d]. In a random sequence, the probabilities that a segment will belong to each of these categories are 2862209/67108864, 8454015/67108864, 14046335/67108864, 8322945/67108864, 8454015/67108864, 13915265/67108864, 8322945/67108864, and 2731135/67108864, respectively. The total number of segments in each category is counted. A chi-squared statistic is calculated. The degrees-of-freedom is 7. */ double FourByteUpDwnStat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt theIndex, theAdjLength; OneByteInt firstByte, secondByte, thirdByte, fourthByte; const OneByteInt *theBytePtr; double theStat, theExpected, theDiff; *theDF = 7.0; if (!(theAdjLength = theDataLength & 0xFFFFFFFC)) return (0.0); theIndex = 8; do gCountArray [--theIndex] = 0; while (theIndex); theIndex = theAdjLength / 4; theBytePtr = theData; do {firstByte = *(theBytePtr++); secondByte = *(theBytePtr++); thirdByte = *(theBytePtr++); fourthByte = *(theBytePtr++); ++gCountArray [(firstByte <= secondByte ? 0 : 1) + (secondByte <= thirdByte ? 0 : 2) + (thirdByte <= fourthByte ? 0 : 4)];} while (--theIndex); theExpected = theAdjLength / 4; theExpected *= 2862209.0 / 67108864.0; theDiff = gCountArray [0]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat = theDiff; theExpected = theAdjLength / 4; theExpected *= 8454015.0 / 67108864.0; theDiff = gCountArray [1]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theExpected = theAdjLength / 4; theExpected *= 14046335.0 / 67108864.0; theDiff = gCountArray [2]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theExpected = theAdjLength / 4; theExpected *= 8322945.0 / 67108864.0; theDiff = gCountArray [3]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theExpected = theAdjLength / 4; theExpected *= 8454015.0 / 67108864.0; theDiff = gCountArray [4]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theExpected = theAdjLength / 4; theExpected *= 13915265.0 / 67108864.0; theDiff = gCountArray [5]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theExpected = theAdjLength / 4; theExpected *= 8322945.0 / 67108864.0; theDiff = gCountArray [6]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theExpected = theAdjLength / 4; theExpected *= 2731135.0 / 67108864.0; theDiff = gCountArray [7]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; return (theStat);} /* The sequence is divided into four-byte segments. For each segment, the absolute difference between the first and third bytes of the segment is added to the absolute difference between the second and fourth bytes of the segment. The result (d) will fall into one of 7 ranges: [d<64], [64<=d<128], [128<=d<192], [192<=d<256], [256<=d<320], [320<=d<384], or [d>=384]. In a random sequence the probabilities of these ranges are 6935371/67108864, 15991947/67108864, 18225611/67108864, 14683915/67108864, 7696309/67108864, 2865781/67108864, and 709930/67108864, respectively. The total number of segments in each of these 7 ranges is counted. A chi-squared statistic is calculated. The degrees-of-freedom is 6. */ double RectDistStat (const void *theData, FourByteInt theDataLength, double *theDF) {FourByteInt theIndex, theAdjLength; OneByteInt firstByte, secondByte, thirdByte, fourthByte; const OneByteInt *theBytePtr; double theStat, theExpected, theDiff; *theDF = 6.0; if (!(theAdjLength = theDataLength & 0xFFFFFFFC)) return (0.0); theIndex = 8; do gCountArray [--theIndex] = 0; while (theIndex); theIndex = theAdjLength / 4; theBytePtr = theData; do {firstByte = *(theBytePtr++); secondByte = *(theBytePtr++); thirdByte = *(theBytePtr++); fourthByte = *(theBytePtr++); ++gCountArray [((int) (firstByte > thirdByte ? firstByte - thirdByte : thirdByte - firstByte) + (secondByte > fourthByte ? secondByte - fourthByte : fourthByte - secondByte)) / 64];} while (--theIndex); theExpected = theAdjLength / 4; theExpected *= 6935371.0 / 67108864.0; theDiff = gCountArray [0]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat = theDiff; theExpected = theAdjLength / 4; theExpected *= 15991947.0 / 67108864.0; theDiff = gCountArray [1]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theExpected = theAdjLength / 4; theExpected *= 18225611.0 / 67108864.0; theDiff = gCountArray [2]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theExpected = theAdjLength / 4; theExpected *= 14683915.0 / 67108864.0; theDiff = gCountArray [3]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theExpected = theAdjLength / 4; theExpected *= 7696309.0 / 67108864.0; theDiff = gCountArray [4]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theExpected = theAdjLength / 4; theExpected *= 2865781.0 / 67108864.0; theDiff = gCountArray [5]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; theExpected = theAdjLength / 4; theExpected *= 709930.0 / 67108864.0; theDiff = gCountArray [6] + gCountArray [7]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; return (theStat);} /* The sequence is divided into segments. The length of the segments, in bytes, is segLength. If any byte value occurs more than once in a segment, a collision is said to have occurred. In a random sequence, the probability that no collisions will occur in a segment is (255/256) * (254/256) * (253/256) * ... * ((257 - segLength)/256). The total number of segments with and without collisions is counted. A chi-squared statistic is calculated. The degrees-of-freedom is 1. */ double ByteCollisionStat (const void *theData, FourByteInt theDataLength, double *theDF, FourByteInt segLength) {FourByteInt theAdjLength, theSegments, segIndex, theTotal; const OneByteInt *theBytePtr; double theStat, theExpected, theDiff; *theDF = 1.0; if (!segLength || segLength >= 256 || segLength > theDataLength) return (0.0); theAdjLength = (theDataLength / segLength) * segLength; theBytePtr = theData; theSegments = theAdjLength / segLength; theTotal = 0; do {segIndex = 256; do gCountArray [--segIndex] = 0; while (segIndex); segIndex = segLength; do {if (++gCountArray [*(theBytePtr++)] > 1) {++theTotal; theBytePtr += --segIndex; break;}} while (--segIndex);} while (--theSegments); theDiff = 256.0; theExpected = 1.0; segIndex = segLength; do {theExpected *= theDiff / 256.0; --theDiff;} while (--segIndex); theDiff = theAdjLength / segLength; theExpected = (1.0 - theExpected) * theDiff; theDiff = theTotal; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat = theDiff; theDiff = theAdjLength / segLength; theExpected = theDiff - theExpected; theDiff -= theTotal; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff; return (theStat);} /* The sequence is divided into pairs of bytes. The offset in the sequence between the bytes in each pair is theOffset. The two bytes in each pair are bitwise XORed together. The total number of resulting bytes of each possible value is counted. In a random sequence the probability of each possible byte value is 1/256. A chi-squared statistic is calculated. The degrees-of-freedom is 255. */ double OffsetByteXORStat (const void *theData, FourByteInt theDataLength, double *theDF, FourByteInt theOffset) {FourByteInt theIndex, gcdIndex, theAdjLength; OneByteInt theByte; const OneByteInt *theBytePtr; double theStat, theExpected, theDiff; *theDF = 255.0; theAdjLength = theDataLength & 0xFFFFFE00; if (!theOffset || theOffset >= theAdjLength || theOffset > 0xFFFFFFFF - theAdjLength) return (0.0); theIndex = 256; do gCountArray [--theIndex] = 0; while (theIndex); theBytePtr = theData; gcdIndex = GCD (theOffset, theAdjLength); do {theIndex = --gcdIndex; do {theByte = theBytePtr [theIndex]; theIndex = (theIndex + theOffset) % theAdjLength; theByte ^= theBytePtr [theIndex]; theIndex = (theIndex + theOffset) % theAdjLength; ++gCountArray [theByte];} while (theIndex != gcdIndex);} while (gcdIndex); theStat = 0.0; theExpected = theAdjLength / 512; theIndex = 256; do {theDiff = gCountArray [--theIndex]; theDiff -= theExpected; theDiff *= theDiff; theStat += theDiff;} while (theIndex); theStat /= theExpected; return (theStat);} /* For each byte of the sequence (x), (x mod theDivisor) is calculated. In a random sequence, the probability of a given remainder (r) is (q+1)/256, if r < (256 mod theDivisor), and q/256 otherwise (where q = 256/theDivisor, rounded down). The total number of byte values that are congruent to each possible remainder (mod theDivisor) is counted. A chi-squared statistic is calculated. The degrees-of-freedom is theDivisor - 1. */ double ModNStat (const void *theData, FourByteInt theDataLength, double *theDF, OneByteInt theDivisor) {FourByteInt theIndex, theAdjLength, theQuot, theRem; const OneByteInt *theBytePtr; double theStat, theExpected, theDiff; theAdjLength = theDataLength & 0xFFFFFF00; if (!theAdjLength || theDivisor < 2) return (*theDF = 1.0, 0.0); *theDF = theDivisor - 1; theIndex = theDivisor; do gCountArray [--theIndex] = 0; while (theIndex); theIndex = theAdjLength; theBytePtr = theData; do ++gCountArray [*(theBytePtr++) % theDivisor]; while (--theIndex); theStat = 0.0; theQuot = 256 / theDivisor; theRem = 256 % theDivisor; theIndex = 0; do {theExpected = (theAdjLength / 256) * (theQuot + (theRem > theIndex ? 1 : 0)); theDiff = gCountArray [theIndex]; theDiff -= theExpected; theDiff *= theDiff; theDiff /= theExpected; theStat += theDiff;} while (++theIndex < theDivisor); return (theStat);} FourByteInt BitSum (OneByteInt theByte) {FourByteInt bitCount, byteValue, theNumberOne; byteValue = theByte; theNumberOne = 1; bitCount = byteValue & theNumberOne; byteValue >>= 1; bitCount += byteValue & theNumberOne; byteValue >>= 1; bitCount += byteValue & theNumberOne; byteValue >>= 1; bitCount += byteValue & theNumberOne; byteValue >>= 1; bitCount += byteValue & theNumberOne; byteValue >>= 1; bitCount += byteValue & theNumberOne; byteValue >>= 1; bitCount += byteValue & theNumberOne; byteValue >>= 1; bitCount += byteValue & theNumberOne; return (bitCount);} FourByteInt GCD (FourByteInt numA, FourByteInt numB) {FourByteInt theA, theB, theGCD; theA = numA; theGCD = numB; do {theB = theGCD; theGCD = theA;} while ((theA = theB % theA)); return (theGCD);} void KeyGenerator (void *theKey, size_t theKeyLength) {OneByteInt *theBytePtr; unsigned long lcgSeed; size_t theIndex; theBytePtr = theKey; theIndex = theKeyLength; lcgSeed = gLCGSeed; do {lcgSeed = (69069UL * lcgSeed) + 1234567UL; *(theBytePtr++) = lcgSeed >> 24;} while (--theIndex); gLCGSeed = lcgSeed; return;} void ZeroDataBytes (void *theData, size_t theDataLength) {OneByteInt *theByte; size_t byteIndex, longIndex, *theLong; theLong = theData; byteIndex = theDataLength; if ((longIndex = byteIndex / sizeof (size_t))) {byteIndex -= longIndex * sizeof (size_t); do *(theLong++) = 0; while (--longIndex);} if (byteIndex) {theByte = (OneByteInt *) theLong; do *(theByte++) = 0; while (--byteIndex);} return;} /* Here begins a small library of math functions. These functions don't call any functions in , so we have to make sure the following macros are appropriately defined. HUGE_VAL is part of the C89 standard. HUGE_VALL and NAN are not part of the C89 standard, but they are in C99. These macros, and other macro definitions below, depend on definitions in . The functions also require and . */ #ifndef HUGE_VAL #define HUGE_VAL (DBL_MAX / DBL_MIN) #endif #ifndef HUGE_VALL #define HUGE_VALL (LDBL_MAX / LDBL_MIN) #endif #ifndef NAN #define NAN (0.0L / 0.0L) #endif /* HexTimesX multiplies the hexadecimal floating-point number specified in theString by theMultiplier. This is a rather computationally expensive way to insure maximum accuracy when multiplying by floating- point literals. It's probably only worth the trouble for literals that can't be precisely represented by the long double type. */ long double HexTimesX (char *theString, long double theMultiplier) {char *strPtr, theChar; long double theProduct, placeValue; size_t numOfDigits, theDigit; strPtr = theString; placeValue = theMultiplier; if (*strPtr == '-') {++strPtr; placeValue = -placeValue;} if (placeValue == HUGE_VALL || placeValue == -HUGE_VALL) return (errno = ERANGE, placeValue); theProduct = 0.0L; numOfDigits = theDigit = 0; while ((theChar = *(strPtr++))) {if (theChar == '.') break; if ((theChar >= 'a' && theChar <= 'f') || (theChar >= 'A' && theChar <= 'F') || (theChar >= '0' && theChar <= '9')) ++numOfDigits; else break;} if (theChar == '.') {while ((theChar = *(strPtr++))) {if ((theChar >= 'a' && theChar <= 'f') || (theChar >= 'A' && theChar <= 'F') || (theChar >= '0' && theChar <= '9')) {++numOfDigits; placeValue /= 16.0L;} else break;}} --strPtr; while (theDigit < numOfDigits) {theChar = *(--strPtr); if (theChar == '.') continue; if (theChar >= 'a' && theChar <= 'f') theChar -= 'a' - 10; else if (theChar >= 'A' && theChar <= 'F') theChar -= 'A' - 10; else theChar -= '0'; theProduct += placeValue * ((long double) theChar); placeValue *= 16.0L; ++theDigit;} if (theProduct == HUGE_VALL || theProduct == -HUGE_VALL) errno = ERANGE; return (theProduct);} /* LDFract returns theX minus the largest integer not greater than theX. E.g., LDFract (2.1L) returns 0.1L, and LDFract (-2.1L) returns 0.9L. */ long double LDFract (long double theX) {long double ldX, theIndex; if (theX == HUGE_VALL || theX == -HUGE_VALL) return (errno = ERANGE, 0.0L); ldX = theX; while (ldX < 0.0L) {theIndex = 2.0L; do {ldX += theIndex; theIndex *= theIndex;} while (-ldX >= theIndex);} while (ldX >= 2.0L) {theIndex = 2.0L; do {ldX -= theIndex; theIndex *= theIndex;} while (ldX >= theIndex);} if (ldX >= 1.0L) --ldX; return (ldX);} /* ATimesLog2OfX calculates theA times the base 2 logarithm of theX and returns the fractional part of this product. The integer part of the product is returned in the variable pointed to by theInt; this will be the largest integer not greater than the product. E.g., if the product is 2.1, ATimesLog2OfX will return 0.1 and 2.0; if the product is -2.1, it will return 0.9 and -3.0. */ long double ATimesLog2OfX (long double theA, long double theX, long double *theInt) {long double fractA, fractB, intA, intB; fractA = Log2OfX (theX, &intA); fractA *= theA; intA *= theA; fractB = LDFract (fractA); intB = fractA - fractB; fractA = LDFract (intA); intA -= fractA; if (fractA + fractB >= 1.0L) {--fractB; ++intB;} fractA += fractB; intA += intB; *theInt = intA; return (fractA);} /* Log2OfX returns the mantissa of the base 2 logarithm of theX. The integer part of the logarithm is returned in the variable pointed to by theInt; this will be the largest integer not greater than the logarithm. E.g., if the logarithm is 2.1, Log2OfX will return 0.1 and 2.0; if the logarithm is -2.1, it will return 0.9 and -3.0. If FLT_RADIX is 2, the mantissa returned should be accurate to LDBL_MANT_DIG bits. The value returned in the variable pointed to by theInt should be precise. */ long double Log2OfX (long double theX, long double *theInt) {long double theSum, lastSum, pSum, qSum; long double ldX, ldY, ldZ, theIndex; if (theX < 0.0L) {errno = EDOM; *theInt = NAN; return (0.0L);} if (theX == 0.0L) {errno = ERANGE; *theInt = -HUGE_VALL; return (0.0L);} if (theX == HUGE_VALL) {errno = ERANGE; *theInt = HUGE_VALL; return (0.0L);} ldX = theX; theSum = pSum = 0.0L; ldZ = 1.0L; while (ldX / ldZ < 1.0L) {theIndex = 1.0L; ldY = 2.0L; do {theSum -= theIndex; ldZ /= ldY; theIndex += theIndex; ldY *= ldY;} while (ldZ / ldX >= ldY);} while (ldX / ldZ >= 2.0L) {theIndex = 1.0L; ldY = 2.0L; do {theSum += theIndex; ldZ *= ldY; theIndex += theIndex; ldY *= ldY;} while (ldX / ldZ >= ldY);} ldX = (ldX - ldZ) / ldZ; if (ldX > 0.0L) {theIndex = 1.0L; do {theIndex /= 2.0L; ldY = ldX * ldX; ldZ = ldX + ldX; ldX = ldY + ldZ;} while (ldX < 1.0L); ldX = (ldY + --ldZ) / 2.0L; qSum = theIndex; if (ldX > 0.0L) {do {lastSum = pSum; do {theIndex /= 2.0L; ldY = ldX * ldX; ldZ = ldX + ldX; ldX = ldY + ldZ;} while (ldX < 1.0L); ldX = (ldY + --ldZ) / 2.0L; pSum += theIndex;} while (ldX > 0.0L && pSum != lastSum);} pSum += qSum;} *theInt = theSum; return (pSum);} /* TwoToTheX returns 2 to the power of a long double argument. If FLT_RADIX is 2, the value returned should be accurate to LDBL_MANT_DIG bits. */ long double TwoToTheX (long double theX) {long double ldX, theIndex, theProd, pProd; long double thePower, prevPower, lastApprox; if (theX == -HUGE_VALL) return (errno = ERANGE, 0.0L); if (theX == HUGE_VALL) return (errno = ERANGE, HUGE_VALL); ldX = theX; theProd = 1.0L; while (ldX < 0.0L) {theIndex = 1.0L; thePower = 2.0L; do {theProd /= thePower; ldX += theIndex; if (theProd == 0.0L) return (errno = ERANGE, theProd); theIndex += theIndex; thePower *= thePower;} while (-ldX >= theIndex);} while (ldX >= 1.0L) {theIndex = 1.0L; thePower = 2.0L; do {theProd *= thePower; ldX -= theIndex; if (theProd == HUGE_VALL) return (errno = ERANGE, theProd); theIndex += theIndex; thePower *= thePower;} while (ldX >= theIndex);} if (ldX > 0.0L) {pProd = 0.0L; thePower = theIndex = 1.0L; do {prevPower = thePower; thePower = 0.0L; do {lastApprox = thePower; thePower += (prevPower - thePower) / (thePower + 1.0L); thePower /= 2.0L;} while (thePower != lastApprox); theIndex /= 2.0L; if (ldX >= theIndex) {pProd += (pProd * thePower) + thePower; ldX -= theIndex;}} while (ldX > 0.0L && (thePower + 1.0L) > 1.0L); theProd += theProd * pProd; if (theProd == HUGE_VALL) errno = ERANGE;} return (theProd);} /* These macros are used in calls to HexTimesX. Cantrell's coefficients are rational numbers, but can't be precisely represented as long doubles. The other constants are irrational. The "S" at the beginning of S_SQUARE_ROOT_OF_PI and S_BASE_2_LOG_OF_E is meant to prevent possible conflicts. MANT_OF_BASE_2_LOG_SQRT_2PI is the mantissa of the base 2 log of the square root of 2pi. */ #if LDBL_DIG < 19 #define MANT_OF_BASE_2_LOG_SQRT_2PI "0.536439A4C6EFBAD86E" #define S_SQUARE_ROOT_OF_PI "1.C5BF891B4EF6AA79C" #define S_BASE_2_LOG_OF_E "1.71547652B82FE1778" #define CANTRELL_COEFFICIENT_EIGHT "1.8F27478402B14D893" #define CANTRELL_COEFFICIENT_SEVEN "0.28E791360B93BC5F5D" #define CANTRELL_COEFFICIENT_SIX "2.0D5E32BC6D9F5D78E" #define CANTRELL_COEFFICIENT_FIVE "0.372642BE5EE71DC552" #define CANTRELL_COEFFICIENT_FOUR "2.FDCDC5670AA34302A" #define CANTRELL_COEFFICIENT_THREE "0.54C0965CAE8514023D" #define CANTRELL_COEFFICIENT_TWO "5.734CE827B8D8E9839" #define CANTRELL_COEFFICIENT_ONE "0.B6679DB440592CF7C4" #else #define MANT_OF_BASE_2_LOG_SQRT_2PI "0.536439A4C6EFBAD86D9727840544713A5C" #define S_SQUARE_ROOT_OF_PI "1.C5BF891B4EF6AA79C3B0520D5DB9383FE" #define S_BASE_2_LOG_OF_E "1.71547652B82FE1777D0FFDA0D23A7D11D" #define CANTRELL_COEFFICIENT_EIGHT "1.8F27478402B14D89312C473E2C34ACB43" #define CANTRELL_COEFFICIENT_SEVEN "0.28E791360B93BC5F5C9D6B186484A0EBEA" #define CANTRELL_COEFFICIENT_SIX "2.0D5E32BC6D9F5D78DF50633852B8005E1" #define CANTRELL_COEFFICIENT_FIVE "0.372642BE5EE71DC5524A4434F49AE415FB" #define CANTRELL_COEFFICIENT_FOUR "2.FDCDC5670AA343029D7681AE3E727CB38" #define CANTRELL_COEFFICIENT_THREE "0.54C0965CAE8514023D78641D8EE67C0D33" #define CANTRELL_COEFFICIENT_TWO "5.734CE827B8D8E983955383CB54E551F73" #define CANTRELL_COEFFICIENT_ONE "0.B6679DB440592CF7C43636F98A326A2641" #endif /* Log2OfGamma returns the mantissa of the base 2 logarithm of the (complete) Gamma function of theX. The integer part of the logarithm is returned in the variable pointed to by theInt; this will be the largest integer not greater than the logarithm. E.g., if the logarithm is 2.1, Log2OfGamma will return 0.1 and 2.0; if the logarithm is -0.1, it will return 0.9 and -1.0. Log2OfGamma uses an approximation developed by David W. Cantrell; thanks to Hugo Pfoertner for pointing me to it. See: http://mathforum.org/epigone/sci.math.num-analysis/zonchahee/ */ long double Log2OfGamma (long double theX, long double *theInt) {long double theValue, ldX, thePochhammer; long double fractA, intA, fractB, intB; if (theX < 0.0L) {errno = EDOM; *theInt = NAN; return (0.0L);} if (theX == 0.0L || theX == HUGE_VALL) {errno = ERANGE; *theInt = HUGE_VALL; return (0.0L);} ldX = theX; /* If it's easy, calculate Gamma the obvious way. */ if (ldX <= 9.0L && LDFract (ldX * 2.0L) == 0.0L) {theValue = 1.0L; while (--ldX > 0.0L) theValue *= ldX; if (ldX < 0.0L) theValue = HexTimesX (S_SQUARE_ROOT_OF_PI, theValue); return (Log2OfX (theValue, theInt));} /* Otherwise, use Cantrell's approximation. */ thePochhammer = 1.0L; while (ldX <= 9.0L) thePochhammer *= ldX++; ldX -= 0.5L; theValue = HexTimesX (CANTRELL_COEFFICIENT_EIGHT, ldX); theValue = HexTimesX (CANTRELL_COEFFICIENT_SEVEN, ldX) + (1.0L / theValue); theValue = HexTimesX (CANTRELL_COEFFICIENT_SIX, ldX) + (1.0L / theValue); theValue = HexTimesX (CANTRELL_COEFFICIENT_FIVE, ldX) + (1.0L / theValue); theValue = HexTimesX (CANTRELL_COEFFICIENT_FOUR, ldX) + (1.0L / theValue); theValue = HexTimesX (CANTRELL_COEFFICIENT_THREE, ldX) + (1.0L / theValue); theValue = HexTimesX (CANTRELL_COEFFICIENT_TWO, ldX) + (1.0L / theValue); theValue = HexTimesX (CANTRELL_COEFFICIENT_ONE, ldX) + (1.0L / theValue); theValue = (24.0L * ldX) + (1.0L / theValue) - 0.5L; theValue = (1.0L / (1.0L + (1.0L / theValue))) / thePochhammer; fractA = Log2OfX (theValue, &intA); fractB = HexTimesX (MANT_OF_BASE_2_LOG_SQRT_2PI, 1.0L); intB = 1.0L; if (fractA + fractB >= 1.0L) {--fractB; ++intB;} fractA += fractB; intA += intB; fractB = ATimesLog2OfX (ldX, ldX, &intB); if (fractA + fractB >= 1.0L) {--fractB; ++intB;} fractA += fractB; intA += intB; intB = HexTimesX (S_BASE_2_LOG_OF_E, ldX); fractB = LDFract (intB); intB -= fractB; if (fractA - fractB < 0.0L) {--fractB; ++intB;} fractA -= fractB; intA -= intB; theValue = fractA; *theInt = intA; if (intA == HUGE_VALL) errno = ERANGE; return (theValue);} /* IncompleteGammaFactor is called by the incomplete gamma functions. Where "^" means "to the power of," and "Gamma" stands for an upper-case Greek letter gamma, this function calculates (x^a)/(Gamma(a)(e^x)). */ long double IncompleteGammaFactor (long double theA, long double theX) {long double fractA, fractB, intA, intB; fractA = ATimesLog2OfX (theA, theX, &intA); fractB = Log2OfGamma (theA, &intB); if (fractA - fractB < 0.0L) {--fractB; ++intB;} fractA -= fractB; intA -= intB; intB = HexTimesX (S_BASE_2_LOG_OF_E, theX); fractB = LDFract (intB); intB -= fractB; if (fractA - fractB < 0.0L) {--fractB; ++intB;} fractA -= fractB; intA -= intB; return (TwoToTheX (fractA) * TwoToTheX (intA));} /* Don't call IncompleteGammaBySeries directly; call IncompleteGamma. */ double IncompleteGammaBySeries (double theA, double theX) {long double theSum, lastSum, aValue, theProd; long double ldA, ldX; ldA = theA; ldX = theX; aValue = ldA; theSum = theProd = 1.0L / aValue; do {lastSum = theSum; theProd /= ++aValue; theProd *= ldX; theSum += theProd;} while (theSum != lastSum); return (1.0L - (theSum * IncompleteGammaFactor (ldA, ldX)));} /* Don't call IncompleteGammaByFraction directly; call IncompleteGamma. */ double IncompleteGammaByFraction (double theA, double theX) {long double aValue, theB, theC, theD, ldA, ldX; long double theIndex, theProd, lastProd; ldA = theA; ldX = theX; theB = 1.0L + (ldX - ldA); theD = theProd = 1.0L / theB; aValue = ldA - 1.0L; theB += 2.0L; theC = theB; theD = 1.0L / (theB + (aValue * theD)); theProd *= theC * theD; theIndex = 1.0L; do {lastProd = theProd; ++theIndex; aValue = theIndex * (ldA - theIndex); theB += 2.0L; theC = theB + (aValue / theC); theD = 1.0L / (theB + (aValue * theD)); theProd *= theC * theD;} while (theProd != lastProd); return (theProd * IncompleteGammaFactor (ldA, ldX));} /* IncompleteGamma calculates the regularized incomplete gamma function. Where "Gamma" stands for an upper-case Greek letter gamma, IncompleteGamma calculates Gamma(a,x)/Gamma(a). W. H. Press, et al (in Numerical Recipes in C) write this function as Q(a,x). It is the complement of P(a,x); i.e. Q(a,x) = 1 - P(a,x). Though the functions here are not from Numerical Recipes in C, they do implement formulae (6.2.5) and (6.2.7) from that book. Cf. formulae (2) and (4) in "Computing the incomplete Gamma function to arbitrary precision" (2003: LNCS 2667) by Serge Winitzki. The continued fraction expansion [Press, et al (6.2.6) and Winitzki (4)] is due to Legendre. */ double IncompleteGamma (double theA, double theX) {if (theA <= 0.0 || theX < 0.0) return (errno = EDOM, 0.0); if (theA == HUGE_VAL) return (errno = ERANGE, theX / HUGE_VAL); if (theX == HUGE_VAL) return (0.0); if (theX == 0.0) return (1.0); if (theX < theA) return (IncompleteGammaBySeries (theA, theX)); return (IncompleteGammaByFraction (theA, theX));} /* ChiSquaredPValue returns the probability of a chi-squared statistic greater than or equal to theChiSquared. */ double ChiSquaredPValue (double degOfFreedom, double theChiSquared) {return (IncompleteGamma (degOfFreedom / 2.0, theChiSquared / 2.0));} /* Here ends this small library of math functions. */ /* Stopwatch functions ------------------- Imagine a stopwatch with three buttons: reset, start, and stop. Time on the watch is displayed in seconds, accurate to 1/CLOCKS_PER_SEC seconds (assuming that clock() and CLOCKS_PER_SEC are accurate). ResetStopwatch must be called before timing. It stops the watch and resets it to 0. StartStopwatch starts the watch moving; it doesn't reset it to 0. If the watch is already moving, StartStopwatch does nothing. StopStopwatch stops the watch; it doesn't reset it to 0. If the watch is already stopped, StopStopwatch does nothing. ReadStopwatch returns the current time on the watch, whether the watch is moving or stopped. Reading will return a wrong result if the time between a start and stop (or between reads) exceeds the range of clock(). A work-around for this problem is to read periodically, even if the time read is ignored. These functions assume that clock_t is an unsigned type. */ void ResetStopwatch (WatchStruct *theWatch) {theWatch->accumTime = 0.0; theWatch->watchMoving = 0; return;} void StartStopwatch (WatchStruct *theWatch) {if (!theWatch->watchMoving) {theWatch->startTime = clock (); theWatch->watchMoving = 1;} return;} void StopStopwatch (WatchStruct *theWatch) {clock_t clockDif; if (theWatch->watchMoving) {clockDif = clock () - theWatch->startTime; theWatch->accumTime += clockDif; theWatch->watchMoving = 0;} return;} double ReadStopwatch (WatchStruct *theWatch) {clock_t clockDif, newStart; if (theWatch->watchMoving) {newStart = clock (); clockDif = newStart - theWatch->startTime; theWatch->accumTime += clockDif; theWatch->startTime = newStart;} return (theWatch->accumTime / CLOCKS_PER_SEC);} void RC4SetKey (RC4EnvStruct *theEnv, const void *theKey, size_t theKeyLength) {size_t theSwap, indexA, indexB, keyEnd, expKey [256]; const OneByteInt *keyPtr; indexA = 0x00; do theEnv->permArray [indexA] = indexA; while (++indexA < 0x100); if (theKeyLength < 0xFF) keyEnd = theKeyLength - 0x01; else keyEnd = 0xFF; indexA = indexB = 0x00; keyPtr = theKey; do {expKey [indexA] = keyPtr [indexB]; if (++indexB > keyEnd) indexB = 0x00;} while (++indexA < 0x100); indexA = indexB = 0x00; do {indexB = (indexB + theEnv->permArray [indexA] + expKey [indexA]) & 0xFF; theSwap = theEnv->permArray [indexA]; theEnv->permArray [indexA] = theEnv->permArray [indexB]; theEnv->permArray [indexB] = theSwap;} while (++indexA < 0x100); theEnv->indexA = 0x01; theEnv->indexB = 0x00; return;} void RC4 (RC4EnvStruct *theEnv, void *theData, size_t theDataLength) {register size_t theIndex; register OneByteInt *theByte; register size_t byteA, byteB, indexA, indexB; theByte = theData; theIndex = theDataLength; indexA = theEnv->indexA; indexB = theEnv->indexB; do {indexA = (indexA + 1) & 0xFF; byteA = theEnv->permArray [indexA]; indexB = (indexB + byteA) & 0xFF; byteB = theEnv->permArray [indexB]; theEnv->permArray [indexA] = byteB; theEnv->permArray [indexB] = byteA; *(theByte++) ^= theEnv->permArray [(byteA + byteB) & 0xFF];} while (--theIndex); theEnv->indexA = indexA; theEnv->indexB = indexB; return;} #if TEST_SNOW2 void Snow2 (void *theData, size_t theDataLength) {size_t theIndex; unsigned long *theLong; theLong = theData; theIndex = theDataLength / sizeof (*theLong); do *(theLong++) ^= snow_keystream (); while (--theIndex); return;} #endif