File indexing completed on 2025-12-17 09:21:21
0001 #include "QAKFParticleTrackPtAsymmetry.h"
0002
0003 #include <fun4all/Fun4AllHistoManager.h>
0004 #include <fun4all/Fun4AllReturnCodes.h>
0005 #include <qautils/QAHistManagerDef.h> // for getHistoManager
0006
0007 #include <trackbase_historic/SvtxTrack.h>
0008 #include <trackbase_historic/SvtxTrackMap.h>
0009
0010 #include <KFParticle.h>
0011
0012 #include <TH1.h>
0013 #include <TH2.h>
0014
0015 #include <cassert>
0016 #include <iostream>
0017
0018 QAKFParticleTrackPtAsymmetry::QAKFParticleTrackPtAsymmetry(const std::string &histo_prefix,
0019 double min_m,
0020 double max_m,
0021 const std::vector<double> &mother_eta_bins,
0022 const std::vector<double> &mother_phi_bins)
0023 : m_prefix(histo_prefix)
0024 , m_min_mass(min_m)
0025 , m_max_mass(max_m)
0026 , m_mother_eta_bins(mother_eta_bins)
0027 , m_mother_phi_bins(mother_phi_bins)
0028 {
0029 if (m_mother_eta_bins.size() < 2)
0030 {
0031 std::cout << __PRETTY_FUNCTION__ << " Error: need at least 2 eta bin edges. Use default" << std::endl;
0032 m_mother_eta_bins = std::vector<double>{-2.0, -1.0, -0.5, 0, 0.5, 1.0, 2.0};
0033 }
0034
0035 if (m_mother_phi_bins.size() < 2)
0036 {
0037 std::cout << __PRETTY_FUNCTION__ << " Error: need at least 2 phi bin edges. Use default" << std::endl;
0038 m_mother_phi_bins = std::vector<double>{-3.15, -2.5, -2.0, -1.5, -1.0, -0.5, 0.0, 0.5, 1.0, 1.5, 2.0, 2.5, 3.15};
0039 }
0040 }
0041
0042 void QAKFParticleTrackPtAsymmetry::bookHistograms(Fun4AllHistoManager *hm)
0043 {
0044 assert(hm);
0045
0046
0047 h_trackPtAsymmetry = new TH1F(TString(m_prefix) + "TrackPtAsymmetry", ";#Delta p_{T}(trk^{+}, trk^{-})/#Sigma p_{T}(trk^{+}, trk^{-});Entries", 100, -1.0, 1.0);
0048 hm->registerHisto(h_trackPtAsymmetry);
0049
0050 h2_trackPtAsymmetry_vs_mass = new TH2F(TString(m_prefix) + "TrackPtAsymmetry_vs_Mass", ";#Delta p_{T}(trk^{+}, trk^{-})/#Sigma p_{T}(trk^{+}, trk^{-});Mass [GeV]", 100, -1.0, 1.0, 100, m_min_mass, m_max_mass);
0051 hm->registerHisto(h2_trackPtAsymmetry_vs_mass);
0052
0053
0054 for (size_t i = 0; i < m_mother_eta_bins.size() - 1; ++i)
0055 {
0056 std::vector<TH2 *> phi_bin_histos;
0057 for (size_t j = 0; j < m_mother_phi_bins.size() - 1; ++j)
0058 {
0059 std::string histo_name = m_prefix + "TrackPtAsymmetry_eta" + std::to_string(m_mother_eta_bins[i]) + "_" + std::to_string(m_mother_eta_bins[i + 1]) + "_phi" + std::to_string(m_mother_phi_bins[j]) + "_" + std::to_string(m_mother_phi_bins[j + 1]);
0060
0061 TH2 *h2 = new TH2F(histo_name.c_str(), ";#Delta p_{T}(trk^{+}, trk^{-})/#Sigma p_{T}(trk^{+}, trk^{-}) ;Mass [GeV]", 100, -1.0, 1.0, 100, m_min_mass, m_max_mass);
0062 hm->registerHisto(h2);
0063 phi_bin_histos.push_back(h2);
0064 }
0065 h2_trackPtAsymmetry_eta_phi_bins.push_back(phi_bin_histos);
0066 }
0067 }
0068
0069 void QAKFParticleTrackPtAsymmetry::analyzeTrackPtAsymmetry(SvtxTrackMap *m_trackMap, KFParticle *mother)
0070 {
0071 if (!m_trackMap)
0072 {
0073 if (m_verbosity > 0)
0074 {
0075 std::cout << __PRETTY_FUNCTION__ << " Error: null SvtxTrackMap pointer" << std::endl;
0076 }
0077 return;
0078 }
0079
0080 if (!mother)
0081 {
0082 if (m_verbosity > 0)
0083 {
0084 std::cout << __PRETTY_FUNCTION__ << " Error: null mother pointer" << std::endl;
0085 }
0086 return;
0087 }
0088
0089 if (mother->DaughterIds().size() != 2)
0090 {
0091 if (m_verbosity > 0)
0092 {
0093 std::cout << __PRETTY_FUNCTION__ << " Error: number of decay products is " << mother->DaughterIds().size() << ", not a 2-body decay" << std::endl;
0094 }
0095 return;
0096 }
0097
0098 const std::vector<int> &track_ids = mother->DaughterIds();
0099 SvtxTrack *trk1 = nullptr;
0100 SvtxTrack *trk2 = nullptr;
0101 auto it1 = m_trackMap->find(track_ids[0]);
0102 if (it1 != m_trackMap->end())
0103 {
0104 trk1 = it1->second;
0105 }
0106 auto it2 = m_trackMap->find(track_ids[1]);
0107 if (it2 != m_trackMap->end())
0108 {
0109 trk2 = it2->second;
0110 }
0111
0112 if (!trk1 || !trk2)
0113 {
0114 if (m_verbosity > 0)
0115 {
0116 std::cout << __PRETTY_FUNCTION__ << " Error: cannot to find one or both daughter tracks in SvtxTrackMap" << std::endl;
0117 }
0118 return;
0119 }
0120
0121
0122 SvtxTrack *pos_trk = nullptr;
0123 SvtxTrack *neg_trk = nullptr;
0124 if (trk1->get_charge() > 0 && trk2->get_charge() < 0)
0125 {
0126 pos_trk = trk1;
0127 neg_trk = trk2;
0128 }
0129 else if (trk1->get_charge() < 0 && trk2->get_charge() > 0)
0130 {
0131 pos_trk = trk2;
0132 neg_trk = trk1;
0133 }
0134 else
0135 {
0136 if (m_verbosity > 0)
0137 {
0138 std::cout << __PRETTY_FUNCTION__ << " Error: both daughter tracks have the same charge" << std::endl;
0139 }
0140 return;
0141 }
0142
0143 float pt_pos = pos_trk->get_pt();
0144 float pt_neg = neg_trk->get_pt();
0145 float pt_sum = pt_pos + pt_neg;
0146 if (pt_sum == 0)
0147 {
0148 if (m_verbosity > 0)
0149 {
0150 std::cout << __PRETTY_FUNCTION__ << " Warning: sum of daughter track pt is zero" << std::endl;
0151 }
0152 return;
0153 }
0154 float pt_asymmetry = (pt_pos - pt_neg) / pt_sum;
0155
0156 if (h_trackPtAsymmetry)
0157 {
0158 h_trackPtAsymmetry->Fill(pt_asymmetry);
0159 }
0160 if (h2_trackPtAsymmetry_vs_mass)
0161 {
0162 h2_trackPtAsymmetry_vs_mass->Fill(pt_asymmetry, mother->GetMass());
0163 }
0164
0165 float mother_eta = 0;
0166 float mother_eta_err = 0;
0167 mother->GetEta(mother_eta, mother_eta_err);
0168 float mother_phi = 0;
0169 float mother_phi_err = 0;
0170 mother->GetPhi(mother_phi, mother_phi_err);
0171 for (size_t i = 0; i < m_mother_eta_bins.size() - 1; ++i)
0172 {
0173 if (mother_eta >= m_mother_eta_bins[i] && mother_eta < m_mother_eta_bins[i + 1])
0174 {
0175 for (size_t j = 0; j < m_mother_phi_bins.size() - 1; ++j)
0176 {
0177 if (mother_phi >= m_mother_phi_bins[j] && mother_phi < m_mother_phi_bins[j + 1])
0178 {
0179 TH2 *h2 = h2_trackPtAsymmetry_eta_phi_bins[i][j];
0180 if (h2)
0181 {
0182 h2->Fill(pt_asymmetry, mother->GetMass());
0183 }
0184 }
0185 }
0186 }
0187 }
0188 }
