📄 main.cpp
字号:
else {
return true;
}
}
}
bool LinearComplexity(int M, int n)
{
int i, ii, j, d, N;
int L, m, N_, parity, sign;
double p_value, T_, mean;
int K = 6;
double pi[7]={0.01047,0.03125,0.12500,0.50000,0.25000,0.06250,0.020833};
double nu[7], chi2;
BitField *T, *P, *B_, *C;
N = (int)floor(n/M);
if ( ((B_ = (BitField*) calloc(M,sizeof(BitField))) == NULL) ||
((C = (BitField*) calloc(M,sizeof(BitField))) == NULL) ||
((P = (BitField*) calloc(M,sizeof(BitField))) == NULL) ||
((T = (BitField*) calloc(M,sizeof(BitField))) == NULL) ) {
printf("Insufficient Memory for Work Space:: Linear Complexity Test\n");
fflush(stdout);
if ( B_!= NULL )
free(B_);
if ( C != NULL )
free(C);
if ( P != NULL )
free(P);
if ( T != NULL )
free(T);
return false;
}
if ( n%M != 0 )
printf("\tNote: %d bits were discarded!\n", n%M);
for( i=0; i<K+1; i++ )
nu[i] = 0.00;
for( ii=0; ii<N; ii++ ) {
for( i=0; i<M; i++ ) {
B_[i].b = 0;
C[i].b = 0;
T[i].b = 0;
P[i].b = 0;
}
L = 0;
m = -1;
d = 0;
C[0].b = 1;
B_[0].b = 1;
/* DETERMINE LINEAR COMPLEXITY */
N_ = 0;
while( N_ < M ) {
d = (int)bit[ii*M+N_];
for( i=1; i<=L; i++ )
d += (int)C[i].b*(int)bit[ii*M+N_-i];
d = d%2;
if ( d == 1 ) {
for( i=0; i<M; i++ ) {
T[i].b = C[i].b;
P[i].b = 0;
}
for( j=0; j<M; j++ )
if ( B_[j].b == 1 )
P[j+N_-m].b = 1;
for( i=0; i<M; i++ )
C[i].b = (C[i].b + P[i].b)%2;
if ( L <= N_/2 ) {
L = N_ + 1 - L;
m = N_;
for( i=0; i<M; i++ )
B_[i].b = T[i].b;
}
}
N_++;
}
if ( (parity = (M+1)%2) == 0 )
sign = -1;
else
sign = 1;
mean = M/2. + (9.+sign)/36. - 1./pow(2,M) * (M/3. + 2./9.);
if ( (parity = M%2) == 0 )
sign = 1;
else
sign = -1;
T_ = sign * (L - mean) + 2./9.;
if ( T_ <= -2.5 )
nu[0]++;
else if ( T_ > -2.5 && T_ <= -1.5 )
nu[1]++;
else if ( T_ > -1.5 && T_ <= -0.5 )
nu[2]++;
else if ( T_ > -0.5 && T_ <= 0.5 )
nu[3]++;
else if ( T_ > 0.5 && T_ <= 1.5 )
nu[4]++;
else if ( T_ > 1.5 && T_ <= 2.5 )
nu[5]++;
else
nu[6]++;
}
chi2 = 0.00;
//for( i=0; i<K+1; i++ )
//fprintf(stats[TESTS_LINEAR_COMPLEXITY], "%4d ",(int)nu[i]);
for( i=0; i<K+1; i++ )
chi2 += pow(nu[i]-N*pi[i],2)/(N*pi[i]);
p_value = igamc(K/2.0,chi2/2.0);
free(B_);
free(P);
free(C);
free(T);
if ( p_value < ALPHA ) {
return false;
}
else {
return true;
}
}
/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
L O N G E S T R U N S T E S T
* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
bool LongestRunOfOnes(int n)
{
double run, v_n_obs, p_value, sum, chi2;
int N, i, j;
//char assignment[7];
#if LONG_RUNS_CASE_8 == 1
unsigned int nu[4] = {0, 0, 0, 0};
double pi[4] = {0.21484375, 0.3671875, 0.23046875, 0.1875};
double k[4] = {1, 2, 3, 8};
double K = 3;
int M = 8;
#elif LONG_RUNS_CASE_128 == 1
double pi[6] = {0.1174035788, 0.242955959, 0.249363483, 0.17517706,
0.1027010710, 0.112398847};
double k[6] = {4, 5, 6, 7, 8, 9};
double K = 5;
unsigned int nu[6] = {0, 0, 0, 0, 0, 0};
int M = 128;
#elif LONG_RUNS_CASE_10000 == 1
double pi[7] = {0.0882, 0.2092, 0.2483, 0.1933, 0.1208, 0.0675, 0.0727};
double K = 6;
unsigned int nu[7] = {0, 0, 0, 0, 0, 0, 0};
int M = 10000;
#endif
/*#ifdef GEN_TIMING_INFO
clock_t start, finish;
FILE *fp;
fp = fopen("LongRuns.txt", "a");
start = clock();
#endif
*/
N = (int)floor(n/M);
for(i = 0; i < N; i++) {
v_n_obs = 0.e0;
run = 0.e0;
for(j = i*M; j < (i+1)*M; j++) {
if ((int)bit[j] == 1) {
run++;
v_n_obs = MAX(v_n_obs, run);
}
else
run = 0.e0;
}
#if LONG_RUNS_CASE_8 == 1
switch((int)v_n_obs) {
case 0:
case 1:
nu[0]++;
break;
case 2:
nu[1]++;
break;
case 3:
nu[2]++;
break;
default:
nu[3]++;
break;
}
#elif LONG_RUNS_CASE_128 == 1
if ( (int)v_n_obs <= 4 )
nu[0]++;
else
switch((int)v_n_obs) {
case 5:
nu[1]++;
break;
case 6:
nu[2]++;
break;
case 7:
nu[3]++;
break;
case 8:
nu[4]++;
break;
default:
nu[5]++;
break;
}
#elif LONG_RUNS_CASE_10000 == 1
if ( (int)v_n_obs <= 10 )
nu[0]++;
else
switch((int)v_n_obs) {
case 11:
nu[1]++;
break;
case 12:
nu[2]++;
break;
case 13:
nu[3]++;
break;
case 14:
nu[4]++;
break;
case 15:
nu[5]++;
break;
default:
nu[6]++;
break;
}
#endif
}
chi2 = 0.0; /* Compute Chi Square */
sum = 0;
for(i = 0; i < K+1; i++) {
chi2 += pow((double)nu[i] - (double)N*pi[i],2)/((double)N*pi[i]);
sum += nu[i];
}
p_value = igamc(K/2.,chi2/2.);
/*p_value = gammq(K/2.,chi2/2.);*/
/* if (LONG_RUNS) {
fprintf(stats[TESTS_LONGEST_RUNS], "\t\t\t LONGEST RUNS OF ONES TEST\n");
fprintf(stats[TESTS_LONGEST_RUNS], "\t\t---------------------------------------------\n");
fprintf(stats[TESTS_LONGEST_RUNS], "\t\tCOMPUTATIONAL INFORMATION:\n");
fprintf(stats[TESTS_LONGEST_RUNS], "\t\t---------------------------------------------\n");
fprintf(stats[TESTS_LONGEST_RUNS], "\t\t(a) N (# of substrings) = %d\n", (int)N);
fprintf(stats[TESTS_LONGEST_RUNS], "\t\t(b) M (Substring Length) = %d\n", (int)M);
fprintf(stats[TESTS_LONGEST_RUNS], "\t\t(c) Chi^2 = %f\n", chi2);
fprintf(stats[TESTS_LONGEST_RUNS], "\t\t---------------------------------------------\n");
fprintf(stats[TESTS_LONGEST_RUNS], "\t\t F R E Q U E N C Y\n");
fprintf(stats[TESTS_LONGEST_RUNS], "\t\t---------------------------------------------\n");
}
if (LONG_RUNS && LONG_RUNS_CASE_8) {
fprintf(stats[TESTS_LONGEST_RUNS], "\t\t <=1 2 3 >=4 P-value Assignment");
fprintf(stats[TESTS_LONGEST_RUNS], "\n\t\t %3d %3d %3d %3d ",nu[0],nu[1],nu[2],nu[3]);
}
else if (LONG_RUNS && LONG_RUNS_CASE_128) {
fprintf(stats[TESTS_LONGEST_RUNS], "\t\t<=4 5 6 7 8 >=9 P-value Assignment");
fprintf(stats[TESTS_LONGEST_RUNS], "\n\t\t %3d %3d %3d %3d %3d %3d ",nu[0],nu[1],nu[2], nu[3], nu[4], nu[5]);
}
else if (LONG_RUNS && LONG_RUNS_CASE_10000) {
fprintf(stats[TESTS_LONGEST_RUNS], "\t\t<=10 11 12 13 14 15 >=16 P-value Assignment");
fprintf(stats[TESTS_LONGEST_RUNS], "\n\t\t %3d %3d %3d %3d %3d %3d %3d ",nu[0],nu[1],nu[2], nu[3], nu[4], nu[5], nu[6]);
}
if ( isNegative(p_value) || isGreaterThanOne(p_value) )
fprintf(stats[TESTS_LONGEST_RUNS], "WARNING: P_VALUE IS OUT OF RANGE.\n");
*/
if ( p_value < ALPHA ) {
return false;
}
else {
return true;
}
/*fprintf(stats[TESTS_LONGEST_RUNS], "%f %s\n\n", p_value, assignment); fflush(stats[TESTS_LONGEST_RUNS]);
fprintf(results[TESTS_LONGEST_RUNS], "%f\n", p_value); fflush(results[TESTS_LONGEST_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++;
#ifdef GEN_TIMING_INFO
finish = clock();
fprintf(fp, "%d\n", finish - start);
fclose(fp);
#endif
return;*/
}
/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
N O N O V E R L A P P I N G T E M P L A T E T E S T
* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
bool NonOverlappingTemplateMatchings(int m, int n) //m是模板的个数
{
int numOfTemplates[22] = { 0, 0, 2, 4, 6, 12, 20,
40, 74, 148, 284, 568, 1116, 2232,
4424, 8848, 17622, 35244,70340, 140680, 281076, 562152};
/*----------------------------------------------------------------------------*
*NOTE: Should additional templates lengths beyond 21 be desired, they must *
*first be constructed, saved into files and then the corresponding *
*number of nonperiodic templates for that file be stored in the m-th *
*position in the numOfTemplates variable. *
*----------------------------------------------------------------------------*/
int SKIP;
unsigned bits;
FILE* fp;
unsigned int W_obs;
double sum, chi2, p_value, lambda;
int i, j, jj, k, match, state;
char directory[100];
int M, N, K = 5;
unsigned int nu[6], *Wj;
double pi[6], varWj /*, U[6] = {1,2,3,4,5,6}*/;
extern double lgam ( double x );
BitField *templates=0;
/*
#ifdef GEN_TIMING_INFO
clock_t start, finish;
FILE *fp1;
fp1 = fopen("NonOverlap.txt", "a");
start = clock();
#endif
*/
N = 8;
M = (int)floor(n/N);
Wj = (unsigned int*)calloc(N,sizeof(unsigned int));
lambda = (M-m+1)/pow(2,m);
varWj = M*(1./pow(2.,m) - (2.*m-1.)/pow(2.,2.*m));
sprintf(directory, "template%d",m);
/*printf("%s\n", directory);*/
if ( ((isNegative(lambda)) || (isZero(lambda))) ||
((fp = fopen(directory, "r")) == NULL) ||
((templates = (BitField*) calloc(m,sizeof(BitField))) == NULL)) {
printf("\tNONOVERLAPPING TEMPLATES TESTS ABORTED DUE ");
printf("TO ONE OF THE FOLLOWING : \n");
printf("\tLambda (%f) not being positive!\n", lambda);
printf("\tTemplate file %s not existing\n", directory);
printf("\tInsufficient memory for required work space.\n");
if ( templates )
free(templates);
state = 0;
exit(-1);
}
else {
/* if ( APERIODIC_TEMPLATES ) {
fprintf(stats[TESTS_NONPERIODIC_TEMPL], "\t\t NONPERIODIC TEMPLATES TEST\n");
fprintf(stats[TESTS_NONPERIODIC_TEMPL], "---------------------------------------------");
fprintf(stats[TESTS_NONPERIODIC_TEMPL], "-----------------------------------\n");
fprintf(stats[TESTS_NONPERIODIC_TEMPL], "\t\t COMPUTATIONAL INFORMATION\n");
fprintf(stats[TESTS_NONPERIODIC_TEMPL], "---------------------------------------------");
fprintf(stats[TESTS_NONPERIODIC_TEMPL], "-----------------------------------\n");
fprintf(stats[TESTS_NONPERIODIC_TEMPL], "\tLAMBDA = %f\tM = %d\tN = %d\tm = %d\tn = %d\n", lambda,M,N,m,n);
fprintf(stats[TESTS_NONPERIODIC_TEMPL], "---------------------------------------------");
fprintf(stats[TESTS_NONPERIODIC_TEMPL], "-----------------------------------\n");
fprintf(stats[TESTS_NONPERIODIC_TEMPL], "\t\tF R E Q U E N C Y\n");
fprintf(stats[TESTS_NONPERIODIC_TEMPL], "Template W_1 W_2 W_3 W_4 W_5 W_6 W_7 W_8");
fprintf(stats[TESTS_NONPERIODIC_TEMPL], " Chi^2 P_value Assignment Index\n");
fprintf(stats[TESTS_NONPERIODIC_TEMPL], "---------------------------------------------");
fprintf(stats[TESTS_NONPERIODIC_TEMPL], "-----------------------------------\n");
}*/
if ( numOfTemplates[m] < MAXNUMOFTEMPLATES )
SKIP = 1;
else
SKIP = (int)(numOfTemplates[m]/MAXNUMOFTEMPLATES);
numOfTemplates[m] = (int)numOfTemplates[m]/SKIP;
sum = 0.0;
for( i=0; i<2; i++ ) { /* Compute Probabilities */
pi[i] = exp(-lambda+i*log(lambda)-lgam(i+1));
sum += pi[i];
}
pi[0] = sum;
for( i=2; i<=K; i++ ) { /* Compute Probabilities */
pi[i-1] = exp(-lambda+i*log(lambda)-lgam(i+1));
sum += pi[i-1];
}
pi[K] = 1 - sum;
for( jj=0; jj<MIN(MAXNUMOFTEMPLATES,numOfTemplates[m]); jj++ ) {
sum = 0;
for(k = 0; k < m; k++) {
fscanf(fp, "%d", &bits);
templates[k].b = bits;
// printf("%d", templates[k].b);
}
//fprintf(stats[TESTS_NONPERIODIC_TEMPL]," ");
for( k=0; k<=K; k++ )
nu[k] = 0;
for( i=0; i<N; i++ ) {
W_obs = 0;
for( j=0; j<M-m+1; j++ ) {
match = 1;
for( k=0; k<m; k++ ) {
if ( (int)templates[k].b != (int)bit[i*M+j+k] ) {
match = 0;
break;
}
}
if ( match == 1 )
W_obs++;
}
Wj[i] = W_obs;
}
sum = 0;
chi2 = 0.0; /* Compute Chi Square */
for( i=0; i<N; i++ ) {
//if ( m == 10 )
// fprintf(stats[TESTS_NONPERIODIC_TEMPL], "%3d ", Wj[i]);
// else
// fprintf(stats[TESTS_NONPERIODIC_TEMPL], "%4d ", Wj[i]);
chi2 += pow(((double)Wj[i] - lambda)/pow(varWj,0.5),2);
}
p_value = igamc(N/2.0,chi2/2.0);
/*p_value = gammq(N/2.0,chi2/2.0);*/
if ( isNegative(p_value) || isGreaterThanOne(p_value) )
printf("\t\tWARNING: P_VALUE IS OUT OF RANGE.\n");
if ( p_value < ALPHA ) {
return false;
}
else {
return true;
}
//fprintf(stats[TESTS_NONPERIODIC_TEMPL], "%9.6f %f %s %3d\n", chi2, p_value, assignment, jj);
fseek(fp, (long)(SKIP-1)*2*m, SEEK_CUR);
// fprintf(results[TESTS_NONPERIODIC_TEMPL], "%f\n", p_value);
// fprintf(grid, "%d", state); fflush(grid);
// ind = jj / ((numOfTemplates[m]+NUMAPERIODICFILES-1)/NUMAPERIODICFILES);
// fprintf(pvalsAperiodic[ind], "%f ", p_value); fflush(pvalsAperiodic[ind]);
// if ( p_value < tp.minimumP )
// tp.minimumP = p_value;
// if ( !_isnan(p_value) )
// tp.lnSum += log(p_value);
// tp.df++;
}
// fprintf(stats[TESTS_NONPERIODIC_TEMPL], "\n");
// for ( ind=0; ind<NUMAPERIODICFILES; ind++ )
// fprintf(pvalsAperiodic[ind], "\n");
free(templates);
}
// fflush(stats[TESTS_NONPERIODIC_TEMPL]);
// fflush(results[TESTS_NONPERIODIC_TEMPL]);
// fflush(grid);
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -