/* file: probability-tables.c (version 0.0)
   compilation: gcc -lm probability-tables.c
   outputs: xxx.table (phan-TeX tables)
   last modification: december 9, 2009
   use: academic purpose?
   author: Anthony Phan
*/
#include <locale.h>
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include "probability-distributions.h"

FILE *output;
double rounding = 1e4;

double rnd(double x){
  return floor(rounding*x+0.5)/rounding;
}

/* Ci-dessous nous trouvons les commandes de mise en forme du fichier
   de sortie. Les s\'eparer du reste a pour avantage de rendre la
   programmation plus claire et plus ais\'ee \`a modifier.
   Notre syntaxe pour les tables est en gros~:
   \table[options]
   \tr\td...\endtr
   .
   .
   .
   \endtable
   c'est fait pour \^etre semblable \`a de l'HTML/TeX simplifi\'e.
*/

table(){
  fprintf(output,
	  "\\table[align=c,width=\\hsize,vrules=t,cellspacing=2pt]\n");
  fprintf(output, "\\noalign{\\hrule}\n");
}

hr(){
  fprintf(output, "\\noalign{\\hrule}\n");
}

tr(){
  fprintf(output, "\\tr");
}

td(char* x){
  fprintf(output, "\\td %s", x);
}
tdshort(double x){
  fprintf(output, "\\td%.2f", x);
}

tdlong(double x){
  fprintf(output, "\\td%.4f", x);
}

endtr(){
  fprintf(output, "\\endtr\n");
}

endtable(){
  fprintf(output,"\\noalign{\\hrule}\n\\endtable");
}

/* Des exemples de tables\dots\ \`A modifier si besoin est\dots  */

normalcdftable(){
  int i, j;
  output = fopen("normalcdf.table", "w");
  table();
  tr(); td("\\enlargedstrut[0.5ex,0.5ex]$z$");
  for (j = 0; j <= 9; j++) tdshort(0.01*j); endtr();
  hr();
  for (i = 0; i <= 29; i++){
    tr();
    if (i % 10 == 0)
      if (i == 0)
	fprintf(output, "\\td\\enlargedstrut[0.5ex,0.0ex]%.1f", 0.1*i);
      else
	fprintf(output, "\\td\\enlargedstrut[1ex,0.0ex]%.1f", 0.1*i);
    else fprintf(output, "\\td%.1f", 0.1*i);
    for  (j = 0; j <= 9; j++)
      tdlong(rnd(normalcdf(0.1*i+0.01*j)));
    endtr();
  }
  endtable();
  fclose(output);
}

normalicdftable(){
  int i, j;
  output = fopen("normalicdf.table", "w");
  table();
  tr(); td("\\enlargedstrut[0.5ex,0.5ex]$\\alpha$");
  for (j = 0; j <= 9; j++) tdshort(0.01*j); endtr();
  hr();
  for (i = 5; i <= 9; i++){
    tr();
    if (i == 5) fprintf(output, "\\td\\enlargedstrut[0.5ex,0ex]%.1f", 0.1*i);
    else
      if (i == 9) fprintf(output, "\\td\\enlargedstrut[0ex,0.5ex]%.1f", 0.1*i);
      else fprintf(output, "\\td%.1f", 0.1*i);
    for  (j = 0; j <= 9; j++)
      tdlong(rnd(normalicdf(0.1*i+0.01*j)));
    endtr();
  }
  fprintf(output,"\\noalign{\\hrule\\smallskip\\hrule}\n");
  tr(); td("\\enlargedstrut[0.5ex,0.5ex]$\\alpha$");
  for (j = 0; j <= 9; j++) fprintf(output, "\\td %.3f", 0.99+0.001*j);
  endtr(); hr();
  tr(); td("\\enlargedstrut[0.5ex,0.5ex]$\\Phi^{-1}(\\alpha)$");
  for  (j = 0; j <= 9; j++)
    tdlong(rnd(normalicdf(0.99+0.001*j)));
  endtr();
  fprintf(output,"\\noalign{\\hrule\\smallskip\\hrule}\n");
  tr(); td("\\enlargedstrut[0.5ex,0.5ex]$\\alpha$");
  for (j = 0; j <= 9; j++) fprintf(output, "\\td %.4f", 0.999+0.0001*j);
  endtr(); hr();
  tr(); td("\\enlargedstrut[0.5ex,0.5ex]$\\Phi^{-1}(\\alpha)$");
  for  (j = 0; j <= 9; j++)
    tdlong(rnd(normalicdf(0.999+0.0001*j)));
  endtr();
  endtable();
  fclose(output);
}

normalicdfbistable(){
  int i, j;
  output = fopen("normalicdfbis.table", "w");
  table();
  tr(); td("\\enlargedstrut[0.5ex,0.5ex]$\\alpha$");
  for (j = 0; j <= 9; j++) tdshort(0.01*j); endtr();
  hr();
  for (i = 0; i <= 9; i++){
    tr();
    if (i == 0) fprintf(output,"\\td\\enlargedstrut[0.5ex,0ex]%.1f",0.1*i);
    else if (i == 9)
      fprintf(output,"\\td\\enlargedstrut[0ex,0.5ex]%.1f",0.1*i);
    else fprintf(output,"\\td%.1f",0.1*i);
    for  (j = 0; j <= 9; j++)
      if (i == 0 && j == 0) td("$\\infty$");
      else tdlong(rnd(normalicdf(1-(0.1*i+0.01*j)/2)));
    endtr();
  }
  endtable();
  fclose(output);
}

normalicdftertable(){
  int i;
  double alpha[7];
  alpha[0] = 1e-3;
  alpha[1] = 1e-4;
  alpha[2] = 1e-5;
  alpha[3] = 1e-6;
  alpha[4] = 1e-7;
  alpha[5] = 1e-8;
  alpha[6] = 1e-9;
  output = fopen("normalicdfter.table", "w");
  table(); tr(); td("$\\alpha$");
  for (i = 0; i <= 6; i++) fprintf(output, "\\td %g", alpha[i]);
  endtr();
  hr();
  tr(); td("$z_{1-\\alpha/2}$");
  for (i = 0; i <= 6; i++) tdlong(rnd(normalicdf_(1-alpha[i]/2)));
  endtr();
  endtable();
  fclose(output);
}

chisquareicdftable(){
  int i, j;
  double alpha[9];
  alpha[0] = 0.99;
  alpha[1] = 0.975;
  alpha[2] = 0.95;
  alpha[3] = 0.90;
  alpha[4] = 0.10;
  alpha[5] = 0.05;
  alpha[6] = 0.025;
  alpha[7] = 0.01;
  alpha[8] = 0.001;
  output = fopen("chisquareicdf.table", "w");
  table();
  tr(); td("\\backslashedcell[25,15]{$\\nu$}{$\\alpha$}");
  for (j = 0; j <= 8; j++) fprintf(output, "\\td %.3f", alpha[j]);
  endtr();
  hr();
  for (i = 1; i <= 30; i++){
    tr();
    if (i % 10 == 1)
      if (i == 1)
	fprintf(output, "\\td\\enlargedstrut[0.5ex,0.0ex]%d", i);
      else
	fprintf(output, "\\td\\enlargedstrut[1ex,0.0ex]%d", i);
    else fprintf(output, "\\td%d", i);
    for  (j = 0; j <= 8; j++)
      tdlong(rnd(chisquareicdf(1.0-alpha[j], i)));
    endtr();
  }
  endtable();
  fclose(output);
}

studenticdftable(){
  int i, j;
  double alpha[9];
  alpha[0] = 0.90;
  alpha[1] = 0.50;
  alpha[2] = 0.30;
  alpha[3] = 0.20;
  alpha[4] = 0.10;
  alpha[5] = 0.05;
  alpha[6] = 0.02;
  alpha[7] = 0.01;
  alpha[8] = 0.001;
  output = fopen("studenticdf.table", "w");
  table();
  tr(); td("\\backslashedcell[25,15]{$\\nu$}{$\\alpha$}");
  for (j = 0; j <= 8; j++) fprintf(output, "\\td %.3f", alpha[j]);
  endtr();
  hr();
  for (i = 1; i <= 30; i++){
    tr();
    if (i % 10 == 1)
      if (i == 1)
	fprintf(output, "\\td\\enlargedstrut[0.5ex,0.0ex]%d", i);
      else
	fprintf(output, "\\td\\enlargedstrut[1ex,0.0ex]%d", i);
    else fprintf(output, "\\td%d", i);
    for  (j = 0; j <= 8; j++)
      tdlong(rnd(studenticdf(1.0-alpha[j]/2, i)));
    endtr();
  }
    tr(); fprintf(output, "\\td\\enlargedstrut[1ex,0.0ex]40");
    for  (j = 0; j <= 8; j++)
      tdlong(rnd(studenticdf(1.0-alpha[j]/2, 40)));
    endtr();
    tr(); fprintf(output, "\\td 60");
    for  (j = 0; j <= 8; j++)
      tdlong(rnd(studenticdf(1.0-alpha[j]/2, 60)));
    endtr();
    tr(); fprintf(output, "\\td 80");
    for  (j = 0; j <= 8; j++)
      tdlong(rnd(studenticdf(1.0-alpha[j]/2, 80)));
    endtr();
    tr(); fprintf(output, "\\td 120");
    for  (j = 0; j <= 8; j++)
      tdlong(rnd(studenticdf(1.0-alpha[j]/2, 120)));
    endtr();
    fprintf(output,"\\noalign{\\hrule\\smallskip\\hrule}\n");
    tr(); fprintf(output, "\\td $\\infty$");
    for  (j = 0; j <= 8; j++)
      tdlong(rnd(normalicdf(1.0-alpha[j]/2)));
    endtr();
  endtable();
  fclose(output);
}

fishericdftable(){
  int i, j; double x, alpha;
  alpha = 0.05;
  output = fopen("fishericdf.table", "w");
  table();
  tr(); td("\\backslashedcell[25,15]{$\\nu_2$}{$\\nu_1$}");
  j = 1;
  while (j <= 30){
    fprintf(output, "\\td %d", j);
    if (j < 6) j = j+1;
    else if (j < 10) j = j+2;
    else if (j < 20) j = j+5; else j = j+10;
  }
  endtr(); hr();
  i = 1;
  while (i <= 100){
    tr(); fprintf(output, "\\td %d", i); j = 1;
    while (j <= 30){
      x = rnd(fishericdf(1.0-alpha, j, i));
      if (x >= 100) fprintf(output, "\\td %.0f", x);
      else if (x >= 10) fprintf(output, "\\td %.1f", x);
      else fprintf(output, "\\td %.2f", x);
      if (j < 6) j = j+1;
      else if (j < 10) j = j+2;
      else if (j < 20) j = j+5; else j = j+10;
    }
    endtr();
    if (i < 20) i = i+1;
    else if (i < 30) i = i+2;
    else if (i < 60) i = i+10; else i = i+20;
  }
  endtable();
  fclose(output);
}

kolmogorovcdftable(int n){
  int i, j;
  output = fopen("kolmogorovcdf.table", "w");
  table();
  tr(); td("\\backslashedcell[25,15]{$n$}{$k$}");
  for (j = 1; j <= 9; j++) fprintf(output, "\\td %.1f", 0.1*j);
  endtr(); hr();
  for (i=1; i <= n; i++){
    fprintf(output, "\\tr\\td%d", i);
    for (j = 1; j <= 9; j++)
      tdlong(rnd(kolmogorovcdf(0.1*j, i)));
    endtr();
  }
  endtable();
  fclose(output);
}

kolmogorovicdfrow(double p, int i, int m, int step){
  int j;
  tr(); tdshort(p);
  for (j = 0; j < m; j++)
    if (i+j*step <= 0) tdlong(1.0000);
    else tdlong(rnd(kolmogorovicdf(1-p, i+j*step)));
  endtr();
}

/* m columns, n rows */

kolmogorovicdftable(){
  int i, j, step;
  output = fopen("kolmogorovicdf.table", "w");
  i = 0; step = 1;
  while (i < 200){
    table();
    tr(); td("\\backslashedcell[25,15]{$\\alpha$}{$n$}\\endtd\n");
    for (j = 0; j < 10; j++) fprintf(output, "\\td %d", i+j*step);
    endtr();
    hr();
    kolmogorovicdfrow(0.01, i, 10, step);
    kolmogorovicdfrow(0.05, i, 10, step);
    kolmogorovicdfrow(0.1, i, 10, step);
    kolmogorovicdfrow(0.15, i, 10, step);
    kolmogorovicdfrow(0.20, i, 10, step);
    endtable();
    fprintf(output, "\\smallbreak");
    i = i+10*step;
    if (i < 40) step = 1;
    else if (i < 60) step = 2;
    else if (i < 80) step = 5;
    else if (i < 200) step = 10;
    else step = 50;
  }
  fclose(output);
}

klmcdftable(int n){
  int i, j;
  output = fopen("klmcdf.table", "w");
  table();
  tr(); td("\\backslashedcell[25,15]{$n$}{$k$}");
  for (j = 1; j <= 9; j++) fprintf(output, "\\td %.1f", 0.1*j);
  endtr(); hr();
  for (i=1; i <= n; i++){
    fprintf(output, "\\tr\\td%d", i);
    for (j = 1; j <= 9; j++)
      tdlong(rnd(klmcdf(0.1*j, i)));
    endtr();
  }
  endtable();
  fclose(output);
}

klmicdfrow(double p, int i, int m, int step){
  int j;
  tr(); tdshort(p);
  for (j = 0; j < m; j++)
    if (i+j*step <= 0) tdlong(1.0000);
    else tdlong(rnd(klmicdf(1-p, i+j*step)));
  endtr();
}

/* m columns, n rows */

klmicdftable(){
  int i, j, step;
  output = fopen("klmicdf.table", "w");
  i = 0; step = 1;
  while (i < 200){
    table();
    tr(); td("\\backslashedcell[25,15]{$\\alpha$}{$n$}\\endtd\n");
    for (j = 0; j < 10; j++) fprintf(output, "\\td %d", i+j*step);
    endtr();
    hr();
    klmicdfrow(0.01, i, 10, step);
    klmicdfrow(0.05, i, 10, step);
    klmicdfrow(0.1, i, 10, step);
    klmicdfrow(0.15, i, 10, step);
    klmicdfrow(0.20, i, 10, step);
    endtable();
    fprintf(output, "\\smallbreak");
    i = i+10*step;
    if (i < 40) step = 1;
    else if (i < 60) step = 2;
    else if (i < 80) step = 5;
    else if (i < 200) step = 10;
    else step = 50;
  }
  fclose(output);
}

kpmicdftable(){
  int i, j;
  double p[11];
  p[0] = 0.01;
  p[1] = 0.05;
  p[2] = 0.10;
  p[3] = 0.20;
  p[4] = 0.25;
  p[5] = 0.50;
  p[6] = 0.75;
  p[7] = 0.80;
  p[8] = 0.90;
  p[9] = 0.95;
  p[10] = 0.99;
  output = fopen("kpmicdf.table", "w");
  table();
  tr(); td("\\backslashedcell[25,15]{$n$}{$\\alpha$}");
  for (j = 0; j <= 10; j++) fprintf(output, "\\td %.2f", 1.0-p[j]);
  endtr();
  hr();
  for (i = 1; i <= 30; i++){
    tr();
    if (i % 10 == 1)
      if (i == 1)
	fprintf(output, "\\td\\enlargedstrut[0.5ex,0ex]%d", i);
	else fprintf(output, "\\td\\enlargedstrut[1ex,0ex]%d", i);
    else if (i == 30) fprintf(output, "\\td\\enlargedstrut[0ex,0.5ex]%d", i);
    else fprintf(output, "\\td%d", i);
    for  (j = 0; j <= 10; j++)
      tdlong(rnd(sqrt(i)*kpmicdf(p[j], i)));
    endtr();
  }
  hr();
  tr(); fprintf(output, "\\td\\enlargedstrut[0.5ex,0.5ex]
$n>30$\\td[colspan=11\\relax]
$y_p-{1\\over6}n^{-1/2}+O(1/n)$, avec $y_p^2={1\\over2}\\ln(1/(1-p))$");
  endtr(); hr();
  tr(); fprintf(output, "\\td\\enlargedstrut[0.5ex,0.5ex]$y_p{}$");
  for  (j = 0; j <= 10; j++)
    tdlong(rnd(sqrt(0.5*log(1.0/(1.0-p[j])))));
  endtr();
  endtable();
  fclose(output);
}

main(){
  setlocale(LC_NUMERIC,"fr_FR");
  normalcdftable();
  normalicdftable();
  normalicdfbistable();
  normalicdftertable();
  chisquareicdftable();
  studenticdftable();
  fishericdftable();
  klmcdftable(60);
  klmicdftable();
  kolmogorovicdftable();
  kpmicdftable();
}