/*@@@@@@@@@@@@@@@@@@@@@@@@ ngop.cpp @@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
 Random Option Pricing
 Summer 2002
 Authors: A. Korobchevsky
          S. Ferrando
 Ryerson University ALL RIGHTS RESERVED
 
 Purpose: Generic Option Calculator

 Input file format: 
		sigma: 0.5
		r:     0.2
		T:     0.3
		N:     40
		S:     20
		M0:    20000
		gap:   300
		steps: 100
 Other parameters are prompted by the user
 Command Line Arg. : random generator to be used (mt, taus, gsfr) 
 @@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@*/

#include <fstream>
#include <typeinfo>
#include <vector>
#include <ctime>
#include <iterator>
#include <iostream>
using namespace std;

#include "../Src/Alg/BasicMC.h"
#include "../Src/Alg/LocalMC.h"
#include "../Src/Alg/ExactValues.h"
#include "../Src/Und/Underlying.h"
#include "../Src/Und/Binom.h"
#include "../Src/Und/Continium.h"
#include "../Src/Rng/GSL_RandomGen.h"


int main(int argc, char **argv)
{ // Let default RNG be mnt
  char defArg[] = "mt";
  char arg[10];
  if (argc == 1)
    strcpy(arg, defArg);
  else if (argc == 2)
    strcpy(arg, argv[1]);
  else
  { cerr << "Extra command line args\n";
    exit(1);
  }

  const char results[20] = "plots.ps ";
  char resCommand[100] = "./extractor ";
  strcat(resCommand, results);

  // Compute binom & continium for now
  ifstream in("input", ios::in);
  string garbage;
  int S, M, Nfinal, M0, gap, steps;
  double sigma;

  Underlying* ptr_binom = new Binom(arg);
  in >> garbage;  in >> sigma;  
  Underlying* ptr_cont = new Continium(arg, sigma);

              
  // Input text file required with fixed format
  in >> garbage;  in >> ptr_binom->r;   ptr_cont->r = ptr_binom->r;   
  in >> garbage;  in >> ptr_binom->T;   ptr_cont->T = ptr_binom->T;  
  in >> garbage;  in >> Nfinal;   
  in >> garbage;  in >> S;     // Starting option value, needed for algorithm  
  in >> garbage;  in >> M0;     // Number of simulations, redundant for Exact Values Algorithm  
  in >> garbage;  in >> gap; 
  in >> garbage;  in >> steps; 

  char temp[20];
  sprintf(temp, "%d ", Nfinal);
  strcat(resCommand, temp);

  M = M0 + gap * steps;  // Compute total number of simulations

  ofstream out("bmc", ios::out);
  out << "##############Option values###################\n";        /* Start of Header */
  out << "Underlying: Binom\n"; 
  out << "Derivative: Strike\n"; 
    
  out << "Algorithm: Basic Monte Carlo\n"; 
  out << "U = " << ptr_binom->U;    
  out << "\nD = " << ptr_binom->D;    
  out << "\nr = " << ptr_binom->r;    
  out << "\nT = " << ptr_binom->T;   
  out << "\nN = " << Nfinal;   
  out << "\nP = " << ptr_binom->P;  
  out << "\nS = " << S;
  out << "\nM = " << M;
  out << "\n#############################################" << "$" << endl;

  ofstream outc("cmc", ios::out);
  outc << "##############Option values###################\n";        /* Start of Header */
  outc << "Underlying: Continium\n"; 
  outc << "Derivative: Strike\n"; 
    
  outc << "Algorithm: Basic Monte Carlo\n"; 
  outc << "U = " << ptr_cont->U;    
  outc << "\nD = " << ptr_cont->D;    
  outc << "\nr = " << ptr_cont->r;    
  outc << "\nT = " << ptr_cont->T;   
  outc << "\nN = " << Nfinal;   
  outc << "\nP = " << ptr_cont->P;  
  outc << "\nS = " << S;
  outc << "\nM = " << M;
  outc << "\n#############################################" << "$" << endl;

  in.close();


  Algorithm* ptr_a = new BasicMC(S, M, M0, gap);

  vector<double> ncv_binom;
  vector<double> nse_binom;
  vector<double> ncv_cont;
  vector<double> nse_cont;

  bool first = true;

  // Compute option value
  for (int N = 1; N < Nfinal; N++)
  { vector<double> cv;
    vector<double> se;
    long start = clock();

  
    ptr_cont->N = ptr_binom->N = N;  
    // Compute probabilities 
    ptr_binom->U = ptr_cont->U = exp(sigma * sqrt(ptr_cont->T / N));
    ptr_binom->D = ptr_cont->D = 1 / ptr_cont->U;
    ptr_binom->P = ptr_cont->P = (ptr_cont->U - exp(ptr_cont->r * ptr_cont->T / N)) / (ptr_cont->U - ptr_cont->D);

    Derivative* ptr_der = new AverageStrikeEuropeanCall(N);
    cout << "N = " << N << ", ptr_binom->U = " << ptr_binom->U << endl;
    ptr_a->calcResults(ptr_binom, ptr_der, cv, se);
    cout << "  Strike Value " << cv.back() << ", standard error " << se.back() <<  ", time = " << (clock() - start)/1000000.0 << "sec\n";
    ncv_binom.push_back(cv.back());    nse_binom.push_back(cv.back());

    if (first)
    { sprintf(temp, "1 %d ", Nfinal + 1);
      strcat(resCommand, temp);
    }
    first = false;

    cv.clear();
    se.clear();
    ptr_a->calcResults(ptr_cont, ptr_der, cv, se);
    cout << "  Strike Value " << cv.back() << ", standard error " << se.back() <<  ", time = " << (clock() - start)/1000000.0 << "sec\n";
    ncv_cont.push_back(cv.back());    nse_cont.push_back(se.back());
  }


  copy(ncv_binom.begin(), ncv_binom.end(), ostream_iterator<double>(out, " "));
  out.close();
  copy(ncv_cont.begin(), ncv_cont.end(), ostream_iterator<double>(outc, " "));
  outc.close();

  // Invoke the plotter
  strcat(resCommand, "bmc cmc");
  system(resCommand);
  system("gv plots.ps &");

  cout << resCommand << endl;
}
