macros/simulation/HF_selections.C

0001 #include <fun4all/Fun4AllServer.h>
0002
0003 #include <calotrigger/TriggerRunInfoReco.h>
0004 #include <globalvertex/GlobalVertexReco.h>
0005
0006 #include <kfparticleqa/QAKFParticle.h>
0007 #include <kfparticleqa/QAKFParticleTrackPtAsymmetry.h>
0008
0009 #pragma GCC diagnostic push
0010
0011 #pragma GCC diagnostic ignored "-Wundefined-internal"
0012
0013 #include <kfparticle_sphenix/KFParticle_sPHENIX.h>
0014
0015 #pragma GCC diagnostic pop
0016
0017 R__LOAD_LIBRARY(libkfparticle_sphenix.so)
0018 R__LOAD_LIBRARY(libcalotrigger.so)
0019 R__LOAD_LIBRARY(libkfparticleqa.so)
0020
0021 namespace HeavyFlavorReco
0022 {
0023   int VERBOSITY_HF = 0;
0024
0025   bool run_pipi_reco = true;
0026   bool run_ppi_reco = true; // set to true if needed
0027   bool run_KK_reco = false; // set to true if needed
0028   bool run_Kpi_reco = true; // set to true if needed
0029
0030   std::string output_dir = "./"; //Top dir of where the output nTuples will be written
0031   std::string kfp_header = "outputKFParticle_";
0032   std::string processing_folder = "inReconstruction/";
0033   std::string trailer = ".root";
0034
0035   // https://wiki.bnl.gov/sPHENIX/index.php/KFParticle
0036   std::string pipi_decay_descriptor = "K_S0 -> pi^+ pi^-"; //See twiki on how to set this
0037   std::string pipi_reconstruction_name = "pipi_reco"; //Used for naming output folder, file and node
0038   std::string pipi_output_reco_file;
0039   std::string pipi_output_dir;
0040
0041   std::string ppi_decay_descriptor = "[Lambda0 -> proton^+ pi^-]cc"; //See twiki on how to set this
0042   std::string ppi_reconstruction_name = "ppi_reco"; //Used for naming output folder, file and node
0043   std::string ppi_output_reco_file;
0044   std::string ppi_output_dir;
0045
0046   std::string KK_decay_descriptor = "phi -> K^+ K^-"; //See twiki on how to set this
0047   std::string KK_reconstruction_name = "KK_reco"; //Used for naming output folder, file and node
0048   std::string KK_output_reco_file;
0049   std::string KK_output_dir;
0050
0051   std::string Kpi_decay_descriptor = "[D0 -> K^- pi^+]cc"; //See twiki on how to set this
0052   std::string Kpi_reconstruction_name = "Kpi_reco"; //Used for naming output folder, file and node
0053   std::string Kpi_output_reco_file;
0054   std::string Kpi_output_dir;
0055
0056   bool save_kfpntuple = true;
0057   bool use_pid = false;
0058   bool save_tracks_to_DST = true;
0059   bool dont_use_global_vertex = true;
0060   bool require_track_and_vertex_match = true;
0061   bool save_all_vtx_info = true;
0062   bool constrain_phi_mass = true;
0063   bool use_2D_matching = false;
0064   bool get_trigger_info = false;
0065   bool get_detector_info = true;
0066   bool get_dEdx_info = true;
0067   float pid_frac = 0.2;
0068   bool constrain_lambda_mass = true;
0069 };  // namespace HeavyFlavorReco'
0070
0071 using namespace HeavyFlavorReco;
0072
0073 void init_kfp_dependencies()
0074 {
0075   //dE/dx needs TRKR_CLUSTER and CYLINDERCELLGEOM_SVTX which need to be in the DST or loaded from a geo file
0076   Fun4AllServer *se = Fun4AllServer::instance();
0077
0078   GlobalVertexReco* gblvertex = new GlobalVertexReco();
0079   gblvertex->Verbosity(VERBOSITY_HF);
0080   se->registerSubsystem(gblvertex);
0081
0082
0083   if (get_trigger_info)
0084   {
0085     TriggerRunInfoReco *triggerruninforeco = new TriggerRunInfoReco();
0086     se->registerSubsystem(triggerruninforeco);
0087   }
0088 }
0089
0090 void create_hf_directories(std::string reconstruction_name, std::string &final_output_dir, std::string &output_reco_file)
0091 {
0092   std::string output_file_name = kfp_header + reconstruction_name + trailer;
0093   final_output_dir = output_dir + reconstruction_name + "/";
0094   std::string output_reco_dir = final_output_dir + processing_folder;
0095   output_reco_file = output_reco_dir + output_file_name;
0096
0097   std::string makeDirectory = "mkdir -p " + output_reco_dir;
0098   system(makeDirectory.c_str());
0099 }
0100
0101 void reconstruct_pipi_mass()
0102 {
0103   Fun4AllServer *se = Fun4AllServer::instance();
0104
0105   KFParticle_sPHENIX *kfparticle = new KFParticle_sPHENIX(pipi_reconstruction_name);
0106   kfparticle->Verbosity(10);
0107
0108   kfparticle->setDecayDescriptor(pipi_decay_descriptor);
0109
0110   //kfparticle->extraolateTracksToSV(false); //To ensure the pT map is accurate
0111
0112   kfparticle->saveOutput(save_kfpntuple);
0113
0114   //kfparticle->useFakePrimaryVertex();
0115   //kfparticle->requireTrackVertexBunchCrossingMatch(false);
0116
0117   kfparticle->usePID(use_pid);
0118   kfparticle->useLocalPIDFile(true);
0119   kfparticle->setLocalPIDFilename("dedx_out.root");
0120   kfparticle->setPIDacceptFraction(pid_frac);
0121   kfparticle->get_dEdx_info();
0122   kfparticle->dontUseGlobalVertex(dont_use_global_vertex);
0123   //kfparticle->requireTrackVertexBunchCrossingMatch(require_track_and_vertex_match);
0124   kfparticle->getAllPVInfo(save_all_vtx_info);
0125   kfparticle->allowZeroMassTracks();
0126   kfparticle->use2Dmatching(use_2D_matching);
0127   kfparticle->getTriggerInfo(get_trigger_info);
0128   kfparticle->getDetectorInfo(get_detector_info);
0129   kfparticle->saveDST(save_tracks_to_DST);
0130   kfparticle->setContainerName(pipi_reconstruction_name);
0131   kfparticle->saveParticleContainer(true);
0132   kfparticle->magFieldFile("FIELDMAP_TRACKING");
0133
0134   //PV to SV cuts
0135   kfparticle->constrainToPrimaryVertex();
0136   //kfparticle->setMotherIPchi2(100);
0137   kfparticle->setMotherIPchi2(FLT_MAX);
0138   kfparticle->setFlightDistancechi2(-1.);
0139   kfparticle->setMinDIRA(0.99);
0140   kfparticle->setMinDIRA_XY(-1.1);
0141   kfparticle->setDecayLengthRange(0.05, FLT_MAX);
0142   kfparticle->setDecayLengthRange_XY(-10., FLT_MAX);
0143   kfparticle->setDecayTimeRange_XY(-10000, FLT_MAX);
0144   kfparticle->setDecayTimeRange(-10000, FLT_MAX);
0145   kfparticle->setMinDecayTimeSignificance(-1e5);
0146   kfparticle->setMinDecayLengthSignificance(-1e5);
0147   kfparticle->setMinDecayLengthSignificance_XY(-1e5);
0148
0149   //Track parameters
0150   kfparticle->setMinimumTrackPT(0.0);
0151   kfparticle->setMinimumTrackIPchi2(-1.);
0152   kfparticle->setMinimumTrackIPchi2_XY(-1.);
0153   kfparticle->setMinimumTrackIP(-1.);
0154   kfparticle->setMinimumTrackIP_XY(0.05);
0155   kfparticle->setMaximumTrackchi2nDOF(300.);
0156   kfparticle->setMinMVTXhits(1);
0157   kfparticle->setMinINTThits(1);
0158   kfparticle->setMinTPChits(20);
0159
0160   //Vertex parameters
0161   kfparticle->setMaximumVertexchi2nDOF(20);
0162   kfparticle->setMaximumDaughterDCA(0.1);
0163   kfparticle->setMaximumDaughterDCA_XY(1);
0164
0165   //Parent parameters
0166   kfparticle->setMotherPT(0);
0167   kfparticle->setMinimumMass(0.40);
0168   kfparticle->setMaximumMass(0.60);
0169   kfparticle->setMaximumMotherVertexVolume(0.1);
0170
0171   kfparticle->setOutputName(pipi_output_reco_file);
0172
0173   se->registerSubsystem(kfparticle);
0174
0175   QAKFParticle *kfpqa = new QAKFParticle("QAKFParticle_K_S0","K_S0",0.4,0.6);
0176   kfpqa->setKFParticleNodeName(pipi_reconstruction_name);
0177   kfpqa->enableTrackPtAsymmetry(true);
0178   kfpqa->Verbosity(VERBOSITY_HF);
0179   se->registerSubsystem(kfpqa);
0180
0181 }
0182
0183 void reconstruct_KK_mass()
0184 {
0185   Fun4AllServer *se = Fun4AllServer::instance();
0186
0187   KFParticle_sPHENIX *kfparticle = new KFParticle_sPHENIX(KK_reconstruction_name);
0188
0189   kfparticle->setDecayDescriptor(KK_decay_descriptor);
0190   kfparticle->saveOutput(save_kfpntuple);
0191
0192   kfparticle->usePID(use_pid);
0193   kfparticle->setPIDacceptFraction(pid_frac);
0194   kfparticle->useLocalPIDFile(true);
0195   kfparticle->setLocalPIDFilename("dedx_out.root");
0196   kfparticle->get_dEdx_info();
0197   kfparticle->dontUseGlobalVertex(dont_use_global_vertex);
0198   kfparticle->requireTrackVertexBunchCrossingMatch(require_track_and_vertex_match);
0199   kfparticle->getAllPVInfo(save_all_vtx_info);
0200   kfparticle->allowZeroMassTracks();
0201   kfparticle->use2Dmatching(use_2D_matching);
0202   kfparticle->getTriggerInfo(get_trigger_info);
0203   kfparticle->getDetectorInfo(get_detector_info);
0204   kfparticle->saveDST(save_tracks_to_DST);
0205   kfparticle->setContainerName(KK_reconstruction_name);
0206   kfparticle->saveParticleContainer(true);
0207   kfparticle->magFieldFile("FIELDMAP_TRACKING");
0208
0209   //PV to SV cuts
0210   kfparticle->constrainToPrimaryVertex();
0211   kfparticle->setMotherIPchi2(100);
0212   kfparticle->setFlightDistancechi2(-1.);
0213
0214   //Track parameters
0215   kfparticle->setMinimumTrackPT(0.1);
0216   kfparticle->setMaximumTrackPT(0.7);
0217   kfparticle->setMaximumTrackchi2nDOF(100.);
0218   kfparticle->setMinTPChits(25);
0219   kfparticle->setMinMVTXhits(1);
0220   kfparticle->setMinINTThits(0);
0221
0222   //Vertex parameters
0223   kfparticle->setMaximumVertexchi2nDOF(20);
0224   kfparticle->setMaximumDaughterDCA(0.05);
0225   kfparticle->setMaximumDaughterDCA_XY(100);
0226
0227   //Parent parameters
0228   kfparticle->setMotherPT(0);
0229   kfparticle->setMinimumMass(0.98);
0230   kfparticle->setMaximumMass(1.1);
0231   kfparticle->setMaximumMotherVertexVolume(0.1);
0232
0233   se->registerSubsystem(kfparticle);
0234
0235   QAKFParticle *kfpqa = new QAKFParticle("QAKFParticle_phi","phi",0.98,1.1);
0236   kfpqa->setKFParticleNodeName(KK_reconstruction_name);
0237   se->registerSubsystem(kfpqa);
0238 }
0239
0240 void reconstruct_ppi_mass()
0241 {
0242   Fun4AllServer *se = Fun4AllServer::instance();
0243
0244   KFParticle_sPHENIX *kfparticle = new KFParticle_sPHENIX(ppi_reconstruction_name);
0245
0246   kfparticle->setDecayDescriptor(ppi_decay_descriptor);
0247   kfparticle->saveOutput(save_kfpntuple);
0248
0249   kfparticle->usePID(use_pid);
0250   kfparticle->setPIDacceptFraction(pid_frac);
0251   kfparticle->useLocalPIDFile(true);
0252   kfparticle->setLocalPIDFilename("dedx_out.root");
0253   kfparticle->get_dEdx_info();
0254   kfparticle->dontUseGlobalVertex(dont_use_global_vertex);
0255   kfparticle->requireTrackVertexBunchCrossingMatch(require_track_and_vertex_match);
0256   kfparticle->getAllPVInfo(save_all_vtx_info);
0257   kfparticle->allowZeroMassTracks();
0258   kfparticle->use2Dmatching(use_2D_matching);
0259   kfparticle->getTriggerInfo(get_trigger_info);
0260   kfparticle->getDetectorInfo(get_detector_info);
0261   kfparticle->saveDST(save_tracks_to_DST);
0262   kfparticle->setContainerName(ppi_reconstruction_name);
0263   kfparticle->saveParticleContainer(true);
0264   kfparticle->magFieldFile("FIELDMAP_TRACKING");
0265
0266   //PV to SV cuts
0267   kfparticle->constrainToPrimaryVertex();
0268   //kfparticle->setMotherIPchi2(100);
0269   kfparticle->setMotherIPchi2(FLT_MAX);
0270   kfparticle->setFlightDistancechi2(-1.);
0271   kfparticle->setMinDIRA(0.99);
0272   kfparticle->setMinDIRA_XY(-1.1);
0273   kfparticle->setDecayLengthRange(0.05, FLT_MAX);
0274   kfparticle->setDecayLengthRange_XY(-10., FLT_MAX);
0275   kfparticle->setDecayTimeRange_XY(-10000, FLT_MAX);
0276   kfparticle->setDecayTimeRange(-10000, FLT_MAX);
0277   kfparticle->setMinDecayTimeSignificance(-1e5);
0278   kfparticle->setMinDecayLengthSignificance(-1e5);
0279   kfparticle->setMinDecayLengthSignificance_XY(-1e5);
0280
0281   //Track parameters
0282   kfparticle->setMinimumTrackPT(0.0);
0283   kfparticle->setMinimumTrackIPchi2(-1.);
0284   kfparticle->setMinimumTrackIPchi2_XY(-1.);
0285   kfparticle->setMinimumTrackIP(-1.);
0286   kfparticle->setMinimumTrackIP_XY(0.05);
0287   kfparticle->setMaximumTrackchi2nDOF(300.);
0288   kfparticle->setMinMVTXhits(1);
0289   kfparticle->setMinINTThits(1);
0290   kfparticle->setMinTPChits(20);
0291
0292   //Vertex parameters
0293   kfparticle->setMaximumVertexchi2nDOF(20);
0294   kfparticle->setMaximumDaughterDCA(0.5);
0295   kfparticle->setMaximumDaughterDCA_XY(100);
0296
0297   //Parent parameters
0298   kfparticle->setMotherPT(0);
0299   kfparticle->setMinimumMass(1.08);
0300   kfparticle->setMaximumMass(1.15);
0301   kfparticle->setMaximumMotherVertexVolume(0.1);
0302   kfparticle->setOutputName(ppi_output_reco_file);
0303
0304   se->registerSubsystem(kfparticle);
0305
0306   QAKFParticle *kfpqa = new QAKFParticle("QAKFParticle_Lambda0","Lambda0",1.08,1.15);
0307   kfpqa->setKFParticleNodeName(ppi_reconstruction_name);
0308   se->registerSubsystem(kfpqa);
0309 }
0310
0311 void reconstruct_Kpi_mass()
0312 {
0313   Fun4AllServer *se = Fun4AllServer::instance();
0314
0315   KFParticle_sPHENIX *kfparticle = new KFParticle_sPHENIX(Kpi_reconstruction_name);
0316   kfparticle->Verbosity(0);
0317
0318   kfparticle->setDecayDescriptor(Kpi_decay_descriptor);
0319   kfparticle->saveOutput(save_kfpntuple);
0320
0321   kfparticle->usePID(use_pid);
0322   kfparticle->setPIDacceptFraction(pid_frac);
0323   kfparticle->useLocalPIDFile(true);
0324   kfparticle->setLocalPIDFilename("dedx_out.root");
0325   kfparticle->get_dEdx_info();
0326   kfparticle->dontUseGlobalVertex(dont_use_global_vertex);
0327   kfparticle->requireTrackVertexBunchCrossingMatch(require_track_and_vertex_match);
0328   kfparticle->getAllPVInfo(save_all_vtx_info);
0329   kfparticle->allowZeroMassTracks();
0330   kfparticle->use2Dmatching(use_2D_matching);
0331   kfparticle->getTriggerInfo(get_trigger_info);
0332   kfparticle->getDetectorInfo(get_detector_info);
0333   kfparticle->saveDST(save_tracks_to_DST);
0334   kfparticle->setContainerName(Kpi_reconstruction_name);
0335   kfparticle->saveParticleContainer(true);
0336   kfparticle->magFieldFile("FIELDMAP_TRACKING");
0337
0338   //PV to SV cuts
0339   kfparticle->constrainToPrimaryVertex();
0340   //kfparticle->setMotherIPchi2(100);
0341   kfparticle->setMotherIPchi2(FLT_MAX);
0342   kfparticle->setFlightDistancechi2(-1.);
0343   kfparticle->setMinDIRA(0.95);
0344   kfparticle->setMinDIRA_XY(-1.1);
0345   kfparticle->setDecayLengthRange(0., FLT_MAX);
0346   kfparticle->setDecayLengthRange_XY(-10., FLT_MAX);
0347   kfparticle->setDecayTimeRange_XY(-10000, FLT_MAX);
0348   kfparticle->setDecayTimeRange(-10000, FLT_MAX);
0349   kfparticle->setMinDecayTimeSignificance(-1e5);
0350   kfparticle->setMinDecayLengthSignificance(-1e5);
0351   kfparticle->setMinDecayLengthSignificance_XY(-1e5);
0352
0353   //Track parameters
0354   kfparticle->setMinimumTrackPT(0.0);
0355   kfparticle->setMinimumTrackIPchi2(-1.);
0356   kfparticle->setMinimumTrackIPchi2_XY(-1.);
0357   kfparticle->setMinimumTrackIP(-1.);
0358   kfparticle->setMinimumTrackIP_XY(0.005);
0359   kfparticle->setMaximumTrackchi2nDOF(300.);
0360   kfparticle->setMinMVTXhits(1);
0361   kfparticle->setMinINTThits(1);
0362   kfparticle->setMinTPChits(20);
0363
0364   //Vertex parameters
0365   kfparticle->setMaximumVertexchi2nDOF(20);
0366   kfparticle->setMaximumDaughterDCA(0.5);
0367   kfparticle->setMaximumDaughterDCA_XY(100);
0368
0369   //Parent parameters
0370   kfparticle->setMotherPT(0);
0371   kfparticle->setMinimumMass(1.75);
0372   kfparticle->setMaximumMass(1.95);
0373   kfparticle->setMaximumMotherVertexVolume(0.1);
0374   kfparticle->setOutputName(Kpi_output_reco_file);
0375
0376   se->registerSubsystem(kfparticle);
0377
0378   QAKFParticle *kfpqa = new QAKFParticle("QAKFParticle_D0","D0",1.75,1.95);
0379   kfpqa->setKFParticleNodeName(Kpi_reconstruction_name);
0380   se->registerSubsystem(kfpqa);
0381 }
0382
0383 void end_kfparticle(std::string full_file_name, std::string final_path)
0384 {
0385   ifstream file(full_file_name.c_str());
0386   if (file.good())
0387   {
0388     string moveOutput = "mv " + full_file_name + " " + final_path;
0389     system(moveOutput.c_str());
0390   }
0391 }