// finance++

/* TimeHalvingTrinom *****************/

struct IP_TimeHalvingTrinom
{
  double S0;
  unsigned int N;
  vector< double > dt;
  double sigma;
  double r;
};

typedef struct IP_TimeHalvingTrinom IP_TimeHalvingTrinom;

class TimeHalvingTrinom
{
 protected:
  static const unsigned int Q = 3;
  static const int K = 0;
  static const int I = 2;
  static const int J = 2;
  bool recomb;

  vector< vector<double> > M;
  vector< vector<double> > P;
  vector< double > disc;

  double S0;
  unsigned int N;
  vector< double > dt;

 public:
  TimeHalvingTrinom(IP_TimeHalvingTrinom inIP)
    {
      if(! (inIP.dt.size() == 1 || inIP.dt.size() == (inIP.N - 1)) )
	{
	  fprintf(stderr, "TimeHalvingTrinom invalid dt size %d, must be 1 or N-1.\n", inIP.dt.size());
	  return;
	}
      
      if(inIP.sigma < 0.0)
	{
	  fprintf(stderr, "TimeHalvingTrinom invalid sigma %f, cannot be less than 0.\n", inIP.sigma);
	  return;
	}

      if(inIP.dt.size() == 1)
	recomb = true;
      else
	recomb = false;

      S0 = inIP.S0;
      N = inIP.N;
      dt = vector<double>(inIP.dt.size());
      
      double p;

      for(unsigned int i = 0; i < inIP.dt.size(); i++)
	{
	  dt[i] = inIP.dt[i];

	  vector<double> vd = vector<double>(3);
      
	  vd[0] = pow(M_E, -2 * inIP.sigma * sqrt(dt[i] / 2));
	  vd[1] = (double)1.0;
	  vd[2] = pow(M_E, 2 * inIP.sigma * sqrt(dt[i] / 2));
      
	  M.push_back(vd);

	  p = (pow(M_E, inIP.r * (dt[i] / 2)) - pow(M_E, -inIP.sigma * sqrt(dt[i] / 2))) /
	    (pow(M_E, inIP.sigma * sqrt(dt[i] / 2)) - pow(M_E, -inIP.sigma * sqrt(dt[i] / 2)));
      
	  vd = vector<double>(3);

	  vd[0] = pow(1 - p, 2);
	  vd[1] = 2 * p * (1 - p);
	  vd[2] = pow(p, 2);
      
	  P.push_back(vd);
	  
	  disc.push_back(pow(M_E, (-inIP.r) * dt[i]));
	}
      
    }
};


