pi0ClusterAna/macros/Fun4All_Pi0Gen_sPHENIX.C

0001 #ifndef MACRO_FUN4ALLG4SPHENIX_C
0002 #define MACRO_FUN4ALLG4SPHENIX_C
0003
0004 #include <../src/pi0ClusterAna.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(libpi0ClusterAna.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(-1.2, 1.2);
0161     INPUTGENERATOR::SimpleEventGenerator[0]->set_phi_range(-M_PI, M_PI);
0162     INPUTGENERATOR::SimpleEventGenerator[0]->set_pt_range(2, 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 = true;
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 = true;
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 = true;
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 = true;
0296   //Enable::MICROMEGAS_CELL = fasle;// = 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 = true;
0301   Enable::TRACKING_EVAL = false;//= Enable::TRACKING_TRACK && true;
0302   Enable::TRACKING_QA = false;//xs= 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   G4BEAMLINE::skin_thickness = 0.5;
0341
0342   Enable::ZDC = false;
0343 //  Enable::ZDC_ABSORBER = true;
0344 //  Enable::ZDC_SUPPORT = true;
0345   Enable::ZDC_TOWER = Enable::ZDC && true;
0346   Enable::ZDC_EVAL = Enable::ZDC_TOWER && true;
0347
0348   //! forward flux return plug door. Out of acceptance and off by default.
0349   //Enable::PLUGDOOR = true;
0350   //Enable::PLUGDOOR_ABSORBER = true;
0351   Enable::PLUGDOOR_BLACKHOLE = true;
0352
0353   Enable::GLOBAL_RECO = true;
0354   //Enable::GLOBAL_FASTSIM = true;
0355
0356   //Enable::KFPARTICLE = true;
0357   //Enable::KFPARTICLE_VERBOSITY = 1;
0358   //Enable::KFPARTICLE_TRUTH_MATCH = true;
0359   //Enable::KFPARTICLE_SAVE_NTUPLE = true;
0360
0361   Enable::CALOTRIGGER = Enable::CEMC_TOWER && Enable::HCALIN_TOWER && Enable::HCALOUT_TOWER && false;
0362
0363   Enable::JETS = false;
0364   Enable::JETS_EVAL = Enable::JETS && true;
0365   Enable::JETS_QA = Enable::JETS && Enable::QA && true;
0366
0367   // HI Jet Reco for p+Au / Au+Au collisions (default is false for
0368   // single particle / p+p-only simulations, or for p+Au / Au+Au
0369   // simulations which don't particularly care about jets)
0370   Enable::HIJETS = false && Enable::JETS && Enable::CEMC_TOWER && Enable::HCALIN_TOWER && Enable::HCALOUT_TOWER;
0371
0372   // 3-D topoCluster reconstruction, potentially in all calorimeter layers
0373   Enable::TOPOCLUSTER = false && Enable::CEMC_TOWER && Enable::HCALIN_TOWER && Enable::HCALOUT_TOWER;
0374   // particle flow jet reconstruction - needs topoClusters!
0375   Enable::PARTICLEFLOW = true && Enable::TOPOCLUSTER;
0376   // centrality reconstruction
0377   Enable::CENTRALITY = true;
0378
0379   // new settings using Enable namespace in GlobalVariables.C
0380   Enable::BLACKHOLE = true;
0381   //Enable::BLACKHOLE_SAVEHITS = false; // turn off saving of bh hits
0382   //BlackHoleGeometry::visible = true;
0383
0384   // run user provided code (from local G4_User.C)
0385   //Enable::USER = true;
0386
0387   //---------------
0388   // World Settings
0389   //---------------
0390   //  G4WORLD::PhysicsList = "FTFP_BERT"; //FTFP_BERT_HP best for calo
0391   //  G4WORLD::WorldMaterial = "G4_AIR"; // set to G4_GALACTIC for material scans
0392
0393   //---------------
0394   // Magnet Settings
0395   //---------------
0396
0397   //  G4MAGNET::magfield =  string(getenv("CALIBRATIONROOT"))+ string("/Field/Map/sphenix3dbigmapxyz.root");  // default map from the calibration database
0398   //  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)
0399 //  G4MAGNET::magfield_rescale = 1.;  // make consistent with expected Babar field strength of 1.4T
0400
0401   //---------------
0402   // Pythia Decayer
0403   //---------------
0404   // list of decay types in
0405   // $OFFLINE_MAIN/include/g4decayer/EDecayType.hh
0406   // default is All:
0407   // G4P6DECAYER::decayType = EDecayType::kAll;
0408
0409   // Initialize the selected subsystems
0410   G4Init();
0411
0412   //---------------------
0413   // GEANT4 Detector description
0414   //---------------------
0415   if (!Input::READHITS)
0416   {
0417     G4Setup();
0418   }
0419
0420   //------------------
0421   // Detector Division
0422   //------------------
0423
0424   if (Enable::BBC || Enable::BBCFAKE) Bbc_Reco();
0425
0426   if (Enable::MVTX_CELL) Mvtx_Cells();
0427   if (Enable::INTT_CELL) Intt_Cells();
0428   if (Enable::TPC_CELL) TPC_Cells();
0429   if (Enable::MICROMEGAS_CELL) Micromegas_Cells();
0430
0431   if (Enable::CEMC_CELL) CEMC_Cells();
0432
0433   if (Enable::HCALIN_CELL) HCALInner_Cells();
0434
0435   if (Enable::HCALOUT_CELL) HCALOuter_Cells();
0436
0437   //-----------------------------
0438   // CEMC towering and clustering
0439   //-----------------------------
0440
0441   if (Enable::CEMC_TOWER) CEMC_Towers();
0442   if (Enable::CEMC_CLUSTER) CEMC_Clusters();
0443
0444   //-----------------------------
0445   // HCAL towering and clustering
0446   //-----------------------------
0447
0448   if (Enable::HCALIN_TOWER) HCALInner_Towers();
0449   if (Enable::HCALIN_CLUSTER) HCALInner_Clusters();
0450
0451   if (Enable::HCALOUT_TOWER) HCALOuter_Towers();
0452   if (Enable::HCALOUT_CLUSTER) HCALOuter_Clusters();
0453
0454   // if enabled, do topoClustering early, upstream of any possible jet reconstruction
0455   if (Enable::TOPOCLUSTER) TopoClusterReco();
0456
0457   //--------------
0458   // SVTX tracking
0459   //--------------
0460   if(Enable::TRACKING_TRACK)
0461     {
0462       TrackingInit();
0463     }
0464   if (Enable::MVTX_CLUSTER) Mvtx_Clustering();
0465   if (Enable::INTT_CLUSTER) Intt_Clustering();
0466   if (Enable::TPC_CLUSTER) TPC_Clustering();
0467   if (Enable::MICROMEGAS_CLUSTER) Micromegas_Clustering();
0468
0469   if (Enable::TRACKING_TRACK)
0470   {
0471     Tracking_Reco();
0472   }
0473   //-----------------
0474   // Global Vertexing
0475   //-----------------
0476
0477   if (Enable::GLOBAL_RECO && Enable::GLOBAL_FASTSIM)
0478   {
0479     cout << "You can only enable Enable::GLOBAL_RECO or Enable::GLOBAL_FASTSIM, not both" << endl;
0480     gSystem->Exit(1);
0481   }
0482   if (Enable::GLOBAL_RECO)
0483   {
0484     Global_Reco();
0485   }
0486   else if (Enable::GLOBAL_FASTSIM)
0487   {
0488     Global_FastSim();
0489   }
0490
0491   //-----------------
0492   // Centrality Determination
0493   //-----------------
0494
0495   if (Enable::CENTRALITY)
0496   {
0497       Centrality();
0498   }
0499
0500   //-----------------
0501   // Calo Trigger Simulation
0502   //-----------------
0503
0504   if (Enable::CALOTRIGGER)
0505   {
0506     CaloTrigger_Sim();
0507   }
0508
0509   //---------
0510   // Jet reco
0511   //---------
0512
0513   if (Enable::JETS) Jet_Reco();
0514   if (Enable::HIJETS) HIJetReco();
0515
0516   if (Enable::PARTICLEFLOW) ParticleFlow();
0517
0518   //----------------------
0519   // Simulation evaluation
0520   //----------------------
0521   string outputroot = outputFile;
0522   string remove_this = ".root";
0523   size_t pos = outputroot.find(remove_this);
0524   if (pos != string::npos)
0525   {
0526     outputroot.erase(pos, remove_this.length());
0527   }
0528
0529   if (Enable::TRACKING_EVAL) Tracking_Eval(outputroot + "_g4svtx_eval.root");
0530
0531   if (Enable::CEMC_EVAL) CEMC_Eval(outputroot + "_g4cemc_eval.root");
0532
0533   if (Enable::HCALIN_EVAL) HCALInner_Eval(outputroot + "_g4hcalin_eval.root");
0534
0535   if (Enable::HCALOUT_EVAL) HCALOuter_Eval(outputroot + "_g4hcalout_eval.root");
0536
0537   if (Enable::JETS_EVAL) Jet_Eval(outputroot + "_g4jet_eval.root");
0538
0539   if (Enable::DSTREADER) G4DSTreader(outputroot + "_DSTReader.root");
0540
0541   if (Enable::USER) UserAnalysisInit();
0542
0543
0544
0545
0546
0547   pi0ClusterAna *eval = new pi0ClusterAna("dummy", outputroot + "_pi0ClusterAna.root");
0548   se->registerSubsystem(eval);
0549
0550
0551   //======================
0552   // Run KFParticle on evt
0553   //======================
0554   if (Enable::KFPARTICLE && Input::UPSILON) KFParticle_Upsilon_Reco();
0555   if (Enable::KFPARTICLE && Input::DZERO) KFParticle_D0_Reco();
0556
0557   //----------------------
0558   // Standard QAs
0559   //----------------------
0560
0561   if (Enable::CEMC_QA) CEMC_QA();
0562   if (Enable::HCALIN_QA) HCALInner_QA();
0563   if (Enable::HCALOUT_QA) HCALOuter_QA();
0564
0565   if (Enable::JETS_QA) Jet_QA();
0566
0567   if (Enable::MVTX_QA) Mvtx_QA();
0568   if (Enable::INTT_QA) Intt_QA();
0569   if (Enable::TPC_QA) TPC_QA();
0570   if (Enable::MICROMEGAS_QA) Micromegas_QA();
0571   if (Enable::TRACKING_QA) Tracking_QA();
0572
0573   if (Enable::TRACKING_QA && Enable::CEMC_QA && Enable::HCALIN_QA && Enable::HCALOUT_QA) QA_G4CaloTracking();
0574
0575   //--------------
0576   // Set up Input Managers
0577   //--------------
0578
0579   InputManagers();
0580
0581   if (Enable::PRODUCTION)
0582   {
0583     Production_CreateOutputDir();
0584   }
0585
0586   if (Enable::DSTOUT)
0587   {
0588     string FullOutFile = DstOut::OutputDir + "/" + DstOut::OutputFile;
0589     Fun4AllDstOutputManager *out = new Fun4AllDstOutputManager("DSTOUT", FullOutFile);
0590     if (Enable::DSTOUT_COMPRESS)
0591     {
0592       ShowerCompress();
0593       DstCompress(out);
0594     }
0595     se->registerOutputManager(out);
0596   }
0597   //-----------------
0598   // Event processing
0599   //-----------------
0600   if (Enable::DISPLAY)
0601   {
0602     DisplayOn();
0603
0604     gROOT->ProcessLine("Fun4AllServer *se = Fun4AllServer::instance();");
0605     gROOT->ProcessLine("PHG4Reco *g4 = (PHG4Reco *) se->getSubsysReco(\"PHG4RECO\");");
0606
0607     cout << "-------------------------------------------------" << endl;
0608     cout << "You are in event display mode. Run one event with" << endl;
0609     cout << "se->run(1)" << endl;
0610     cout << "Run Geant4 command with following examples" << endl;
0611     gROOT->ProcessLine("displaycmd()");
0612
0613     return 0;
0614   }
0615
0616   // if we use a negative number of events we go back to the command line here
0617   if (nEvents < 0)
0618   {
0619     return 0;
0620   }
0621   // if we run the particle generator and use 0 it'll run forever
0622   // for embedding it runs forever if the repeat flag is set
0623   if (nEvents == 0 && !Input::HEPMC && !Input::READHITS && INPUTEMBED::REPEAT)
0624   {
0625     cout << "using 0 for number of events is a bad idea when using particle generators" << endl;
0626     cout << "it will run forever, so I just return without running anything" << endl;
0627     return 0;
0628   }
0629
0630   se->skip(skip);
0631   se->run(nEvents);
0632
0633   //-----
0634   // QA output
0635   //-----
0636
0637   if (Enable::QA) QA_Output(outputroot + "_qa.root");
0638
0639   //-----
0640   // Exit
0641   //-----
0642
0643   se->End();
0644   std::cout << "All done" << std::endl;
0645   delete se;
0646   if (Enable::PRODUCTION)
0647   {
0648     Production_MoveOutput();
0649   }
0650
0651   gSystem->Exit(0);
0652   return 0;
0653 }
0654 #endif