macros/gen/Fun4All_Pi0Gen_sPHENIX.C

0001 #ifndef MACRO_FUN4ALLG4SPHENIX_C
0002 #define MACRO_FUN4ALLG4SPHENIX_C
0003
0004 #include <../../pi0Efficiency.h>
0005 #include <GlobalVariables.C>
0006
0007 #include <DisplayOn.C>
0008 #include <G4Setup_sPHENIX.C>
0009 #include <G4_Bbc.C>
0010 #include <G4_CaloTrigger.C>
0011 #include <G4_Centrality.C>
0012 #include <G4_DSTReader.C>
0013 #include <G4_Global.C>
0014 #include <G4_HIJetReco.C>
0015 #include <G4_Input.C>
0016 #include <G4_Jets.C>
0017 #include <G4_KFParticle.C>
0018 #include <G4_ParticleFlow.C>
0019 #include <G4_Production.C>
0020 #include <G4_TopoClusterReco.C>
0021 #include <G4_Tracking.C>
0022 #include <G4_User.C>
0023 #include <QA.C>
0024
0025 #include <fun4all/Fun4AllDstOutputManager.h>
0026 #include <fun4all/Fun4AllOutputManager.h>
0027 #include <fun4all/Fun4AllServer.h>
0028
0029 #include <phool/PHRandomSeed.h>
0030 #include <phool/recoConsts.h>
0031
0032 R__LOAD_LIBRARY(libfun4all.so)
0033 R__LOAD_LIBRARY(libpi0Efficiency.so)
0034 // For HepMC Hijing
0035 // try inputFile = /sphenix/sim/sim01/sphnxpro/sHijing_HepMC/sHijing_0-12fm.dat
0036
0037 int Fun4All_Pi0Gen_sPHENIX(
0038     const int nEvents = 1,
0039     const int pid = 111,
0040     const string &outputFile = "G4sPHENIX.root",
0041     const string &inputFile = "https://www.phenix.bnl.gov/WWW/publish/phnxbld/sPHENIX/files/sPHENIX_G4Hits_sHijing_9-11fm_00000_00010.root",
0042     const string &embed_input_file = "https://www.phenix.bnl.gov/WWW/publish/phnxbld/sPHENIX/files/sPHENIX_G4Hits_sHijing_9-11fm_00000_00010.root",
0043     const int skip = 0,
0044     const string &outdir = ".")
0045 {
0046   Fun4AllServer *se = Fun4AllServer::instance();
0047   se->Verbosity(0);
0048
0049   //Opt to print all random seed used for debugging reproducibility. Comment out to reduce stdout prints.
0050   PHRandomSeed::Verbosity(0);
0051
0052   // just if we set some flags somewhere in this macro
0053   recoConsts *rc = recoConsts::instance();
0054   // By default every random number generator uses
0055   // PHRandomSeed() which reads /dev/urandom to get its seed
0056   // if the RANDOMSEED flag is set its value is taken as seed
0057   // You can either set this to a random value using PHRandomSeed()
0058   // which will make all seeds identical (not sure what the point of
0059   // this would be:
0060   //  rc->set_IntFlag("RANDOMSEED",PHRandomSeed());
0061   // or set it to a fixed value so you can debug your code
0062   //  rc->set_IntFlag("RANDOMSEED", 12345);
0063
0064   //===============
0065   // conditions DB flags
0066   //===============
0067   // tag
0068   rc->set_StringFlag("XPLOAD_TAG","sPHENIX_ExampleGT_1");
0069   // database config
0070   rc->set_StringFlag("XPLOAD_CONFIG","sPHENIX_cdb");
0071   // 64 bit timestamp
0072   rc->set_uint64Flag("TIMESTAMP",12345678912345);
0073
0074
0075   //===============
0076   // Input options
0077   //===============
0078   // verbosity setting (applies to all input managers)
0079   Input::VERBOSITY = 0;
0080   // First enable the input generators
0081   // Either:
0082   // read previously generated g4-hits files, in this case it opens a DST and skips
0083   // the simulations step completely. The G4Setup macro is only loaded to get information
0084   // about the number of layers used for the cell reco code
0085   //  Input::READHITS = true;
0086   INPUTREADHITS::filename[0] = inputFile;
0087   // if you use a filelist
0088   // INPUTREADHITS::listfile[0] = inputFile;
0089   // Or:
0090   // Use particle generator
0091   // And
0092   // Further choose to embed newly simulated events to a previous simulation. Not compatible with `readhits = true`
0093   // In case embedding into a production output, please double check your G4Setup_sPHENIX.C and G4_*.C consistent with those in the production macro folder
0094   // E.g. /sphenix/sim//sim01/production/2016-07-21/single_particle/spacal2d/
0095   //  Input::EMBED = true;
0096   INPUTEMBED::filename[0] = embed_input_file;
0097   // if you use a filelist
0098   //INPUTEMBED::listfile[0] = embed_input_file;
0099
0100   Input::SIMPLE = true;
0101   // Input::SIMPLE_NUMBER = 2; // if you need 2 of them
0102   Input::SIMPLE_VERBOSITY = 0;
0103
0104   //  Input::PYTHIA6 = true;
0105
0106   // Input::PYTHIA8 = true;
0107
0108   //  Input::GUN = true;
0109   //  Input::GUN_NUMBER = 3; // if you need 3 of them
0110   // Input::GUN_VERBOSITY = 1;
0111
0112   //D0 generator
0113   //Input::DZERO = false;
0114   //Input::DZERO_VERBOSITY = 0;
0115   //Lambda_c generator //Not ready yet
0116   //Input::LAMBDAC = false;
0117   //Input::LAMBDAC_VERBOSITY = 0;
0118   // Upsilon generator
0119   //Input::UPSILON = true;
0120   //Input::UPSILON_NUMBER = 3; // if you need 3 of them
0121   //Input::UPSILON_VERBOSITY = 0;
0122
0123   //  Input::HEPMC = true;
0124   INPUTHEPMC::filename = inputFile;
0125
0126   // Event pile up simulation with collision rate in Hz MB collisions.
0127   //Input::PILEUPRATE = 100e3;
0128
0129   //-----------------
0130   // Initialize the selected Input/Event generation
0131   //-----------------
0132   // This creates the input generator(s)
0133   InputInit();
0134
0135   //--------------
0136   // Set generator specific options
0137   //--------------
0138   // can only be set after InputInit() is called
0139
0140   // Simple Input generator:
0141   // if you run more than one of these Input::SIMPLE_NUMBER > 1
0142   // add the settings for other with [1], next with [2]...
0143   if (Input::SIMPLE)
0144   {
0145     INPUTGENERATOR::SimpleEventGenerator[0]->add_particles(pid, 1);
0146
0147     if (Input::HEPMC || Input::EMBED)
0148     {
0149       INPUTGENERATOR::SimpleEventGenerator[0]->set_reuse_existing_vertex(true);
0150       INPUTGENERATOR::SimpleEventGenerator[0]->set_existing_vertex_offset_vector(0.0, 0.0, 0.0);
0151     }
0152     else
0153     {
0154       INPUTGENERATOR::SimpleEventGenerator[0]->set_vertex_distribution_function(PHG4SimpleEventGenerator::Gaus,
0155                                                                                 PHG4SimpleEventGenerator::Gaus,
0156                                                                                 PHG4SimpleEventGenerator::Gaus);
0157       INPUTGENERATOR::SimpleEventGenerator[0]->set_vertex_distribution_mean(0., 0., 0.);
0158       INPUTGENERATOR::SimpleEventGenerator[0]->set_vertex_distribution_width(0.01, 0.01, 5.);
0159     }
0160     INPUTGENERATOR::SimpleEventGenerator[0]->set_eta_range(-0.3, 0.3);
0161     INPUTGENERATOR::SimpleEventGenerator[0]->set_phi_range(-M_PI, M_PI);
0162     INPUTGENERATOR::SimpleEventGenerator[0]->set_pt_range(0.5, 20.);
0163   }
0164   // Upsilons
0165   // if you run more than one of these Input::UPSILON_NUMBER > 1
0166   // add the settings for other with [1], next with [2]...
0167   if (Input::UPSILON)
0168   {
0169     INPUTGENERATOR::VectorMesonGenerator[0]->add_decay_particles("e", 0);
0170     INPUTGENERATOR::VectorMesonGenerator[0]->set_rapidity_range(-1, 1);
0171     INPUTGENERATOR::VectorMesonGenerator[0]->set_pt_range(0., 10.);
0172     // Y species - select only one, last one wins
0173     INPUTGENERATOR::VectorMesonGenerator[0]->set_upsilon_1s();
0174     if (Input::HEPMC || Input::EMBED)
0175     {
0176       INPUTGENERATOR::VectorMesonGenerator[0]->set_reuse_existing_vertex(true);
0177       INPUTGENERATOR::VectorMesonGenerator[0]->set_existing_vertex_offset_vector(0.0, 0.0, 0.0);
0178     }
0179   }
0180   // particle gun
0181   // if you run more than one of these Input::GUN_NUMBER > 1
0182   // add the settings for other with [1], next with [2]...
0183   if (Input::GUN)
0184   {
0185     INPUTGENERATOR::Gun[0]->AddParticle("pi-", 0, 1, 0);
0186     INPUTGENERATOR::Gun[0]->set_vtx(0, 0, 0);
0187   }
0188
0189   // pythia6
0190   if (Input::PYTHIA6)
0191   {
0192     //! apply sPHENIX nominal beam parameter with 2mrad crossing as defined in sPH-TRG-2020-001
0193     Input::ApplysPHENIXBeamParameter(INPUTGENERATOR::Pythia6);
0194   }
0195   // pythia8
0196   if (Input::PYTHIA8)
0197   {
0198     //! apply sPHENIX nominal beam parameter with 2mrad crossing as defined in sPH-TRG-2020-001
0199     Input::ApplysPHENIXBeamParameter(INPUTGENERATOR::Pythia8);
0200   }
0201
0202   //--------------
0203   // Set Input Manager specific options
0204   //--------------
0205   // can only be set after InputInit() is called
0206
0207   if (Input::HEPMC)
0208   {
0209     //! apply sPHENIX nominal beam parameter with 2mrad crossing as defined in sPH-TRG-2020-001
0210     Input::ApplysPHENIXBeamParameter(INPUTMANAGER::HepMCInputManager);
0211
0212     // optional overriding beam parameters
0213     //INPUTMANAGER::HepMCInputManager->set_vertex_distribution_width(100e-4, 100e-4, 8, 0);  //optional collision smear in space, time
0214     //    INPUTMANAGER::HepMCInputManager->set_vertex_distribution_mean(0,0,0,0);//optional collision central position shift in space, time
0215     // //optional choice of vertex distribution function in space, time
0216     //INPUTMANAGER::HepMCInputManager->set_vertex_distribution_function(PHHepMCGenHelper::Gaus, PHHepMCGenHelper::Gaus, PHHepMCGenHelper::Gaus, PHHepMCGenHelper::Gaus);
0217     //! embedding ID for the event
0218     //! positive ID is the embedded event of interest, e.g. jetty event from pythia
0219     //! negative IDs are backgrounds, .e.g out of time pile up collisions
0220     //! Usually, ID = 0 means the primary Au+Au collision background
0221     //INPUTMANAGER::HepMCInputManager->set_embedding_id(Input::EmbedID);
0222     if (Input::PILEUPRATE > 0)
0223     {
0224       // Copy vertex settings from foreground hepmc input
0225       INPUTMANAGER::HepMCPileupInputManager->CopyHelperSettings(INPUTMANAGER::HepMCInputManager);
0226       // and then modify the ones you want to be different
0227       // INPUTMANAGER::HepMCPileupInputManager->set_vertex_distribution_width(100e-4,100e-4,8,0);
0228     }
0229   }
0230   if (Input::PILEUPRATE > 0)
0231   {
0232     //! apply sPHENIX nominal beam parameter with 2mrad crossing as defined in sPH-TRG-2020-001
0233     Input::ApplysPHENIXBeamParameter(INPUTMANAGER::HepMCPileupInputManager);
0234   }
0235   // register all input generators with Fun4All
0236   InputRegister();
0237
0238   // set up production relatedstuff
0239   //   Enable::PRODUCTION = true;
0240
0241   //======================
0242   // Write the DST
0243   //======================
0244
0245   //Enable::DSTOUT = true;
0246   Enable::DSTOUT_COMPRESS = false;
0247   DstOut::OutputDir = outdir;
0248   DstOut::OutputFile = outputFile;
0249
0250   //Option to convert DST to human command readable TTree for quick poke around the outputs
0251   //  Enable::DSTREADER = true;
0252
0253   // turn the display on (default off)
0254    //Enable::DISPLAY = true;
0255
0256   //======================
0257   // What to run
0258   //======================
0259
0260   // QA, main switch
0261   Enable::QA = false;
0262
0263   // Global options (enabled for all enables subsystems - if implemented)
0264   //  Enable::ABSORBER = true;
0265   //  Enable::OVERLAPCHECK = true;
0266   //  Enable::VERBOSITY = 1;
0267
0268   // Enable::BBC = true;
0269   // Enable::BBC_SUPPORT = true; // save hist in bbc support structure
0270   Enable::BBCFAKE = true;  // Smeared vtx and t0, use if you don't want real BBC in simulation
0271
0272   Enable::PIPE = true;
0273   Enable::PIPE_ABSORBER = true;
0274
0275   // central tracking
0276   Enable::MVTX = false;
0277   Enable::MVTX_CELL = Enable::MVTX && true;
0278   Enable::MVTX_CLUSTER = Enable::MVTX_CELL && true;
0279   Enable::MVTX_QA = Enable::MVTX_CLUSTER && Enable::QA && true;
0280   Enable::TrackingService = false;
0281
0282   Enable::INTT = false;
0283 //  Enable::INTT_ABSORBER = true; // enables layerwise support structure readout
0284 //  Enable::INTT_SUPPORT = true; // enable global support structure readout
0285   Enable::INTT_CELL = Enable::INTT && true;
0286   Enable::INTT_CLUSTER = Enable::INTT_CELL && true;
0287   Enable::INTT_QA = Enable::INTT_CLUSTER && Enable::QA && true;
0288
0289   Enable::TPC = false;
0290   Enable::TPC_ABSORBER = true;
0291   Enable::TPC_CELL = Enable::TPC && true;
0292   Enable::TPC_CLUSTER = Enable::TPC_CELL && true;
0293   Enable::TPC_QA = Enable::TPC_CLUSTER && Enable::QA && true;
0294
0295   Enable::MICROMEGAS = false;
0296   Enable::MICROMEGAS_CELL = Enable::MICROMEGAS && true;
0297   Enable::MICROMEGAS_CLUSTER = Enable::MICROMEGAS_CELL && true;
0298   Enable::MICROMEGAS_QA = Enable::MICROMEGAS_CLUSTER && Enable::QA && true;
0299
0300   Enable::TRACKING_TRACK = false;
0301   Enable::TRACKING_EVAL = Enable::TRACKING_TRACK && true;
0302   Enable::TRACKING_QA = Enable::TRACKING_TRACK && Enable::QA && true;
0303
0304   //  cemc electronics + thin layer of W-epoxy to get albedo from cemc
0305   //  into the tracking, cannot run together with CEMC
0306   //  Enable::CEMCALBEDO = true;
0307
0308   Enable::CEMC = true;
0309   Enable::CEMC_ABSORBER = true;
0310   Enable::CEMC_CELL = Enable::CEMC && true;
0311   Enable::CEMC_TOWER = Enable::CEMC_CELL && true;
0312   Enable::CEMC_CLUSTER = Enable::CEMC_TOWER && true;
0313   Enable::CEMC_EVAL = false;//Enable::CEMC_CLUSTER && true;
0314   Enable::CEMC_QA = false;//Enable::CEMC_CLUSTER && Enable::QA && true;
0315
0316   Enable::HCALIN =false;
0317   Enable::HCALIN_ABSORBER = true;
0318   Enable::HCALIN_CELL = Enable::HCALIN && true;
0319   Enable::HCALIN_TOWER = Enable::HCALIN_CELL && true;
0320   Enable::HCALIN_CLUSTER = Enable::HCALIN_TOWER && true;
0321   Enable::HCALIN_EVAL = Enable::HCALIN_CLUSTER && true;
0322   Enable::HCALIN_QA = Enable::HCALIN_CLUSTER && Enable::QA && true;
0323
0324   Enable::MAGNET = false;
0325   Enable::MAGNET_ABSORBER = false;
0326
0327   Enable::HCALOUT = false;
0328   Enable::HCALOUT_ABSORBER = true;
0329   Enable::HCALOUT_CELL = Enable::HCALOUT && true;
0330   Enable::HCALOUT_TOWER = Enable::HCALOUT_CELL && true;
0331   Enable::HCALOUT_CLUSTER = Enable::HCALOUT_TOWER && true;
0332   Enable::HCALOUT_EVAL = Enable::HCALOUT_CLUSTER && true;
0333   Enable::HCALOUT_QA = Enable::HCALOUT_CLUSTER && Enable::QA && true;
0334
0335   Enable::EPD = false;
0336
0337   Enable::BEAMLINE = true;
0338 //  Enable::BEAMLINE_ABSORBER = true;  // makes the beam line magnets sensitive volumes
0339 //  Enable::BEAMLINE_BLACKHOLE = true; // turns the beamline magnets into black holes
0340   Enable::ZDC = false;
0341 //  Enable::ZDC_ABSORBER = true;
0342 //  Enable::ZDC_SUPPORT = true;
0343   Enable::ZDC_TOWER = Enable::ZDC && true;
0344   Enable::ZDC_EVAL = Enable::ZDC_TOWER && true;
0345
0346   //! forward flux return plug door. Out of acceptance and off by default.
0347   //Enable::PLUGDOOR = true;
0348   Enable::PLUGDOOR_ABSORBER = true;
0349
0350   Enable::GLOBAL_RECO = true;
0351   //Enable::GLOBAL_FASTSIM = true;
0352
0353   //Enable::KFPARTICLE = true;
0354   //Enable::KFPARTICLE_VERBOSITY = 1;
0355   //Enable::KFPARTICLE_TRUTH_MATCH = true;
0356   //Enable::KFPARTICLE_SAVE_NTUPLE = true;
0357
0358   Enable::CALOTRIGGER = Enable::CEMC_TOWER && Enable::HCALIN_TOWER && Enable::HCALOUT_TOWER && false;
0359
0360   Enable::JETS = false;
0361   Enable::JETS_EVAL = Enable::JETS && true;
0362   Enable::JETS_QA = Enable::JETS && Enable::QA && true;
0363
0364   // HI Jet Reco for p+Au / Au+Au collisions (default is false for
0365   // single particle / p+p-only simulations, or for p+Au / Au+Au
0366   // simulations which don't particularly care about jets)
0367   Enable::HIJETS = false && Enable::JETS && Enable::CEMC_TOWER && Enable::HCALIN_TOWER && Enable::HCALOUT_TOWER;
0368
0369   // 3-D topoCluster reconstruction, potentially in all calorimeter layers
0370   Enable::TOPOCLUSTER = false && Enable::CEMC_TOWER && Enable::HCALIN_TOWER && Enable::HCALOUT_TOWER;
0371   // particle flow jet reconstruction - needs topoClusters!
0372   Enable::PARTICLEFLOW = true && Enable::TOPOCLUSTER;
0373   // centrality reconstruction
0374   Enable::CENTRALITY = true;
0375
0376   // new settings using Enable namespace in GlobalVariables.C
0377   Enable::BLACKHOLE = true;
0378   //Enable::BLACKHOLE_SAVEHITS = false; // turn off saving of bh hits
0379   //BlackHoleGeometry::visible = true;
0380
0381   // run user provided code (from local G4_User.C)
0382   //Enable::USER = true;
0383
0384   //---------------
0385   // World Settings
0386   //---------------
0387   //  G4WORLD::PhysicsList = "FTFP_BERT"; //FTFP_BERT_HP best for calo
0388   //  G4WORLD::WorldMaterial = "G4_AIR"; // set to G4_GALACTIC for material scans
0389
0390   //---------------
0391   // Magnet Settings
0392   //---------------
0393
0394   //  G4MAGNET::magfield =  string(getenv("CALIBRATIONROOT"))+ string("/Field/Map/sphenix3dbigmapxyz.root");  // default map from the calibration database
0395   //  G4MAGNET::magfield = "1.5"; // alternatively to specify a constant magnetic field, give a float number, which will be translated to solenoidal field in T, if string use as fieldmap name (including path)
0396 //  G4MAGNET::magfield_rescale = 1.;  // make consistent with expected Babar field strength of 1.4T
0397
0398   //---------------
0399   // Pythia Decayer
0400   //---------------
0401   // list of decay types in
0402   // $OFFLINE_MAIN/include/g4decayer/EDecayType.hh
0403   // default is All:
0404   // G4P6DECAYER::decayType = EDecayType::kAll;
0405
0406   // Initialize the selected subsystems
0407   G4Init();
0408
0409   //---------------------
0410   // GEANT4 Detector description
0411   //---------------------
0412   if (!Input::READHITS)
0413   {
0414     G4Setup();
0415   }
0416
0417   //------------------
0418   // Detector Division
0419   //------------------
0420
0421   if (Enable::BBC || Enable::BBCFAKE) Bbc_Reco();
0422
0423   if (Enable::MVTX_CELL) Mvtx_Cells();
0424   if (Enable::INTT_CELL) Intt_Cells();
0425   if (Enable::TPC_CELL) TPC_Cells();
0426   if (Enable::MICROMEGAS_CELL) Micromegas_Cells();
0427
0428   if (Enable::CEMC_CELL) CEMC_Cells();
0429
0430   if (Enable::HCALIN_CELL) HCALInner_Cells();
0431
0432   if (Enable::HCALOUT_CELL) HCALOuter_Cells();
0433
0434   //-----------------------------
0435   // CEMC towering and clustering
0436   //-----------------------------
0437
0438   if (Enable::CEMC_TOWER) CEMC_Towers();
0439   if (Enable::CEMC_CLUSTER) CEMC_Clusters();
0440
0441   //-----------------------------
0442   // HCAL towering and clustering
0443   //-----------------------------
0444
0445   if (Enable::HCALIN_TOWER) HCALInner_Towers();
0446   if (Enable::HCALIN_CLUSTER) HCALInner_Clusters();
0447
0448   if (Enable::HCALOUT_TOWER) HCALOuter_Towers();
0449   if (Enable::HCALOUT_CLUSTER) HCALOuter_Clusters();
0450
0451   // if enabled, do topoClustering early, upstream of any possible jet reconstruction
0452   if (Enable::TOPOCLUSTER) TopoClusterReco();
0453
0454   //--------------
0455   // SVTX tracking
0456   //--------------
0457   if(Enable::TRACKING_TRACK)
0458     {
0459       TrackingInit();
0460     }
0461   if (Enable::MVTX_CLUSTER) Mvtx_Clustering();
0462   if (Enable::INTT_CLUSTER) Intt_Clustering();
0463   if (Enable::TPC_CLUSTER) TPC_Clustering();
0464   if (Enable::MICROMEGAS_CLUSTER) Micromegas_Clustering();
0465
0466   if (Enable::TRACKING_TRACK)
0467   {
0468     Tracking_Reco();
0469   }
0470   //-----------------
0471   // Global Vertexing
0472   //-----------------
0473
0474   if (Enable::GLOBAL_RECO && Enable::GLOBAL_FASTSIM)
0475   {
0476     cout << "You can only enable Enable::GLOBAL_RECO or Enable::GLOBAL_FASTSIM, not both" << endl;
0477     gSystem->Exit(1);
0478   }
0479   if (Enable::GLOBAL_RECO)
0480   {
0481     Global_Reco();
0482   }
0483   else if (Enable::GLOBAL_FASTSIM)
0484   {
0485     Global_FastSim();
0486   }
0487
0488   //-----------------
0489   // Centrality Determination
0490   //-----------------
0491
0492   if (Enable::CENTRALITY)
0493   {
0494       Centrality();
0495   }
0496
0497   //-----------------
0498   // Calo Trigger Simulation
0499   //-----------------
0500
0501   if (Enable::CALOTRIGGER)
0502   {
0503     CaloTrigger_Sim();
0504   }
0505
0506   //---------
0507   // Jet reco
0508   //---------
0509
0510   if (Enable::JETS) Jet_Reco();
0511   if (Enable::HIJETS) HIJetReco();
0512
0513   if (Enable::PARTICLEFLOW) ParticleFlow();
0514
0515   //----------------------
0516   // Simulation evaluation
0517   //----------------------
0518   string outputroot = outputFile;
0519   string remove_this = ".root";
0520   size_t pos = outputroot.find(remove_this);
0521   if (pos != string::npos)
0522   {
0523     outputroot.erase(pos, remove_this.length());
0524   }
0525
0526   if (Enable::TRACKING_EVAL) Tracking_Eval(outputroot + "_g4svtx_eval.root");
0527
0528   if (Enable::CEMC_EVAL) CEMC_Eval(outputroot + "_g4cemc_eval.root");
0529
0530   if (Enable::HCALIN_EVAL) HCALInner_Eval(outputroot + "_g4hcalin_eval.root");
0531
0532   if (Enable::HCALOUT_EVAL) HCALOuter_Eval(outputroot + "_g4hcalout_eval.root");
0533
0534   if (Enable::JETS_EVAL) Jet_Eval(outputroot + "_g4jet_eval.root");
0535
0536   if (Enable::DSTREADER) G4DSTreader(outputroot + "_DSTReader.root");
0537
0538   if (Enable::USER) UserAnalysisInit();
0539
0540
0541
0542
0543
0544   pi0Efficiency *eval = new pi0Efficiency("dummy", outputroot + "_pi0Efficiency.root");
0545   se->registerSubsystem(eval);
0546
0547
0548   //======================
0549   // Run KFParticle on evt
0550   //======================
0551   if (Enable::KFPARTICLE && Input::UPSILON) KFParticle_Upsilon_Reco();
0552   if (Enable::KFPARTICLE && Input::DZERO) KFParticle_D0_Reco();
0553
0554   //----------------------
0555   // Standard QAs
0556   //----------------------
0557
0558   if (Enable::CEMC_QA) CEMC_QA();
0559   if (Enable::HCALIN_QA) HCALInner_QA();
0560   if (Enable::HCALOUT_QA) HCALOuter_QA();
0561
0562   if (Enable::JETS_QA) Jet_QA();
0563
0564   if (Enable::MVTX_QA) Mvtx_QA();
0565   if (Enable::INTT_QA) Intt_QA();
0566   if (Enable::TPC_QA) TPC_QA();
0567   if (Enable::MICROMEGAS_QA) Micromegas_QA();
0568   if (Enable::TRACKING_QA) Tracking_QA();
0569
0570   if (Enable::TRACKING_QA && Enable::CEMC_QA && Enable::HCALIN_QA && Enable::HCALOUT_QA) QA_G4CaloTracking();
0571
0572   //--------------
0573   // Set up Input Managers
0574   //--------------
0575
0576   InputManagers();
0577
0578   if (Enable::PRODUCTION)
0579   {
0580     Production_CreateOutputDir();
0581   }
0582
0583   if (Enable::DSTOUT)
0584   {
0585     string FullOutFile = DstOut::OutputDir + "/" + DstOut::OutputFile;
0586     Fun4AllDstOutputManager *out = new Fun4AllDstOutputManager("DSTOUT", FullOutFile);
0587     if (Enable::DSTOUT_COMPRESS)
0588     {
0589       ShowerCompress();
0590       DstCompress(out);
0591     }
0592     se->registerOutputManager(out);
0593   }
0594   //-----------------
0595   // Event processing
0596   //-----------------
0597   if (Enable::DISPLAY)
0598   {
0599     DisplayOn();
0600
0601     gROOT->ProcessLine("Fun4AllServer *se = Fun4AllServer::instance();");
0602     gROOT->ProcessLine("PHG4Reco *g4 = (PHG4Reco *) se->getSubsysReco(\"PHG4RECO\");");
0603
0604     cout << "-------------------------------------------------" << endl;
0605     cout << "You are in event display mode. Run one event with" << endl;
0606     cout << "se->run(1)" << endl;
0607     cout << "Run Geant4 command with following examples" << endl;
0608     gROOT->ProcessLine("displaycmd()");
0609
0610     return 0;
0611   }
0612
0613   // if we use a negative number of events we go back to the command line here
0614   if (nEvents < 0)
0615   {
0616     return 0;
0617   }
0618   // if we run the particle generator and use 0 it'll run forever
0619   // for embedding it runs forever if the repeat flag is set
0620   if (nEvents == 0 && !Input::HEPMC && !Input::READHITS && INPUTEMBED::REPEAT)
0621   {
0622     cout << "using 0 for number of events is a bad idea when using particle generators" << endl;
0623     cout << "it will run forever, so I just return without running anything" << endl;
0624     return 0;
0625   }
0626
0627   se->skip(skip);
0628   se->run(nEvents);
0629
0630   //-----
0631   // QA output
0632   //-----
0633
0634   if (Enable::QA) QA_Output(outputroot + "_qa.root");
0635
0636   //-----
0637   // Exit
0638   //-----
0639
0640   se->End();
0641   std::cout << "All done" << std::endl;
0642   delete se;
0643   if (Enable::PRODUCTION)
0644   {
0645     Production_MoveOutput();
0646   }
0647
0648   gSystem->Exit(0);
0649   return 0;
0650 }
0651 #endif