// ----------------------------------------------------------------------------
// $Id: root_GoGenMC.C,v 1.8 2001/06/19 19:45:40 masako Exp $
// ----------------------------------------------------------------------------
// $Log: root_GoGenMC.C,v $
// Revision 1.8  2001/06/19 19:45:40  masako
// Remove cross section.
//
// Revision 1.7  2001/05/16 17:46:08  masako
// Add a comment to show the unit of cross section (in fb)
//
// Revision 1.6  2001/05/11 21:52:34  toshi
// Remove .so suffix from gSystem->Load("...") in order to call macros from
// Windows system.
//
// Revision 1.5  2001/05/09 21:59:52  toshi
// Add $Id: root_GoGenMC.C,v 1.8 2001/06/19 19:45:40 masako Exp $.
//
// Revision 1.4  2001/05/04 23:25:07  masako
// Remove Comment out to store the LCDBeam/LCDEventHeader...
//
// Revision 1.2  2001/05/01 20:32:01  masako
// Bug fix.
// Correct typo miss.
//
//
//  Example file for event generation interactively with Root.
//  
//  In this file, 1) event generation (e+e- -> ZH) with pandora_pythia
//                2) out put a root file
//        
//                               Apr. 26 2001  Toshi
//
//   :: How to run ::
// 
//   1) type root
//   2)
//      root [0]  .x root_GoGenMC.C
//
// For example, if you like to generate 20 events, type like
//      root [0]  .x root_GoGenMC.C(20)

//__________________________________________________________________________
Int_t root_GoGenMC(Int_t   nEvent  = 100){

  gROOT->Reset();

  // In order to refer to shareable libraries in this simple form,
  // include their paths in a .rootrc file.

  gSystem->Load("libPandora");
  gSystem->Load("libPandoraProc");
  gSystem->Load("libPythia6");
  gSystem->Load("libEG");
  gSystem->Load("libEGPythia6");
  gSystem->Load("libTauola");
  gSystem->Load("libPandPyth");
  gSystem->Load("libLCDEvent");
  //gSystem->Load("libLCDRootAppsUtil");
  //gSystem->Load("libLCDFastMC");
  //gSystem->Load("libLCDPhUtil");
  gSystem->Load("libLCDGenUtil");


  //------------ Out put file (Root)
  Char_t* outputfilename="test.root";
  TFile* ofile = new TFile(outputfilename,"RECREATE","LCD Event ROOT file");

  LCDEvent* event=new LCDEvent();
  Int_t comp    = 1;       // by default file is compressed
  Int_t split   = 1;       // by default, split Event in sub branches
  Int_t bufsize = 64000;
  ofile->SetCompressionLevel(comp);
  TTree* tree = new TTree("T","LCD ROOT tree");
  tree->SetAutoSave(1000000000);  // autosave when 1 Gbyte written      
  tree->Branch("LCDEvent","LCDEvent", &event, bufsize,split);

  //Pyjet -> McPart converter
  LCDPyj2McPart py2mc;


  //------------< Setup Generator (Pandora-Pythia) >------------------------
  cout<<"Start Generator Setup!"<<endl;

#include "pandpythconsts.h"

  Int_t iseed_pan = 7;  // default pandora initial seed

  /*    define the pandora event selection in this space     */
  
  double ECM  = 500.;
  double PolE = -0.8; // Polarization Electron
  double PolP =  0.0; // Polarization Positron

  ebeam b1(ECM/2.0, PolE, electron, electron);
  b1.setup(NLC500);
  b1.ISRoff();            // Turn off initial-state radiation
  b1.beamstrahlungoff();  // Turn off beamstrahlung

  ebeam b2(ECM/2.0, PolP, positron, positron);
  b2.setup(NLC500);
  b2.ISRoff();
  b2.beamstrahlungoff();

  double Hmass = 120.; // Higgs Mass
//double Hmass = 140.; // Higgs Mass
  eetoZHiggs pr(Hmass);

  //For example:: Z0 decays to anything. H0 decays to bb-bar only.
  pr.onlyDecay(allStates,bOnly); 

  pandora P(b1,b2,pr);
  
  pandorarun  PR(&P, epluseminus, ECM, nEvent, iseed_pan);
  
  //----- For PYTHIA ------
  //Or you can also change PYTHIA seed here.
  int iseed_pyth = 22846; PR.newseed(iseed_pyth); // Change Seed

  //PR.partonshoweroff();   // Turn off parton shower
  //PR.fragmentationoff();  // Turn off fragmentation
  //PR.FSRoff();            // Turn off final-state QCD QED radiation

  PR.StdHepIOoff();  // Don't use stdhep IO

  //       end of definition

  PR.initialize(0);
  
  // If you want/don't want to decay some specific particle(s),
  // you can call PartDecayOn/PartDecayOff here.   
  // For example: Do not allow to decay K0S and Lambda
  //int ipart1 = 310;  // K0S
  //int ipart2 = 3122; // Lambda
  //PR.PartDecayOff(ipart1);
  //PR.PartDecayOff(ipart2);

  cout<<"End of Generator Setup!"<<endl;
  //------------< End Generator (Pandora-Pythia) setup>-------------------

  gBenchmark->Start("gen_eetoZH");

  Int_t iEvent;
  for (iEvent=0 ; iEvent < nEvent ; iEvent++) {
    PR->getevent();//Event Generation

    event->Clean();

    // Make McPart
    py2mc.MakeMcPart(event->MCparticles());

    // Store Event Header
    event->EventHeader()->SetEvtNum(iEvent);
    event->EventHeader()->SetProcessName(pr->name);

    // Store beam information
    event->BEAM()->SetEcm(ECM);
    event->BEAM()->SetPolElec(PolE);
    event->BEAM()->SetPolPosi(PolP);

    tree->Fill();      // Fill information into output root file
  }

  gBenchmark->Stop("gen_eetoZH");
  gBenchmark->Print("gen_eetoZH");

  PR.terminate();

  ofile->Write();
  tree->Print();
  ofile->Close();

  return iEvent;
}