// Generate simple multivariate data according to a specified
// hypothesis (signal or backround) or as a mixture.
// Glen Cowan, RHUL, Physics, November 2007

#include <iostream>
#include <fstream>
#include <cstdlib>
#include <cmath>
#include <TRandom3.h>
#include "Event.h"

using namespace std;

int main(){

// Open output files

  const int numFiles = 2;           // one for signal, one for background
  string fileName[numFiles];
  fileName[0] = "background.txt";
  fileName[1] = "signal.txt";

// Set up to generate random values

// Set up to generate random values

  int seed = 12345;
  TRandom3* ran = new TRandom3(seed);
  const int numValues = 10000;
  const double r0 = 2.0;
  // historical defaults
  const double sigma = 0.5;
  const double a = 0.7;

  // try for more overlap
  // const double sigma = 0.7;
  // const double a = 0.9;

  for (int j=0; j<numFiles; j++){      // j = 0 bkg, j = 1 sig

    ofstream outFile;
    outFile.setf(ios::fixed);
    outFile.precision(8);

// generate  events

    if ( j == 0 ) {           // background

      outFile.open(fileName[j].c_str());
      for (int i=0; i<numValues; ++i){
        double x = sigma*ran->Gaus();
        double y = sigma*ran->Gaus();
        double r = ran->Rndm();
        double z = ran->Rndm();
        outFile << x << "\t" << y << "\t" << z << endl;
      }
      outFile.close();

    }
    else if ( j == 1 ) {

      outFile.open(fileName[j].c_str());
      for (int i=0; i<numValues; ++i){
        double r = ran->Rndm();
        double theta = r*a*M_PI;
        double x0 = r0*cos(theta);
        double y0 = r0*sin(theta);
        double x = x0 + sigma*ran->Gaus();
        double y = y0 + sigma*ran->Gaus();
        double u = ran->Rndm();
        double z = pow(u, 1./(1.+2.*theta));
        outFile << x << "\t" << y << "\t" << z << endl;
      }
      outFile.close();

    }

  }           // loop over j = 0, 1  (bkg, sig)

  return 0;

}