/*2kadai02-5.c*/

/*ran1の定義*/
#define IA 16807
#define IM 2147483647
#define AM (1.0/IM)
#define IQ 127773
#define IR 2836
#define NTAB 32
#define NDIV (1+(IM-1)/NTAB)
#define EPS 1.2e-7
#define RNMX (1.0-EPS)

double ran1(long *idum)
{
int j;
long k;
static long iy=0;
static long iv[NTAB];
double temp;

if (*idum <= 0 || !iy){
if (-(*idum) < 1) *idum=1;
else *idum = -(*idum);
for (j=NTAB+7;j>=0;j--) {
k=(*idum)/IQ;
*idum=IA*(*idum-k*IQ)-IR*k;
if (*idum <0) *idum += IM;
if (j <NTAB) iv[j] = *idum;
}
iy=iv[0];
}
k=(*idum)/IQ;
*idum=IA*(*idum-k*IQ)-IR*k;
if (*idum < 0) *idum += IM;
j=iy/NDIV;
iy=iv[j];
iv[j] = *idum;
if ((temp=AM*iy) > RNMX) return RNMX;
else return temp;
}


/*gaussの定義*/
#include<math.h>
double gasdev(long *idum)
{
double ran1(long *idum);
static int iset=0;
static double gset;
double fac,rsq,v1,v2;

if (iset ==0){
do{
v1=2.0*ran1(idum)-1.0;
v2=2.0*ran1(idum)-1.0;
rsq=v1*v1+v2*v2;
}while (rsq >= 1.0 || rsq == 0.0);
fac=sqrt(-2.0*log(rsq)/rsq);
gset=v1*fac;
iset=1;
return v2*fac;
}else{
iset=0;
return gset;
}
}

#define IB1 1
#define IB2 2
#define IB5 16
#define IB18 131072
#define MASK (IB1+IB2+IB5)

int irbit2(unsigned long *iseed)
{
if(*iseed & IB18){
*iseed=((*iseed ^ MASK)<<1) | IB1;
return 1;
}else{
*iseed <<= 1;
return 0;
}
}


#include<stdio.h>
#include<stdlib.h>
#include<math.h>
main(){
double SN,a,b,B,S,A,r,P,p,m,d,ctr,k,total;

unsigned long *iseed;
long *idum;
printf("SN比を入力してください:");
scanf("%lf",&SN);
a=SN/10;
printf("SN/10=%lf\n",a);
b=pow(10,a);
printf("真値=%lf\n",b);

B=sqrt(b);
printf("B=%lf\n",B);


total=1000.0;
A=1.0;
*iseed=189;
*idum=-1;
for(k=0;k<total;k++){
/*送信機*/
m=irbit2(iseed);
if(m==1){
S=A;
}else{
S=-A;
}
/*通信路*/
r=B*S+gasdev(idum);

/*受信機*/
if(r>0.0){
d=1;
}else{
d=0;}
/*カウント*/
if(m==d)
{
ctr=ctr;
}else{
ctr=ctr+1;
}

}
printf("誤った個数は%lfです。\n",ctr);
P=ctr/total; /*誤り率*/
printf("誤り率は%lfです。\n",P);
p=P*100.0;
printf("誤り率は%lf％です。\n",p);


}