Files
rd_sqngen/sqngen.c
T

789 lines
28 KiB
C

//to compile use the command g++ -Wall sqngen.c -lgsl -lgslcblas -o sqngen
//syntax:
// sqngen num_of_samples length gmin_max loss_4state p13 p31 [p32 p23 [p14]]
//OR sqngen num_of_samples length gmin_max loss_GI [ploss [burst_length [density [pisol [good_burst_length]]]]
//OR sqngen num_of_samples length gmin_max loss_bern_4s p
//OR sqngen num_of_samples length gmin_max loss_gilb_4s p r [1-h]
//OR sqngen num_of_samples length gmin_max loss_gilb_ell_4s p r [1-h [1-k]]
//OR sqngen num_of_samples length gmin_max loss_gilb p r [1-h]
//OR sqngen num_of_samples length gmin_max loss_gilb_ell p r [1-h [1-k]]
//p13 p31 p32 p23 p14 ploss density pisol p r 1-k 1-h are percentage
//num_of_samples length gmin_max burst_length good_burst_length are integer
//
//This program is free software; you can redistribute it and/or
//modify it under the terms of the GNU General Public License
//s published by the Free Software Foundation; either version
//2 of the License, or (at your option) any later version.
//
//Author: Fabio Ludovici <fabio.ludovici at yahoo.it>
//Version 2.4
//January 2010
//
//This program was used to test the loss generation algorithms implemented in Netem2
//available at http://netgroup.uniroma2.it/twiki/bin/view.cgi/Main/NetEm2
#include <unistd.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <time.h>
#include <math.h>
#include <stdbool.h>
#include <gsl/gsl_math.h>
#include <gsl/gsl_rng.h>
#include <gsl/gsl_randist.h>
#ifdef WIN32
#include <windows.h>
#define sleep(a) Sleep(a * 1000)
#endif
#define ONE 4294967295u
#include <iostream>
using std::fixed;
#include <fstream>
using namespace std;
double p31;
double p13;
double p23=1;
double p32=0;
double p14=0;
double p;
double r;
double h;
double k;
double burst_length;
double good_burst_length;
double rho;
double rho_min;
double pisol;
double ploss;
void arg_error() {
printf("Error: too few arguments\n");
printf("Syntax is: sqngen num_of_samples length gmin_max loss_4state p13 p31 [p32 p23 p14]\n");
printf("or sqngen num_of_samples length gmin_max loss_GI ploss [burst_length [density [pisol [good_burst_length]]]]\n");
printf("or sqngen num_of_samples length gmin_max loss_bern_4s p\n");
printf("or sqngen num_of_samples length gmin_max loss_gilb_4s p r [1-h]\n");
printf("or sqngen num_of_samples length gmin_max loss_gilb_ell_4s p r [1-h [1-k]]\n");
printf("or sqngen num_of_samples length gmin_max loss_gilb p r [1-h]\n");
printf("or sqngen num_of_samples length gmin_max loss_gilb_ell p r [1-h [1-k]]\n");
printf("p13 p31 p32 p23 p14 ploss density pisol p r 1-k 1-h are percentage\n");
printf("num_of_samples length gmin_man burst_length good_burst_length are integer\n");
}
void sqn_gen(gsl_rng * q, FILE *out_file, char* sequence, int length, double p13, double p31, double p32, double p23, double p14) {
double GI_random;
int state=1;
int tags=0;
//if(length>0) printf("The sequence is:\n");
gsl_rng_set(q, time(NULL));
for (int i=0;i<(length+tags);i++) {
GI_random=gsl_rng_get(q);
GI_random=GI_random/ONE;
switch (state) {
case 1:
sequence[i]=0;
//printf("0");
fprintf(out_file, "1\n");
if((0<GI_random)&&(GI_random<p13)){
state=3;
sequence[i+1]='S';
//printf("S");
tags++;
i++;
break;
} else if((p13<GI_random)&&(GI_random<(1-p14))){
state=1;
break;
} else if(((1-p14)<GI_random) && (GI_random<1)){
state=4;
break;
}
case 2:
sequence[i]=0;
//printf("0");
fprintf(out_file, "1\n");
if((0<GI_random)&&(GI_random<p23)){
state=3;
break;
} else {
state=2;
break;
}
case 3:
sequence[i]=1;
//printf("1");
fprintf(out_file, "0\n");
if((0<GI_random)&&(GI_random<p32)){
state=2;
break;
} else if((p32<GI_random)&&(GI_random<(p31+p32))){
state=1;
sequence[i+1]='E';
//printf("E");
i++;
tags++;
break;
} else if((((p31+p32))<GI_random) && (GI_random<1)){
state=3;
break;
}
case 4:
sequence[i]=1;
//printf("1");
fprintf(out_file, "0\n");
state=1;
break;
}
}
}
void gilb_ell_gen(gsl_rng * q, char* sequence, int length, double p, double r, double k, double h) {
double gilb_ell_random;
int state=1;
int tags=0;
if(length>0) printf("The sequence is:\n");
gsl_rng_set(q, time(NULL));
for (int i=0;i<(length+tags);i++) {
switch (state) {
case 1:
gilb_ell_random=gsl_rng_get(q);
gilb_ell_random=gilb_ell_random/ONE;
if(gilb_ell_random<(1-k)) {
sequence[i]=1;
printf("1");
} else {
sequence[i]=0;
printf("0");
}
gilb_ell_random=gsl_rng_get(q);
gilb_ell_random=gilb_ell_random/ONE;
if(gilb_ell_random<p) {
sequence[i+1]='S';
printf("S");
i++;
tags++;
state=2;
}
break;
case 2:
gilb_ell_random=gsl_rng_get(q);
gilb_ell_random=gilb_ell_random/ONE;
if(gilb_ell_random<(1-h)) {
sequence[i]=1;
printf("1");
} else {
sequence[i]=0;
printf("0");
}
gilb_ell_random=gsl_rng_get(q);
gilb_ell_random=gilb_ell_random/ONE;
if(gilb_ell_random<r) {
sequence[i+1]='E';
printf("E");
i++;
tags++;
state=1;
}
break;
}
};
}
void internal_check(int length, char sequence[], int sample, double *ploss_bis, double *mean_length_of_burst_bis, double *density_bis, double *pisolated_loss_bis, double *mean_length_of_good_burst_bis) {
bool gilb_ell_burst=false;
bool real_burst=false;
bool good_burst=false;
int num_of_bursts_bis=0;
int length_of_all_bursts_bis=0;
int num_of_good_bursts_bis=0;
int length_of_all_good_bursts_bis=0;
int num_of_drops_bis=0;
int num_of_burst_drops_bis=0;
int num_of_isolated_drops_bis=0;
int initial_zeros=0;
int final_zeros=0;
int ones=0;
int tags=0;
for (int k=0;k<(length+tags);k++) {
if (sequence[k]=='S') {
gilb_ell_burst=true;
num_of_bursts_bis++;
tags++;
}
if (sequence[k]=='E') {
gilb_ell_burst=false;
real_burst=false;
good_burst=false;
tags++;
if (ones==0) {
num_of_bursts_bis--;
}
length_of_all_bursts_bis=length_of_all_bursts_bis-initial_zeros-final_zeros;
length_of_all_good_bursts_bis=length_of_all_good_bursts_bis-initial_zeros-final_zeros;
ones=0;
initial_zeros=0;
final_zeros=0;
}
if ((gilb_ell_burst==true) & (ones==0) & (sequence[k]==1)) {
real_burst=true;
}
if ((gilb_ell_burst==true) & ((sequence[k]==1) || (sequence[k]==0))) {
length_of_all_bursts_bis++;
if (sequence[k]==0) {
if (real_burst==false)
initial_zeros++;
else
final_zeros++;
if (real_burst==true) {
length_of_all_good_bursts_bis++;
if (good_burst==false) {
good_burst=true;
num_of_good_bursts_bis++;
}
}
}
if (sequence[k]==1) {
if (real_burst==false) real_burst=true;
if (good_burst==true) good_burst=false;
ones++;
final_zeros=0;
num_of_drops_bis++;
num_of_burst_drops_bis++;
}
}
if ((gilb_ell_burst==false) & (sequence[k]==1)) {
num_of_drops_bis++;
num_of_isolated_drops_bis++;
}
}
printf("\nThe generated sequence has the following characteristics (from the internal check):\n");
ploss_bis[sample]=100*((double)num_of_drops_bis)/((double)length);
pisolated_loss_bis[sample]=100*((double)num_of_isolated_drops_bis)/((double)length);
if (num_of_bursts_bis>0) mean_length_of_burst_bis[sample]=(double) length_of_all_bursts_bis/num_of_bursts_bis;
else mean_length_of_burst_bis[sample]=0;
if (num_of_good_bursts_bis>0) mean_length_of_good_burst_bis[sample]=(double) length_of_all_good_bursts_bis/num_of_good_bursts_bis;
else mean_length_of_good_burst_bis[sample]=0;
density_bis[sample]=100*(double) num_of_burst_drops_bis/length_of_all_bursts_bis;
printf("P_loss %.3f%% Mean_length_of_burst %.3f Density %.3f%% P_isol %.3f%% Mean_length_of_good_burst %.3f\n", ploss_bis[sample], mean_length_of_burst_bis[sample], density_bis[sample], pisolated_loss_bis[sample], mean_length_of_good_burst_bis[sample]);
fstream data;
data.open("data.dat", ios::out | ios::app);
data<<fixed<<ploss_bis[sample]<<";"<<mean_length_of_burst_bis[sample]<<";"<<density_bis[sample]<<";"<<pisolated_loss_bis[sample]<<";"<<mean_length_of_good_burst_bis[sample]<<";";
data.close();
}
void external_blind_check(int length, int gmin_max, char sequence[], int sample, double *ploss_ter, double *mean_length_of_burst_ter, double *density_ter, double *pisolated_loss_ter, double *mean_length_of_good_burst_ter) {
bool this_is_burst=false;
int num_of_bursts_ter=0;
int length_of_all_bursts_ter=0;
int num_of_good_bursts_ter=0;
int length_of_all_good_bursts_ter=0;
int num_of_drops_ter=0;
int num_of_burst_drops_ter=0;
int num_of_isolated_drops_ter=0;
int tags=0;
for (int k=0;k<(length+tags);k++) {
if((sequence[k]=='S') || (sequence[k]=='E')) {
sequence[k]=sequence[k+1];
tags++;
}
}
for(int gmin=2;gmin<gmin_max*2;gmin=gmin*2) {
this_is_burst=false;
num_of_bursts_ter=0;
num_of_good_bursts_ter=0;
length_of_all_bursts_ter=0;
length_of_all_good_bursts_ter=0;
num_of_drops_ter=0;
num_of_burst_drops_ter=0;
num_of_isolated_drops_ter=0;
for (int k=0;k<length;k++) {
if (this_is_burst==true) { /* BURST ON this_is_burst=true */
if ((sequence[k-1]==1) & (sequence[k]==0)) num_of_good_bursts_ter++;
if (sequence[k]==0) {
length_of_all_bursts_ter++;
length_of_all_good_bursts_ter++;
} else if (sequence[k]==1) {
num_of_drops_ter++;
num_of_burst_drops_ter++;
length_of_all_bursts_ter++;
}
if ((k<(length-gmin)) & (sequence[k]==1) & (sequence[k+1]==0)) {
int zeros=1;
for(int l=1;l<gmin;l++) {
if (sequence[k+1+l]==0) zeros++;
}
if (zeros==gmin) this_is_burst=false;
}
} else { /* BURST OFF this_is_burst=false */
if (k==(length-1)) {
if(sequence[k]==1) {
num_of_drops_ter++;
num_of_isolated_drops_ter++;
}
} else if ((sequence[k]==1) & (sequence[k+1]==1)) {
this_is_burst=true;
num_of_drops_ter++;
num_of_bursts_ter++;
num_of_burst_drops_ter++;
length_of_all_bursts_ter++;
} else if ((sequence[k]==1) & (sequence[k+1]==0)) {
num_of_drops_ter++;
num_of_isolated_drops_ter++;
}
}
}
printf("\nFor gmin=%d the generated sequence has the following characteristics (from sequence analysis):\n",gmin);
int i=log2(gmin/2)+1;
int j=log2(gmin_max/2)+1;
ploss_ter[sample*j+i-1]=100*((double)num_of_drops_ter)/((double)length);
pisolated_loss_ter[sample*j+i-1]=100*((double)num_of_isolated_drops_ter)/((double)length);
if (num_of_bursts_ter>0) mean_length_of_burst_ter[sample*j+i-1]=(double) length_of_all_bursts_ter/num_of_bursts_ter;
else mean_length_of_burst_ter[sample*j+i-1]=0;
if (num_of_good_bursts_ter>0) mean_length_of_good_burst_ter[sample*j+i-1]=(double) length_of_all_good_bursts_ter/num_of_good_bursts_ter;
else mean_length_of_good_burst_ter[sample*j+i-1]=0;
density_ter[sample*j+i-1]=100*(((double)num_of_burst_drops_ter)/((double)length_of_all_bursts_ter));
printf("P_loss %.3f%% Mean_length_of_burst %.3f Density %.3f%% P_isol %.3f%% Mean_length_of_good_burst %.3f\n", ploss_ter[sample*j+i-1], mean_length_of_burst_ter[sample*j+i-1], density_ter[sample*j+i-1], pisolated_loss_ter[sample*j+i-1], mean_length_of_good_burst_ter[sample*j+i-1]);
fstream data;
data.open("data.dat", ios::out | ios::app);
data<<gmin<<";"<<fixed<<ploss_ter[sample*j+i-1]<<";"<<mean_length_of_burst_ter[sample*j+i-1]<<";"<<density_ter[sample*j+i-1]<<";"<<pisolated_loss_ter[sample*j+i-1]<<";"<<mean_length_of_good_burst_ter[sample*j+i-1]<<";";
data.close();
}
}
void statistics(char *filename, int num_of_samples, double ploss[], double mean_length_of_burst[], double density[], double pisolated_loss[], double mean_length_of_good_burst[], int gmin, int gmin_max) {
double ploss_mean=0;
double mean_length_of_burst_mean=0;
double density_mean=0;
double pisolated_loss_mean=0;
double mean_length_of_good_burst_mean=0;
double ploss_variance=0;
double mean_length_of_burst_variance=0;
double density_variance=0;
double pisolated_loss_variance=0;
double mean_length_of_good_burst_variance=0;
double ploss_stddev=0;
double mean_length_of_burst_stddev=0;
double density_stddev=0;
double pisolated_loss_stddev=0;
double mean_length_of_good_burst_stddev=0;
double ploss_cx=0;
double mean_length_of_burst_cx=0;
double density_cx=0;
double pisolated_loss_cx=0;
double mean_length_of_good_burst_cx=0;
int i=log2(gmin/2)+1;
int j=log2(gmin_max/2)+1;
for(int sample=0;sample<num_of_samples;sample++) {
ploss_mean=ploss_mean+ploss[sample*j+i-1];
mean_length_of_burst_mean=mean_length_of_burst_mean+mean_length_of_burst[sample*j+i-1];
density_mean=density_mean+density[sample*j+i-1];
pisolated_loss_mean=pisolated_loss_mean+pisolated_loss[sample*j+i-1];
mean_length_of_good_burst_mean=mean_length_of_good_burst_mean+mean_length_of_good_burst[sample*j+i-1];
}
printf("\nThese are the mean values:\n");
ploss_mean=ploss_mean/num_of_samples;
mean_length_of_burst_mean=mean_length_of_burst_mean/num_of_samples;
density_mean=density_mean/num_of_samples;
pisolated_loss_mean=pisolated_loss_mean/num_of_samples;
mean_length_of_good_burst_mean=mean_length_of_good_burst_mean/num_of_samples;
printf("ploss %.3f%% mean_length_of_burst %.3f density %.3f%% pisol %.3f%% mean_length_of_good_burst %.3f\n", ploss_mean, mean_length_of_burst_mean, density_mean, pisolated_loss_mean, mean_length_of_good_burst_mean);
for(int sample=0;sample<num_of_samples;sample++) {
ploss_variance=ploss_variance+(ploss[sample*j+i-1]-ploss_mean)*(ploss[sample*j+i-1]-ploss_mean);
mean_length_of_burst_variance=mean_length_of_burst_variance+(mean_length_of_burst[sample*j+i-1]-mean_length_of_burst_mean)*(mean_length_of_burst[sample*j+i-1]-mean_length_of_burst_mean);
density_variance=density_variance+(density[sample*j+i-1]-density_mean)*(density[sample*j+i-1]-density_mean);
pisolated_loss_variance=pisolated_loss_variance+(pisolated_loss[sample*j+i-1]-pisolated_loss_mean)*(pisolated_loss[sample*j+i-1]-pisolated_loss_mean);
mean_length_of_good_burst_variance=mean_length_of_good_burst_variance+(mean_length_of_good_burst[sample*j+i-1]-mean_length_of_good_burst_mean)*(mean_length_of_good_burst[sample*j+i-1]-mean_length_of_good_burst_mean);
}
printf("\nThese are the variance values:\n");
ploss_variance=ploss_variance/(num_of_samples-1);
mean_length_of_burst_variance=mean_length_of_burst_variance/(num_of_samples-1);
density_variance=density_variance/(num_of_samples-1);
pisolated_loss_variance=pisolated_loss_variance/(num_of_samples-1);
mean_length_of_good_burst_variance=mean_length_of_good_burst_variance/(num_of_samples-1);
printf("ploss %.3f mean_length_of_burst %.3f density %.3f pisol %.3f mean_length_of_good_burst %.3f\n", ploss_variance, mean_length_of_burst_variance, density_variance, pisolated_loss_variance, mean_length_of_good_burst_variance);
printf("\nThese are the standard deviation values:\n");
ploss_stddev=sqrt(ploss_variance);
mean_length_of_burst_stddev=sqrt(mean_length_of_burst_variance);
density_stddev=sqrt(density_variance);
pisolated_loss_stddev=sqrt(pisolated_loss_variance);
mean_length_of_good_burst_stddev=sqrt(mean_length_of_good_burst_variance);
printf("ploss %.3f%% mean_length_of_burst %.3f density %.3f%% pisol %.3f%% mean_length_of_good_burst %.3f\n", ploss_stddev, mean_length_of_burst_stddev, density_stddev, pisolated_loss_stddev, mean_length_of_good_burst_stddev);
printf("\nThese are the variation coefficients:\n");
ploss_cx=ploss_stddev/ploss_mean;
mean_length_of_burst_cx=mean_length_of_burst_stddev/mean_length_of_burst_mean;
density_cx=density_stddev/density_mean;
if(pisolated_loss_mean>0)
pisolated_loss_cx=pisolated_loss_stddev/pisolated_loss_mean;
else
pisolated_loss_cx=0;
mean_length_of_good_burst_cx=mean_length_of_good_burst_stddev/mean_length_of_good_burst_mean;
printf("ploss %.3f mean_length_on_burst %.3f density %.3f pisol %.3f mean_length_of_good_burst %.3f\n", ploss_cx, mean_length_of_burst_cx, density_cx, pisolated_loss_cx, mean_length_of_good_burst_cx);
fstream data;
data.open(filename, ios::out | ios::app);
data <<"Average Prob. Loss; Average Burst Length; Average Density; Average Prob. Loss Isol.; Average Good Burst Length;\n"; data<<fixed<<ploss_mean<<";"<<mean_length_of_burst_mean<<";"<<density_mean<<";"<<pisolated_loss_mean<<";"<<mean_length_of_good_burst_mean<<";\n";
data <<"Variance Prob. Loss; Variance Burst Length; Variance Density; Variance Prob. Loss Isol.; Variance Good Burst length;\n"; data<<ploss_variance<<";"<<mean_length_of_burst_variance<<";"<<density_variance<<";"<<pisolated_loss_variance<<";"<<mean_length_of_good_burst_variance<<";\n";
data <<"Standard Dev. Prob. Loss; Standard Dev. Burst Length; Standard Dev. Density; Standard Dev. Prob. Loss Isol.; Standard Dev. Good Burst length;\n";
data<<ploss_stddev<<";"<<mean_length_of_burst_stddev<<";"<<density_stddev<<";"<<pisolated_loss_stddev<<";"<<mean_length_of_good_burst_stddev<<";\n";
data <<"Variation Prob. Loss; Variation Burst Length; Variation Density; Variation Prob. Loss Isol.; Variation Good Burst length;\n";
data<<ploss_cx<<";"<<mean_length_of_burst_cx<<";"<<density_cx<<";"<<pisolated_loss_cx<<";"<<mean_length_of_good_burst_cx<<";\n";
data.close();
}
int main (int argc, char *argv[]) {
int num_of_samples;
int length;
int gmin_max;
fstream data;
char *filename_common;
fprintf(stdout,"********************** sqngen ***********************\n");
if (argc<7) {
arg_error();
return 0;
}
for (int i=0;i<argc;i++)
{
if (strcmp (argv[i],"-o") == 0)
{
filename_common = argv[i+1];
argc = argc - 2;
break;
}
}
if (argv[1]) num_of_samples=atoi(argv[1]);
if (argv[2]) length=atoi(argv[2]);
gmin_max=atoi(argv[3]);
if (strcmp (argv[4],"loss_4state") == 0) {
if (argc<7) {
arg_error();
return 0;
}
if (argc>7 && (atof(argv[6])+atof(argv[7])>100)) {printf("\nError: p31+p32>100\n"); return 0;}
if (argc>9 && (atof(argv[5])+atof(argv[9])>100)) {printf("\nError: p13+p14>100\n"); return 0;}
p13=atof(argv[5])/100;
p31=atof(argv[6])/100;
if (argc>7)
p32=atof(argv[7])/100;
else
p32=0;
if (argc>8)
p23=atof(argv[8])/100;
else
p23=1;
if (argc>9)
p14=atof(argv[9])/100;
else
p14=0;
} else if (strcmp (argv[4],"loss_GI") == 0) {
if (argc<6) {
arg_error();
return 0;
}
/* Change here to get a 0.1% resolution */
//ploss=(double) atof(argv[5])/100;
ploss=(double) atof(argv[5])/1000;
if (ploss==0) {
printf("\nError: ploss=0\n");
return 0;
}
if (argc>8 && (atof(argv[8])>atof(argv[5]))) {
printf("\nError: pisol>ploss\n");
return 0;
}
if (argc>9 && (atof(argv[9])>atof(argv[6]))) {
printf("\nError: good_burst_length>burst_length\n");
return 0;
}
if (argc>6)
burst_length=atof(argv[6]);
else
burst_length=1/(1-ploss);
if (argc>7)
rho=atof(argv[7])/100;
else
rho=1;
if (argc>8)
pisol=atof(argv[8])/100;
if (argc>9)
good_burst_length=atoi(argv[9]);
else
good_burst_length=0;
p13=(pisol-ploss)/(burst_length*(pisol*(rho-1)-rho+ploss));
p31=1/(burst_length*rho);
if (rho<1) {
if (good_burst_length>0) {
p23=(1/good_burst_length);
p32=(1-rho)/(rho*good_burst_length);
}
else {
p32=(1+rho*rho*burst_length-rho-burst_length*rho)/(rho-rho*burst_length);
p23=(burst_length*rho-1)/(burst_length-1);
}
} else {
p32=0;
p23=1;
}
if (pisol>0)
p14=(double) pisol/(1-pisol);
else
p14=0;
} else if (strcmp (argv[4],"loss_gilb_4s") == 0) {
if (argc<7) {
arg_error();
return 0;
}
p=(atof(argv[5]))/100;
r=(atof(argv[6]))/100;
if(argv[7])
h=1-((atof(argv[7]))/100);
else
h=0;
ploss=(p-h*p)/(p+r);
burst_length=(1-h+r*r*h)/(r*(1-h+r*h));
if(argv[7])
rho=((1-h+r*h)*(1-h+r*h)) / (1-h+r*r*h);
else
rho=1;
p13=(-ploss)/(burst_length*(-rho+ploss));
p31=1/(burst_length*rho);
if (rho<1) {
p32=(1+rho*rho*burst_length-rho-burst_length*rho)/(rho-rho*burst_length);
p23=(burst_length*rho-1)/(burst_length-1);
} else {
p32=0;
p23=1;
}
p14=0;
} else if (strcmp (argv[4],"loss_gilb_ell_4s") == 0) {
if (argc<7) {
arg_error();
return 0;
}
p=(atof(argv[5]))/100;
r=(atof(argv[6]))/100;
if(argv[7])
h=1-((atof(argv[7]))/100);
else
h=0;
if(argv[8])
k=1-((atof(argv[8]))/100);
else
k=1;
ploss=(r-k*r+p-h*p)/(p+r);
burst_length=(1-h+r*r*h)/(r*(1-h+r*h));
if(argv[7])
rho=((1-h+r*h)*(1-h+r*h)) / (1-h+r*r*h);
else
rho=1;
if(argv[8])
pisol=(r-k*r)/(p+r);
else
pisol=0;
p13=(pisol-ploss)/(burst_length*(pisol*(rho-1)-rho+ploss));
p31=1/(burst_length*rho);
if (argv[7]) {
p32=(1+rho*rho*burst_length-rho-burst_length*rho) / (rho-rho*burst_length);
p23=(burst_length*rho-1) / (burst_length-1);
} else
p32=0;
if (pisol>0)
p14=pisol/(1-pisol);
else
p14=0;
} else if (strcmp (argv[4],"loss_bern") == 0) {
p=atof(argv[5])/100;
r=1-p;
h=0;
k=1;
} else if (strcmp (argv[4],"loss_gilb") == 0) {
if (argc<7) {
arg_error();
return 0;
}
p=atof(argv[5])/100;
r=atof(argv[6])/100;
if(argc>7) {
h=1-((atof(argv[7]))/100);
k=0;
} else {
h=0;
k=1;
}
} else if (strcmp (argv[4],"loss_gilb_ell") == 0) {
if (argc<7) {
arg_error();
return 0;
}
p=atof(argv[5])/100;
r=atof(argv[6])/100;
if(argc>7){
h=1-((atof(argv[7]))/100);
k=0;
} else {
h=0;
k=1;
}
if(argc>8) {
h=1-((atof(argv[7]))/100);
k=1-((atof(argv[8]))/100);
}
else {
k=1;
}
}
char sequence[2*length];
double ploss_bis[num_of_samples];
double pisolated_loss_bis[num_of_samples];
double mean_length_of_burst_bis[num_of_samples];
double density_bis[num_of_samples];
double mean_length_of_good_burst_bis[num_of_samples];// good burst length
double ploss_ter[num_of_samples*gmin_max];
double pisolated_loss_ter[num_of_samples*gmin_max];
double mean_length_of_burst_ter[num_of_samples*gmin_max];
double density_ter[num_of_samples*gmin_max];
double mean_length_of_good_burst_ter[num_of_samples]; // good burst length
gsl_rng * q=gsl_rng_alloc (gsl_rng_taus2);
data.open("data.dat", ios::out);
data << "Sample; Num of Samples; Legth; Model; PLoss; BurstLegth";
data.close();
for (int sample=0;sample<num_of_samples;sample++) {
FILE *out_file = NULL;
char filename[200];
sprintf(filename,"%s_P%d_%d.dat", filename_common, (int)(ploss*1000), sample);
out_file = fopen (filename, "w");
if (out_file == NULL)
{
printf("Error opening the error pattern file\n");
return 0;
}
//fprintf(out_file, "nPackets: %d\nVector:\n",length);
data.open("data.dat", ios::out | ios::app);
data <<"\n"<<(sample+1)<<";";
data <<num_of_samples<<";";
data <<length<<";";
data <<gmin_max<<";";
data <<argv[4]<<";"; //chosen model
data <<atof(argv[5])<<";"; //saving input values
if(argv[6]) data <<atof(argv[6])<<";"; else data <<"N/A"<<";";
if((argv[6]) && (argv[7])) data <<atof(argv[7])<<";"; else data <<"100"<<";";
if((argv[6]) && (argv[7]) && (argv[8])) data <<atof(argv[8])<<";"; else data <<"0"<<";";
if((argv[6]) && (argv[7]) && (argv[8]) && (argv[9])) data <<atof(argv[9])<<";"; else data <<"N/A"<<";";
data.close();
sleep (1);
if ((strcmp (argv[4],"loss_bern") == 0) || (strcmp (argv[4],"loss_gilb") == 0) || (strcmp (argv[4],"loss_gilb_ell") == 0))
gilb_ell_gen(q, sequence, length, p, r, k, h);
else
sqn_gen(q, out_file, sequence, length, p13, p31, p32, p23, p14);
fclose (out_file);
printf ("The selected generation option is %s with the following input parameters:\n", argv[4]);
if (strcmp (argv[4],"loss_4state") == 0)
printf("p13=%.3f%% \np31=%.3f%% \np32=%.3f%% \np23=%.3f%% \np14=%.3f%%\n", (double)100*p13/1, (double)100*p31/1, (double)100*p32/1, (double)100*p23/1, (double)100*p14/1);
else if ((strcmp (argv[4],"loss_GI") == 0) & (good_burst_length==0))
printf("P_loss %.3f%% Mean_length_of_burst %.3f%% Density %.3f%% P_isol %.3f%% Mean_length_of_good_burst N/A\n", 100*ploss, burst_length, 100*rho, 100*pisol);
else if ((strcmp (argv[4],"loss_GI") == 0) & (good_burst_length==0))
printf("P_loss %.3f%% Mean_length_of_burst %.3f%% Density %.3f%% P_isol %.3f%% Mean_length_of_good_burst N/A\n", 100*ploss, burst_length, 100*rho, 100*pisol);
else if ((strcmp (argv[4],"loss_GI") == 0) & (good_burst_length!=0))
printf("P_loss %.3f%% Mean_length_of_burst %.3f%% Density %.3f%% P_isol %.3f%% Mean_length_of_good_burst %.3f%%\n", 100*ploss, burst_length, 100*rho, 100*pisol, good_burst_length);
else if (strcmp (argv[4],"loss_gilb_4s") == 0)
printf("p=%.3f%% r=%.3f%% 1-h=%.3f%%\n", 100*p, 100*r, 100*(1-h)/1);
else if (strcmp (argv[4],"loss_gilb_ell_4s") == 0)
printf("p=%.3f%% r=%.3f%% 1-h=%.3f%% 1-k=%.3f%%\n", 100*p, 100*r, 100*(1-h), 100*(1-k));
else if (strcmp (argv[4],"loss_bern") == 0)
printf("p=%.3f%%\n", (double)100*p/1);
else if (strcmp (argv[4],"loss_gilb") == 0)
printf("p=%.3f%% r=%.3f%% 1-h=%.3f%%\n", (double)100*p/1, (double)100*r/1, (double) 100*(1-h)/1);
else if (strcmp (argv[4],"loss_gilb_ell") == 0)
printf("p=%.3f%% r=%.3f%% 1-h=%.3f%% 1-k=%.3f%%\n", (double)100*p, (double)100*r, (double) 100*(1-h), (double) 100*(1-k));
if ((strcmp (argv[4],"loss_GI") == 0) || (strcmp (argv[4],"loss_bern_4s") == 0) || (strcmp (argv[4],"loss_gilb_4s") == 0) || (strcmp (argv[4],"loss_gilb_ell_4s") == 0)) {
printf ("The correspondent transition probabilities are:\n");
printf("p13=%.3f%% \np31=%.3f%% \np32=%.3f%% \np23=%.3f%% \np14=%.3f%%\n", (double)100*p13/1, (double)100*p31/1, (double)100*p32/1, (double)100*p23/1, (double)100*p14/1);
}
if ((strcmp (argv[4],"loss_GI") != 0) & (strcmp (argv[4],"loss_gilb") != 0) & (strcmp (argv[4],"loss_gilb_ell") != 0)) {
printf ("The correspondent GI (General and Intuitive) parameters are:");
printf ("P_loss=%.3f%% Mean_length_of_burst=%.3f Density=%.3f%% P_isol=%.3f%% Mean_length_of_good_burst %.3f \n", (double) 100*(p13*p23+p14*p23*p31) / (p13*p23+p23*p31+p14*p23*p31+p13*p32), (double) (p32+p23) / (p23*p31), (double) 100*(p13*p23) / (p13*p23+p13*p32), (double) 100*(p14*p23*p31) / (p14*p23*p31+p23*p31), (double) 1/p23 );
}
if(length>0) {
internal_check(length, sequence, sample, ploss_bis, mean_length_of_burst_bis, density_bis, pisolated_loss_bis, mean_length_of_good_burst_bis );
if (gmin_max>0)
external_blind_check (length, gmin_max, sequence, sample, ploss_ter, mean_length_of_burst_ter, density_ter, pisolated_loss_ter, mean_length_of_good_burst_ter);
}
}
if (num_of_samples>1) {
char filename[100];
sprintf(filename,"%s_P%d_B%d_EPStatistics.csv", filename_common, atoi(argv[5]), atoi(argv[6]) );
printf("INTERNAL CHECK STATISTICS\n");
data.open(filename, ios::out);
data << "Internal Check Statistics\n";
data << "N Samples; Length of vector; Model;";
if (strcmp (argv[4],"loss_GI") == 0){
data << "PLoss; Burst Legth;\n";
}
data <<num_of_samples<<";";
data <<length<<";";
data <<argv[4]<<";"; //chosen model
data <<fixed<<atof(argv[5])<<";"; //saving internal check values
if(argv[6])
data <<atof(argv[6])<<";";
else
data <<"N/A"<<";";
if((argv[6]) && (argv[7]))
data <<atof(argv[7])<<";";
else
data <<"100"<<";";
if((argv[6]) && (argv[7]) && (argv[8]))
data <<atof(argv[8])<<";";
else
data <<"0"<<";";
if((argv[6]) && (argv[7]) && (argv[8]) && (argv[9]))
data <<atof(argv[9])<<";";
else
data <<"N/A"<<";";
data <<"\n";
data.close();
statistics(filename, num_of_samples, ploss_bis, mean_length_of_burst_bis, density_bis, pisolated_loss_bis, mean_length_of_good_burst_bis, 2, 2);
}
if ((gmin_max>0) && (num_of_samples>1)) {
char filename[100];
sprintf(filename,"%s_P%d_B%d_EPStatistics.csv", filename_common, atoi(argv[5]), atoi(argv[6]) );
if (gmin_max>0)
printf("\nEXTERNAL CHECK STATISTICS\n");
for (int gmin=2;gmin<gmin_max*2;gmin=gmin*2) {
printf("\nFor gmin=%d",gmin);
data.open("statistics.dat", ios::out | ios::app);
data <<"external"<<";"<<gmin<<";";
data <<num_of_samples<<";";
data <<length<<";";
data <<argv[4]<<";"; //chosen model
data <<fixed<<atof(argv[5])<<";";
if(argv[6])
data <<atof(argv[6])<<";";
else data <<"N/A"<<";";
if((argv[6]) && (argv[7]))
data <<atof(argv[7])<<";";
else
data <<"N/A"<<";";
if((argv[6]) && (argv[7]) && (argv[8]))
data <<atof(argv[8])<<";";
else
data <<"N/A"<<";";
if((argv[6]) && (argv[7]) && (argv[8]) && (argv[9]))
data <<atof(argv[9])<<";";
else
data <<"N/A"<<";";
data.close();
statistics(filename, num_of_samples, ploss_ter, mean_length_of_burst_ter, density_ter, pisolated_loss_ter, mean_length_of_good_burst_ter, gmin, gmin_max);
}
}
}