STrackCutStudy/macros/Fun4All_G4_sPHENIX_ForTrackCutStudy.C

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