#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <time.h>

#define EVENPC 0
#define ODDPC 1

typedef struct _node {
   double p;
   struct _node *L, *R;
} node;

typedef node *bintree;

double *genprobs ( int n )
{
   double *H;
   int i;

   if (n <= 0) return NULL;
   H = (double *)malloc(n * sizeof(double));
   if (n == 4) {
      H[0] = 0.1; H[1] = 0.2; H[2] = 0.3; H[3] = 0.4;
   } else {
      for (i=0; i<n; ++i) H[i] = (double)(100 + rand() % 801) / (double)1000;
   }
   printf("\n\n");
   return H;
}

bintree buildtree ( double H[], int n, int nh, int nt, double nodep, int parity )
{
   node *T;

   if (nh + nt > n) return NULL;
   T = (node *)malloc(sizeof(node));
   T -> p = nodep;
   if (  (nh + nt == n) ||
         ((parity == EVENPC) && (nh == nt) && (nh > 0)) ||
         ((parity == ODDPC) && (nh > nt)) ) {
      T -> L = T -> R = NULL;
   } else {
      T -> L = buildtree(H,n,nh+1,nt,nodep*H[nh+nt],parity);
      T -> R = buildtree(H,n,nh,nt+1,nodep*(1-H[nh+nt]),parity);
   }
   return T;
}

/* Helper function for allevents */
void allevrec ( bintree T, int nh, int nt, char E[], int i, int parity )
{
   if (T == NULL) return;
   if ( (T -> L == NULL) && (T -> R == NULL) ) {
      E[i] = '\0';
      printf("    Event: %-10s    Probability = %15.12lf    ", E, T -> p);
      if ( ((parity == EVENPC) && (nh == nt) && (nh > 0)) ||
           ((parity == ODDPC) && (nh > nt)) )
         printf("[Successful termination]\n");
      else
         printf("[Unsuccessful termination]\n");
   }
   E[i] = 'H'; allevrec(T->L,nh+1,nt,E,i+1,parity);
   E[i] = 'T'; allevrec(T->R,nh,nt+1,E,i+1,parity);
}

void allevents ( bintree T, int n, int parity )
{
   char *E;

   printf("\n+++ All termination events\n");
   E = (char *)malloc((n + 1) * sizeof(char));
   allevrec(T,0,0,E,0,parity);
   free(E);
}

/* Helper function for extremeevents */
void extevrec( bintree T, int nh, int nt, char E[], int i,
               char ESmin[], char ESmax[], char EUmin[], char EUmax[],
               double *spmin, double *spmax, double *upmin, double *upmax, int parity )
{
   if (T == NULL) return;
   if ( (T -> L == NULL) && (T -> R == NULL) ) {
      E[i] = '\0';
      if ( ((parity == EVENPC) && (nh == nt) && (nh > 0)) ||
           ((parity == ODDPC) && (nh > nt)) ) {
         if (T -> p < *spmin) { *spmin = T -> p; strcpy(ESmin, E); }
         if (T -> p > *spmax) { *spmax = T -> p; strcpy(ESmax, E); }
      } else {
         if (T -> p < *upmin) { *upmin = T -> p; strcpy(EUmin, E); }
         if (T -> p > *upmax) { *upmax = T -> p; strcpy(EUmax, E); }
      }
   }
   E[i] = 'H'; extevrec(T->L,nh+1,nt,E,i+1,ESmin,ESmax,EUmin,EUmax,spmin,spmax,upmin,upmax,parity);
   E[i] = 'T'; extevrec(T->R,nh,nt+1,E,i+1,ESmin,ESmax,EUmin,EUmax,spmin,spmax,upmin,upmax,parity);
}

void extremeevents ( bintree T, int n, int parity )
{
   char *E, *ESmin, *ESmax, *EUmin, *EUmax;
   double spmin, spmax, upmin, upmax;

   printf("\n+++ Extreme termination events\n");
   E = (char *)malloc((n + 1) * sizeof(char));
   ESmin = (char *)malloc((n + 1) * sizeof(char));
   ESmax = (char *)malloc((n + 1) * sizeof(char));
   EUmin = (char *)malloc((n + 1) * sizeof(char));
   EUmax = (char *)malloc((n + 1) * sizeof(char));
   spmin = upmin = 1; spmax = upmax = 0;
   extevrec(T,0,0,E,0,ESmin,ESmax,EUmin,EUmax,&spmin,&spmax,&upmin,&upmax,parity);
   printf("    Most likely successful termination event     : %s\n", ESmax);
   printf("    Most unlikely successful termination event   : %s\n", ESmin);
   printf("    Most likely unsuccessful termination event   : %s\n", EUmax);
   printf("    Most unlikely unsuccessful termination event : %s\n", EUmin);
   free(E); free(ESmin); free(ESmax); free(EUmin); free(EUmax);
}

/* Helper function for successprob */
void calcprobrec ( bintree T, int nh, int nt, double *succp, double *failp, int parity )
{
   if (T == NULL) return;
   if ( (T -> L == NULL) && (T -> R == NULL) ) {
      if ( ((parity == EVENPC) && (nh == nt) && (nh > 0)) ||
           ((parity == ODDPC) && (nh > nt)) )
         *succp += T -> p;
      else
         *failp += T -> p;
   }
   calcprobrec(T->L,nh+1,nt,succp,failp,parity);
   calcprobrec(T->R,nh,nt+1,succp,failp,parity);
}

void successprob ( bintree T, int parity )
{
   double succp, failp;

   succp = failp = 0;
   calcprobrec(T,0,0,&succp,&failp,parity);
   printf("\n+++ Total probabilities\n");
   printf("    Probability of successful termination   : %15.12lf\n", succp);
   printf("    Probability of unsuccessful termination : %15.12lf\n", failp);
}

void equalhead ( double H[], int n )
{
   double *SP, *P, *SN, *N, succp, failp;
   int i, t;

   /* The S arrays store probabilities of surplus of heads over tails so far.
      Only the surplus count matters, not how it is achieved.
      SP - previous iteration. SN - next iteration.
      These arrays are logically indexed by P and N in the range [-n,n]. */

   SP = (double *)malloc((2*n+1) * sizeof(double));
   P = SP + n;
   SN = (double *)malloc((2*n+1) * sizeof(double));
   N = SN + n;
   for (i=-n; i<=n; ++i) N[i] = 0;
   N[0] = 1; succp = 0;
   for (t=0; t<n; ++t) {
      for (i=-t; i<=t; ++i) { P[i] = N[i]; N[i] = 0; }
      for (i=-t; i<=t; ++i) {
         N[i-1] += P[i] * (1 - H[t]);
         N[i+1] += P[i] * H[t];
      }
      succp += N[0];
      N[0] = 0;
   }
   /* failp = 1 - succp. But I compute failp more stylistically in O(n) time.
      We have already spent O(n^2) time anyway. */
   failp = 0; for (i=-n; i<=n; ++i) failp += N[i];
   printf("    Probability of successful termination   : %15.12lf\n", succp);
   printf("    Probability of unsuccessful termination : %15.12lf\n", failp);
   free(SP); free(SN);
}

void morehead ( double H[], int n )
{
   double *P, *N, succp, failp;
   int i, t;

   /* The double arrays store probabilities of surplus of tails over heads so far.
      Only the surplus count matters, not how it is achieved.
      P - previous iteration. N - next iteration. */
   P = (double *)malloc((n + 1) * sizeof(double));
   N = (double *)malloc((n + 1) * sizeof(double));
   succp = 0;
   N[0] = 1; for (i=1; i<=n; ++i) N[i] = 0;
   for (t=0; t<n; ++t) {
      succp += N[0] * H[t];
      for (i=0; i<=n; ++i) { P[i] = N[i]; N[i] = 0; }
      for (i=1; i<=t; ++i) N[i-1] += P[i] * H[t];
      for (i=0; i<=t; ++i) N[i+1] += P[i] * (1 - H[t]);
   }
   /* failp = 1 - succp. But I compute failp more stylistically in O(n) time.
      We have already spent O(n^2) time anyway. */
   failp = 0; for (i=0; i<=n; ++i) failp += N[i];
   printf("    Probability of successful termination   : %15.12lf\n", succp);
   printf("    Probability of unsuccessful termination : %15.12lf\n", failp);
   free(P); free(N);
}

void notreecomp ( double H[], int n, int parity )
{
   printf("\n+++ Computation without the tree\n");
   if (parity == EVENPC) equalhead(H,n);
   if (parity == ODDPC) morehead(H,n);
}

int main ( int argc, char *argv[] )
{
   int n, i;
   double *H;
   bintree T0, T1;

   srand((unsigned int)time(NULL));
   if (argc > 1) n = atoi(argv[1]); else scanf("%d", &n);

   H = genprobs(n);

   printf("*************************************************************************************\n");
   printf("*********************************   E V E N   P C   *********************************\n");
   printf("*************************************************************************************\n");
   printf("n = %d\n", n);
   printf("Head probabilities:"); for (i=0; i<n; ++i) printf(" %5.3lf", H[i]); printf("\n");
   T0 = buildtree(H,n,0,0,1,EVENPC);
   allevents(T0,n,EVENPC);
   extremeevents(T0,n,EVENPC);
   successprob(T0,EVENPC);
   notreecomp(H,n,EVENPC);
   printf("\n");

   printf("*************************************************************************************\n");
   printf("**********************************   O D D   P C   **********************************\n");
   printf("*************************************************************************************\n");
   printf("n = %d\n", n);
   printf("Head probabilities:"); for (i=0; i<n; ++i) printf(" %5.3lf", H[i]); printf("\n");
   T1 = buildtree(H,n,0,0,1,ODDPC);
   allevents(T1,n,ODDPC);
   extremeevents(T1,n,ODDPC);
   successprob(T1,ODDPC);
   notreecomp(H,n,ODDPC);
   printf("\n");

   exit(0);
}