// ----------------------------------------------------------------------------
// $Id: pythia_GoGenMC.C,v 1.1 2001/06/29 21:40:25 toshi Exp $
// ----------------------------------------------------------------------------
// $Log: pythia_GoGenMC.C,v $
// Revision 1.1  2001/06/29 21:40:25  toshi
// Initial release.
//
//
//  Example file to generate mc events interactively with Root.
//
//  In this file, 1) event generation with Pythia (TPythia6)
//
//  For example, generate e+e- --> qq-bar events
//
//                               June 29 2001  T.Abe
//
//   :: How to run ::
// 
//   1) type root
//   2)
//      root [0]  .x pythia_GoGenMC.C
//
// For example, if you like to run 20 events, type like
//      root [0]  .x pythia_GoGenMC.C(1000)
//    


Int_t pythia_GoGenMC(Int_t nEvent=300){
  
  gROOT->Reset();
  
  // In order to refer to shareable libraries in this simple 
  // form, include their paths in a .rootrc file.
  
  if (strcmp(gSystem->GetName(),"WinNT")) {
    gSystem->Load("libPythia6");
  } else {
    gSystem->Load("libPandora");
    gSystem->Load("libPandPyth");
  }
  gSystem->Load("libEG");
  gSystem->Load("libEGPythia6");
  gSystem->Load("libLCDEvent");
  gSystem->Load("libLCDRootAppsUtil");
  gSystem->Load("libLCDFastMC");
  gSystem->Load("libLCDPhUtil");
  gSystem->Load("libLCDGenUtil");

  //------------< Setup Generator (TPythia6) >------------------------
  printf("Start Generator Setup!\n");

  double ECM  = 91.26;
  int    iseed_pyth =  19518503;// pythia initial seed

  TPythia6 pyth6;

  //...t mass.
  pyth6.SetPMAS(6,1,175.);

  // Select gamma*/Z0 production process
  pyth6.SetMSEL(0);
  pyth6.SetMSUB(1,1);

  // Switch off all Z decay modes except qq-bar.
  int MDCY_23_2 = pyth6.GetMDCY(23,2);
  int MDCY_23_3 = pyth6.GetMDCY(23,3);
  for (int IDCY = MDCY_23_2; IDCY< MDCY_23_2+MDCY_23_3; IDCY++){
    int partID = pyth6.GetKFDP(IDCY,1);
    if (TMath::Abs(partID) > 5 ) pyth6.SetMDME(IDCY,1,0); 
    // if (TMath::Abs(partID) != 5 ) pyth6.SetMDME(IDCY,1,0); 
  }

  pyth6.Pyinit("CMS","E-","E+",ECM);
 
  // Change pythia seed  
  pyth6.SetMRPY(1,iseed_pyth); 

  //Pyjet -> McPart converter
  LCDPyj2McPart py2mc;

  printf("End of Generator Setup!\n");
  //------------< End Generator setup>-------------------

  LCDEvent *event=new LCDEvent();

  //
  // Open root file.
  //
  Int_t comp=2;
  Char_t* rootFile = "eetobbZ0.root";
  TFile* hfile = new TFile(rootFile,"RECREATE","LCD Event ROOT file");
  hfile->SetCompressionLevel(comp);

  //
  // Create a tree with one superbranch.
  //
  TTree* tree = new TTree( "T", "LCD ROOT tree" );
  tree->SetAutoSave( 1000000000 );
  tree->Branch( "LCDEvent", "LCDEvent", &event);

  // Event loop
  Int_t iEvent;
  for (iEvent = 0; iEvent < nEvent; iEvent++) {

    pyth6.GenerateEvent();                //Event Generation with pythia
    //if (iEvent < 3) pyth6.Pylist(1);      //Print Event list for example

    event.MCparticles()->Delete();        //prepare for making McPart
    py2mc.MakeMcPart(event.MCparticles());//Make McPart
 
    if (iEvent < 5) pyth6.Pylist(1);

    tree->Fill();
  }

  //
  // Clean up the ROOT file.
  //
  hfile->Write();
  tree->Print();
  hfile->Close();

  pyth6.Pystat(1); // Pythia summary

  return iEvent;
}