📄 main.cpp
字号:
r = 32; /* COMPUTE PROBABILITIES */
product = 1;
for( i=0; i<=r-1; i++ )
product *= ((1.e0-pow(2,i-32))*(1.e0-pow(2,i-32)))/(1.e0-pow(2,i-r));
p_32 = pow(2,r*(32+32-r)-32*32) * product;
r = 31;
product = 1;
for( i=0; i<=r-1; i++ )
product *= ((1.e0-pow(2,i-32))*(1.e0-pow(2,i-32)))/(1.e0-pow(2,i-r));
p_31 = pow(2,r*(32+32-r)-32*32) * product;
p_30 = 1 - (p_32+p_31);
F_32 = 0;
F_31 = 0;
for( k=0; k<N; k++ ) { /* FOR EACH 32x32 MATRIX */
def_matrix(32, 32, matrix, k);
if ( MATRICES )
display_matrix(32, 32, matrix);
R = computeRank(32, 32, matrix);
if ( R == 32 )
F_32++; /* DETERMINE FREQUENCIES */
if ( R == 31 )
F_31++;
}
F_30 = (double)N - (F_32+F_31);
chi_squared =(pow(F_32 - N*p_32,2)/(double)(N*p_32) +
pow(F_31 - N*p_31,2)/(double)(N*p_31) +
pow(F_30 - N*p_30,2)/(double)(N*p_30));
arg1 = -chi_squared/2.e0;
/*if ( RANK ) {
fprintf(stats[TESTS_RANK], "\t\t\t\tRANK TEST\n");
fprintf(stats[TESTS_RANK], "\t\t---------------------------------------------\n");
fprintf(stats[TESTS_RANK], "\t\tCOMPUTATIONAL INFORMATION:\n");
fprintf(stats[TESTS_RANK], "\t\t---------------------------------------------\n");
fprintf(stats[TESTS_RANK], "\t\t(a) Probability P_%d = %f\n", 32,p_32);
fprintf(stats[TESTS_RANK], "\t\t(b) P_%d = %f\n", 31,p_31);
fprintf(stats[TESTS_RANK], "\t\t(c) P_%d = %f\n", 30,p_30);
fprintf(stats[TESTS_RANK], "\t\t(d) Frequency F_%d = %d\n", 32,(int)F_32);
fprintf(stats[TESTS_RANK], "\t\t(e) F_%d = %d\n", 31,(int)F_31);
fprintf(stats[TESTS_RANK], "\t\t(f) F_%d = %d\n", 30,(int)F_30);
fprintf(stats[TESTS_RANK], "\t\t(g) # of matrices = %d\n", N);
fprintf(stats[TESTS_RANK], "\t\t(h) Chi^2 = %f\n", chi_squared);
if ( n%(32*32) != 0 )
fprintf(stats[TESTS_RANK], "\t\t(i) NOTE: %d BITS WERE DISCARDED.\n", n%(32*32));
fprintf(stats[TESTS_RANK], "\t\t---------------------------------------------\n");
}*/
p_value = exp(arg1);
if ( isNegative(p_value) || isGreaterThanOne(p_value) )
printf("WARNING: P_VALUE IS OUT OF RANGE.\n");
for( i=0; i<32; i++ ) /* DEALLOCATE MATRIX */
free(matrix[i]);
free(matrix);
if ( p_value < ALPHA ) {
return false;
}
else {
return true;
}
}
/*fprintf(stats[TESTS_RANK], "%s\t\tp_value = %f\n\n", assignment, p_value); fflush(stats[TESTS_RANK]);
fprintf(results[TESTS_RANK], "%f\n", p_value); fflush(results[TESTS_RANK]);
fprintf(grid, "%d", state); fflush(grid);
fprintf(pvals, "%f ", p_value); fflush(pvals);
if ( p_value < tp.minimumP )
tp.minimumP = p_value;
if ( !_isnan(p_value) )
tp.lnSum += log(p_value);
tp.df++;
#ifdef GEN_TIMING_INFO
finish = clock();
fprintf(fp, "%d\n", finish - start);
fclose(fp);
#endif
return;*/
}
/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
R U N S T E S T
* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
bool Runs(int n)
{
int i, *r;
double argument, pi, V_n_obs, tau;
double p_value, product, sum;
// char assignment[7];
/*#ifdef GEN_TIMING_INFO
clock_t start, finish;
FILE *fp;
fp = fopen("Runs.txt", "a");
start = clock();
#endif
*/
if ( (r = (int *)calloc(n, sizeof(int))) == NULL ) {
printf("\t\tRUNS TEST: Insufficient space for work array,");
printf(" test aborted!\n");
//printf("0.000000\n"); fflush(results[TESTS_RUNS]);
//printf(pvals, "0.000000 "); fflush(pvals);
}
else {
sum = 0.0;
for(i = 0; i < n; i++)
sum += (int)bit[i];
pi = sum/n;
tau = 2.0/sqrt(n);
if (fabs(pi - 0.5) < tau) {
for(i = 0; i < n-1; i++) {
if ((int)bit[i] == (int)bit[i+1])
r[i] = 0;
else
r[i] = 1;
}
V_n_obs = 0;
for(i = 0; i < n-1; i++)
V_n_obs += r[i];
V_n_obs++;
product = pi * (1.e0 - pi);
argument = fabs(V_n_obs - 2.e0*n*product)/(2.e0*sqrt(2.e0*n)*product);
p_value = erfc(argument);
/* if ( RUNS ) {
fprintf(stats[TESTS_RUNS], "\t\t\t\tRUNS TEST\n");
fprintf(stats[TESTS_RUNS], "\t\t------------------------------------------\n");
fprintf(stats[TESTS_RUNS], "\t\tCOMPUTATIONAL INFORMATION:\n");
fprintf(stats[TESTS_RUNS], "\t\t------------------------------------------\n");
fprintf(stats[TESTS_RUNS], "\t\t(a) Pi = %f\n", pi);
fprintf(stats[TESTS_RUNS], "\t\t(b) V_n_obs (Total # of runs) = %d\n", (int)V_n_obs);
fprintf(stats[TESTS_RUNS], "\t\t(c) V_n_obs - 2 n pi (1-pi)\n");
fprintf(stats[TESTS_RUNS], "\t\t ----------------------- = %f\n", argument);
fprintf(stats[TESTS_RUNS], "\t\t 2 sqrt(2n) pi (1-pi)\n");
fprintf(stats[TESTS_RUNS], "\t\t------------------------------------------\n");
}
*/
if ( isNegative(p_value) || isGreaterThanOne(p_value) )
printf("WARNING: P_VALUE IS OUT OF RANGE.\n");
free(r);
if ( p_value < ALPHA ) {
return false;
}
else {
return true;
}
// fprintf(stats[TESTS_RUNS],"%s\t\tp_value = %f\n\n", assignment, p_value);
}
else {
if (RUNS) {
printf("\t\t\t\tRUNS TEST\n");
printf("\t\t------------------------------------------\n");
printf("\t\tPI ESTIMATOR CRITERIA NOT MET! PI = %1e\n", pi);
}
//strcpy(assignment,"REJECTION");
//p_value = 0.0;
//state = 0;
return false;
}
//fprintf(results[TESTS_RUNS], "%f\n", p_value); fflush(results[TESTS_RUNS]);
//fprintf(grid, "%d", state); fflush(grid);
//fprintf(pvals, "%f ", p_value); fflush(pvals);
/*if ( p_value < tp.minimumP )
tp.minimumP = p_value;
if ( !_isnan(p_value) )
tp.lnSum += log(p_value);
tp.df++;
fflush(stats[TESTS_RUNS]);
*/
}
/*#ifdef GEN_TIMING_INFO
finish = clock();
fprintf(fp, "%d\n", finish - start);
fclose(fp);
#endif
return;*/
}
/*******************************************************************************
串行测试
××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××××8*****/
double psi2(int m, int n);
bool Serial(int m, int n)
{
//int state;
//char assignment[7];
double p_value1, p_value2, psim0, psim1, psim2, del1, del2;
/*#ifdef GEN_TIMING_INFO
clock_t start, finish;
FILE *fp;
fp = fopen("Serial.txt", "a");
start = clock();
#endif
*/
psim0 = psi2(m, n);
psim1 = psi2(m-1, n);
psim2 = psi2(m-2, n);
del1 = psim0 - psim1;
del2 = psim0 - 2.0*psim1 + psim2;
p_value1 = igamc(pow(2,m-1)/2,del1/2.0);
p_value2 = igamc(pow(2,m-2)/2,del2/2.0);
/* fprintf(stats[TESTS_SERIAL], "\t\t\t SERIAL TEST\n");
fprintf(stats[TESTS_SERIAL], "\t\t---------------------------------------------\n");
fprintf(stats[TESTS_SERIAL], "\t\t COMPUTATIONAL INFORMATION: \n");
fprintf(stats[TESTS_SERIAL], "\t\t---------------------------------------------\n");
fprintf(stats[TESTS_SERIAL], "\t\t(a) Block length (m) = %d\n", m);
fprintf(stats[TESTS_SERIAL], "\t\t(b) Sequence length (n) = %d\n", n);
fprintf(stats[TESTS_SERIAL], "\t\t(c) Psi_m = %f\n", psim0);
fprintf(stats[TESTS_SERIAL], "\t\t(d) Psi_m-1 = %f\n", psim1);
fprintf(stats[TESTS_SERIAL], "\t\t(e) Psi_m-2 = %f\n", psim2);
fprintf(stats[TESTS_SERIAL], "\t\t(f) Del_1 = %f\n", del1);
fprintf(stats[TESTS_SERIAL], "\t\t(g) Del_2 = %f\n", del2);
fprintf(stats[TESTS_SERIAL], "\t\t---------------------------------------------\n");
*/
if ( (p_value1 < ALPHA) || (p_value2 < ALPHA) ) {
return false;
}
else {
return true;
}
/* fprintf(stats[TESTS_SERIAL], "%s\t\tp_value1 = %f\n", assignment, p_value1);
fprintf(results[TESTS_SERIAL], "%f\n", p_value1); fflush(results[TESTS_SERIAL]);
fprintf(grid, "%d", state); fflush(grid);
fprintf(pvals, "%f ", p_value1); fflush(pvals);
*/
/* if ( p_value1 < tp.minimumP )
tp.minimumP = p_value1;
if ( !_isnan(p_value1) )
tp.lnSum += log(p_value1);
tp.df++;*/
/* fprintf(stats[TESTS_SERIAL], "%s\t\tp_value2 = %f\n\n", assignment, p_value2); fflush(stats[TESTS_SERIAL]);
fprintf(results[TESTS_SERIAL], "%f\n", p_value2); fflush(results[TESTS_SERIAL]);
fprintf(grid, "%d", state); fflush(grid);
fprintf(pvals, "%f ", p_value2); fflush(pvals);
if ( p_value2 < tp.minimumP )
tp.minimumP = p_value2;
if ( !_isnan(p_value2) )
tp.lnSum += log(p_value2);
tp.df++;
#ifdef GEN_TIMING_INFO
finish = clock();
fprintf(fp, "%d\n", finish - start);
fclose(fp);
#endif
return;*/
}
double psi2(int m, int n)
{
int i, j, k, powLen;
double sum, numOfBlocks;
unsigned int* P;
if ( (m == 0) || (m == -1) )
return 0.0;
numOfBlocks = n;
powLen = (int)pow(2,m+1)-1;
if ( (P = (unsigned int*)calloc(powLen,sizeof(unsigned int)))== NULL ) {
printf("Serial Test: Insufficient memory available.\n");
//fflush(stats[TESTS_SERIAL]);
return 0;
}
for( i=1; i<powLen-1; i++ )
P[i] = 0; /* INITIALIZE NODES */
for( i=0; i<numOfBlocks; i++ ) { /* COMPUTE FREQUENCY */
k = 1;
for( j=0; j<m; j++ ) {
if ( bit[(i+j)%n] == 0 )
k *= 2;
else if ( bit[(i+j)%n] == 1 )
k = 2*k+1;
}
P[k-1]++;
}
/* DISPLAY FREQUENCY */
sum = 0.0;
for( i=(int)pow(2,m)-1; i<(int)pow(2,m+1)-1; i++ )
sum += pow(P[i],2);
sum = (sum * pow(2,m)/(double)n) - (double)n;
free(P);
return sum;
}
/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
U N I V E R S A L T E S T
* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
bool Universal(int L, int Q, int n)
{
int i, j, p, K;
double arg, sqrt2, sigma, phi, sum, p_value, c;
long* T, decRep;
// char assignment[7];
double expected_value[17] = {0, 0, 0, 0, 0, 0, 5.2177052, 6.1962507,
7.1836656, 8.1764248, 9.1723243, 10.170032,
11.168765, 12.168070, 13.167693, 14.167488,
15.167379};
double variance[17] = {0, 0, 0, 0, 0, 0, 2.954, 3.125, 3.238, 3.311,
3.356, 3.384, 3.401, 3.410, 3.416, 3.419, 3.421};
/*#ifdef GEN_TIMING_INFO
clock_t start, finish;
FILE *fp;
fp = fopen("Universal.txt", "a");
start = clock();
#endif
*/
/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
* THE FOLLOWING REDEFINES L, SHOULD THE CONDITION: n >= 1010*2^L*L *
* NOT BE MET, FOR THE BLOCK LENGTH L. *
* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
if ( n >= 387840 ) L = 6;
if ( n >= 904960 ) L = 7;
if ( n >= 2068480 ) L = 8;
if ( n >= 4654080 ) L = 9;
if ( n >= 10342400 ) L = 10;
if ( n >= 22753280 ) L = 11;
if ( n >= 49643520 ) L = 12;
if ( n >= 107560960 ) L = 13;
if ( n >= 231669760 ) L = 14;
if ( n >= 496435200 ) L = 15;
if ( n >= 1059061760 ) L = 16;
Q = 10*(int)pow(2,L);
K = (int) (floor(n/L) - (double)Q); /* BLOCKS TO TEST */
if ( (L < 6) || (L > 16) || ((double)Q < 10*pow(2,L)) ) {
printf("\t\tUNIVERSAL STATISTICAL TEST\n");
printf("\t\t---------------------------------------------\n");
if ( (L < 6) || (L > 16) )
printf("\t\tERROR: L IS OUT OF RANGE.\n");
else
printf("\t\tERROR: Q IS LESS THAN %f.\n", 10*pow(2,L));
//printf(pvals, "0.000000 "); fflush(pvals);
//fprintf(results[TESTS_UNIVERSAL], "0.000000\n"); fflush(results[TESTS_UNIVERSAL]);
}
else {
/* COMPUTE THE EXPECTED: Formula 16, in Marsaglia's Paper */
c = 0.7 - 0.8/(double)L + (4 + 32/(double)L)*pow(K,-3/(double)L)/15;
sigma = c * sqrt(variance[L]/(double)K);
sqrt2 = sqrt(2);
sum = 0.0;
p = (int)pow(2,L);
T = (long*) calloc(p, sizeof(long));
for( i=0; i<p; i++ )
T[i] = 0;
for( i=1; i<=Q; i++ ) { /* INITIALIZE TABLE */
decRep = 0;
for( j=0; j<L; j++ )
decRep += bit[(i-1)*L+j] * (long)pow(2,L-1-j);
T[decRep] = i;
}
for( i=Q+1; i<=Q+K; i++ ) { /* PROCESS BLOCKS */
decRep = 0;
for( j=0; j<L; j++ )
decRep += bit[(i-1)*L+j] * (long)pow(2,L-1-j);
sum += log(i - T[decRep])/log(2);
T[decRep] = i;
}
phi = (double)(sum/(double)K);
/*if ( UNIVERSAL ) {
fprintf(stats[TESTS_UNIVERSAL], "\t\tUNIVERSAL STATISTICAL TEST\n");
fprintf(stats[TESTS_UNIVERSAL], "\t\t--------------------------------------------\n");
fprintf(stats[TESTS_UNIVERSAL], "\t\tCOMPUTATIONAL INFORMATION:\n");
fprintf(stats[TESTS_UNIVERSAL], "\t\t--------------------------------------------\n");
fprintf(stats[TESTS_UNIVERSAL], "\t\t(a) L = %d\n", L);
fprintf(stats[TESTS_UNIVERSAL], "\t\t(b) Q = %d\n", Q);
fprintf(stats[TESTS_UNIVERSAL], "\t\t(c) K = %d\n", K);
fprintf(stats[TESTS_UNIVERSAL], "\t\t(d) sum = %f\n", sum);
fprintf(stats[TESTS_UNIVERSAL], "\t\t(e) sigma = %f\n", sigma);
fprintf(stats[TESTS_UNIVERSAL], "\t\t(f) variance = %f\n", variance[L]);
fprintf(stats[TESTS_UNIVERSAL], "\t\t(g) exp_value = %f\n", expected_value[L]);
fprintf(stats[TESTS_UNIVERSAL], "\t\t(h) phi = %f\n", phi);
if ( ((Q+K)*L) < n )
fprintf(stats[TESTS_UNIVERSAL], "\t\t(i) WARNING: %d bits were discarded.\n", n-(Q+K)*L);
fprintf(stats[TESTS_UNIVERSAL], "\t\t-----------------------------------------\n");
}*/
if ( ((Q+K)*L) < n )
printf("\t\t(i) WARNING: %d bits were discarded.\n", n-(Q+K)*L);
arg = fabs(phi-expected_value[L])/(sqrt2 * sigma);
p_value = erfc(arg);
if ( isNegative(p_value) || isGreaterThanOne(p_value) )
printf("\t\tWARNING: P_VALUE IS OUT OF RANGE\n");
free(T);
if ( p_value < ALPHA ) { /* INTERPRETATION */
return false;
}
else {
return true;
}
//fprintf(stats[TESTS_UNIVERSAL], "%s\t\tp_value = %f\n\n", assignment, p_value); fflush(stats[TESTS_UNIVERSAL]);
//fprintf(results[TESTS_UNIVERSAL], "%f\n", p_value); fflush(results[TESTS_UNIVERSAL]);
//fprintf(grid, "%d", state); fflush(grid);
//fprintf(pvals, "%f ", p_value); fflush(pvals);
/*if ( p_value < tp.minimumP )
tp.minimumP = p_value;
if ( !_isnan(p_value) )
tp.lnSum += log(p_value);
tp.df++;*/
}
/*#ifdef GEN_TIMING_INFO
finish = clock();
fprintf(fp, "%d\n", finish - start);
fclose(fp);
#endif
return;*/
}
int main(int argc, char* argv[])
{
int max,num[6];
bool bTest;
double x;
printf("测试名称后面的1表示通过测试，0表示出错，或者没有通过测试。\n");
Getdata();
printf("--------------------------------------------------------------------");
printf("\t\tMonobitTest:\n");
bTest=MonobitTest(num[0]);
printf("MonobitTest:%d\t%d\n", bTest, num[0]);
printf("--------------------------------------------------------------------");
printf("\t\tPokerTest:\n");
bTest=PokerTest(x);
printf("PokerTest:%d\t%f\n", bTest, x);
printf("--------------------------------------------------------------------");
printf("\t\tRunsTest:\n");
bTest=RunsTest(num);
printf("RunsTest:%d\t", bTest);
for(int i=0; i<6; i++)
printf("%d\t", num[i]);
printf("\n");
printf("-------------------------------------⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -