/* C includes */
#include <stdio.h>
#include <stdlib.h>


/* local includes */
#include "stmc.h"


/* *** parameters *** */
/* seeds */
stmc_int seed1=1;
stmc_int seed2=0;
/* file to store the state of the random number object */
stmc_char ranmar_file[]="ranmar.d";

/* polls */
#define NCH 3 /* number of candidates */
#define NCP1 (NCH+1) /* number of candidates + undecided */
stmc_int Npoll=10000; /* number of poll repetitions */
stmc_int Ndat=700; /* number of poll entries */
stmc_real08 prob[NCH]={0.42,0.38,0.16}; /* initial poll probabilities */


int main(void)
{
  stmc_rannum *r; /* random number generator object */
  stmc_int i,idat,ipoll;
  stmc_real08 cp[NCP1],score[NCP1];
  stmc_real08 **result;
  stmc_real08 xr,bush,clinton;

  printf("Program POLLS:\n\n");
  
  printf("Probabilities:\n\n");
  for(i=0;i<NCH;i++) printf("%10.3f",prob[i]);
  printf("\n\n");
  
  printf("Cumulative probabilities:\n\n");
  cp[0]=prob[0];
  for(i=1;i<NCH;i++) cp[i]=cp[i-1]+prob[i];
  cp[NCP1-1]=1.0;
  for(i=0;i<NCP1;i++) printf("%10.3f",cp[i]);
  printf("\n\n");
  
  /* allocate 2D array for results */
  result=(stmc_real08**)malloc((size_t)NCP1*sizeof(stmc_real08*));
  for(i=0;i<NCP1;i++)
    result[i]=(stmc_real08*)malloc((size_t)Npoll*sizeof(stmc_real08));
  
  /* initialize random number object */
  r=stmc_rma_set(seed1, seed2, ranmar_file);

  bush=0;
  clinton=0;
  
  for(ipoll=0;ipoll<Npoll;ipoll++)
    {
      /* zero out array with scores */
      for(i=0;i<NCP1;i++) score[i]=0.0;
      
      /* simulate number of people asked */
      for(idat=0;idat<Ndat;idat++)
        {
          xr=stmc_ranmar(r);
          i=0;
          
          while((i<NCH)&&(xr>=cp[i]))
            {
              i++;
            }
          score[i]++;
        }
      
      if(score[0]>score[1]) bush++;
      if(score[0]<score[1]) clinton++;

      /* normalize results */
      for(i=0;i<NCP1;i++)
        {
          result[i][ipoll]=score[i]/Ndat;
        }    
    }
    
  bush=bush/Npoll;
  clinton=clinton/Npoll;
  
  printf("\nBush    = %f\n",bush);
  printf("Clinton = %f\n\n",clinton);
  
  /* sort and build peaked probability distribution
     for each candidate */
  for(i=0;i<NCP1;i++)
    {
      stmc_heapsort(Npoll,result[i]);
      stmc_qdf_gnu(Npoll,result[i]);
    }
  
  /* free random number object */
  stmc_rma_free(r);
  
  /* free results array */
  for(i=0;i<NCP1;i++)
    free(result[i]);
  free(result);
  
  return 0;
}


/* include all used routines */
#ifndef STMC_C_COMPILED_AS_LIBRARY
#include "stmc_rma_set.c"
#include "stmc_ranmar.c"
#include "stmc_rma_free.c"
#include "stmc_heapsort.c"
#include "stmc_qdf_gnu.c"
#include "stmc_min.c"
#endif