2 #include <qautils/QAUtil.h>
3 #include <qautils/QAHistManagerDef.h>
55 auto geom_container = findNode::getClass<PHG4CylinderGeomContainer>(topNode,
"CYLINDERGEOM_MVTX");
58 std::cout <<
PHWHERE <<
" unable to find DST node CYLINDERGEOM_MVTX" << std::endl;
63 const auto range = geom_container->get_begin_end();
64 for (
auto iter = range.first; iter != range.second; ++iter)
79 auto h =
new TH1F(Form(
"%sdrphi_%i",
get_histo_prefix().c_str(),
layer), Form(
"MVTX r#Delta#phi_{cluster-truth} layer_%i",
layer), 100, -2
e-3, 2
e-3);
80 h->GetXaxis()->SetTitle(
"r#Delta#phi_{cluster-truth} (cm)");
86 auto h =
new TH1F(Form(
"%srphi_error_%i",
get_histo_prefix().c_str(),
layer), Form(
"MVTX r#Delta#phi error layer_%i",
layer), 100, 0, 2
e-3);
87 h->GetXaxis()->SetTitle(
"r#Delta#phi error (cm)");
93 auto h =
new TH1F(Form(
"%sphi_pulls_%i",
get_histo_prefix().c_str(),
layer), Form(
"MVTX #Delta#phi_{cluster-truth}/#sigma#phi layer_%i",
layer), 100, -3, 3);
94 h->GetXaxis()->SetTitle(
"#Delta#phi_{cluster-truth}/#sigma#phi");
100 auto h =
new TH1F(Form(
"%sdz_%i",
get_histo_prefix().c_str(),
layer), Form(
"MVTX #Deltaz_{cluster-truth} layer_%i",
layer), 100, -3
e-3, 3
e-3);
101 h->GetXaxis()->SetTitle(
"#Deltaz_{cluster-truth} (cm)");
102 hm->registerHisto(
h);
108 h->GetXaxis()->SetTitle(
"z error (cm)");
109 hm->registerHisto(
h);
114 auto h =
new TH1F(Form(
"%sz_pulls_%i",
get_histo_prefix().c_str(),
layer), Form(
"MVTX #Deltaz_{cluster-truth}/#sigmaz layer_%i",
layer), 100, -3, 3);
115 h->GetXaxis()->SetTitle(
"#Deltaz_{cluster-truth}/#sigmaz");
116 hm->registerHisto(
h);
121 auto h =
new TH1F(Form(
"%sclus_size_%i",
get_histo_prefix().c_str(),
layer), Form(
"MVTX cluster size layer_%i",
layer), 20, 0, 20);
122 h->GetXaxis()->SetTitle(
"csize");
123 hm->registerHisto(
h);
128 auto h =
new TH1F(Form(
"%sclus_size_phi_%i",
get_histo_prefix().c_str(),
layer), Form(
"MVTX cluster size (#phi) layer_%i",
layer), 10, 0, 10);
129 h->GetXaxis()->SetTitle(
"csize_{#phi}");
130 hm->registerHisto(
h);
135 auto h =
new TH1F(Form(
"%sclus_size_z_%i",
get_histo_prefix().c_str(),
layer), Form(
"MVTX cluster size (z) layer_%i",
layer), 10, 0, 10);
136 h->GetXaxis()->SetTitle(
"csize_{z}");
137 hm->registerHisto(
h);
164 m_tGeometry = findNode::getClass<ActsGeometry>(topNode,
"ActsGeometry");
167 std::cout <<
PHWHERE <<
"No acts tracking geometry, exiting."
172 m_cluster_map = findNode::getClass<TrkrClusterContainer>(topNode,
"TRKR_CLUSTER");
175 std::cout <<
PHWHERE <<
" unable to find DST node TRKR_CLUSTER" << std::endl;
179 m_cluster_hit_map = findNode::getClass<TrkrClusterHitAssoc>(topNode,
"TRKR_CLUSTERHITASSOC");
182 std::cout <<
PHWHERE <<
" unable to find DST node TRKR_CLUSTERHITASSOC" << std::endl;
186 m_hit_truth_map = findNode::getClass<TrkrHitTruthAssoc>(topNode,
"TRKR_HITTRUTHASSOC");
189 std::cout <<
PHWHERE <<
" unable to find DST node TRKR_HITTRUTHASSOC" << std::endl;
193 m_g4hits_mvtx = findNode::getClass<PHG4HitContainer>(topNode,
"G4HIT_MVTX");
196 std::cout <<
PHWHERE <<
" unable to find DST node G4HIT_MVTX" << std::endl;
213 TH1* drphi =
nullptr;
214 TH1* rphi_error =
nullptr;
215 TH1* phi_pulls =
nullptr;
218 TH1* z_error =
nullptr;
219 TH1* z_pulls =
nullptr;
221 TH1* csize =
nullptr;
222 TH1* csize_phi =
nullptr;
223 TH1* csize_z =
nullptr;
226 using HistogramMap = std::map<int, HistogramList>;
227 HistogramMap histograms;
233 h.rphi_error =
dynamic_cast<TH1*
>(hm->getHisto(Form(
"%srphi_error_%i",
get_histo_prefix().c_str(),
layer)));
234 h.phi_pulls =
dynamic_cast<TH1*
>(hm->getHisto(Form(
"%sphi_pulls_%i",
get_histo_prefix().c_str(),
layer)));
237 h.z_error =
dynamic_cast<TH1*
>(hm->getHisto(Form(
"%sz_error_%i",
get_histo_prefix().c_str(),
layer)));
238 h.z_pulls =
dynamic_cast<TH1*
>(hm->getHisto(Form(
"%sz_pulls_%i",
get_histo_prefix().c_str(),
layer)));
240 h.csize =
dynamic_cast<TH1*
>(hm->getHisto(Form(
"%sclus_size_%i",
get_histo_prefix().c_str(),
layer)));
241 h.csize_phi =
dynamic_cast<TH1*
>(hm->getHisto(Form(
"%sclus_size_phi_%i",
get_histo_prefix().c_str(),
layer)));
242 h.csize_z =
dynamic_cast<TH1*
>(hm->getHisto(Form(
"%sclus_size_z_%i",
get_histo_prefix().c_str(),
layer)));
244 histograms.insert(std::make_pair(
layer, h));
250 for (
auto clusterIter = range.first; clusterIter != range.second; ++clusterIter)
253 const auto& key = clusterIter->first;
255 const auto& cluster = clusterIter->second;
260 const auto z_cluster = global(2);
261 const auto phi_cluster = (float) std::atan2(global(1), global(0));
263 double phi_error = cluster->getRPhiError() / r_cluster;
264 double z_error = cluster->getZError();
270 const auto x_truth = QAG4Util::interpolate<&PHG4Hit::get_x>(g4hits, r_cluster);
271 const auto y_truth = QAG4Util::interpolate<&PHG4Hit::get_y>(g4hits, r_cluster);
272 const auto z_truth = QAG4Util::interpolate<&PHG4Hit::get_z>(g4hits, r_cluster);
273 const auto phi_truth = std::atan2(y_truth, x_truth);
276 const auto dz = z_cluster - z_truth;
280 const auto hiter = histograms.find(
layer);
281 if (hiter == histograms.end())
continue;
284 auto fill = [](TH1*
h,
float value)
285 {
if( h ) h->Fill(
value ); };
286 fill(hiter->second.drphi, r_cluster *
dphi);
287 fill(hiter->second.rphi_error, r_cluster * phi_error);
288 fill(hiter->second.phi_pulls, dphi / phi_error);
291 fill(hiter->second.dz,
dz);
292 fill(hiter->second.z_error, z_error);
293 fill(hiter->second.z_pulls,
dz / z_error);
298 fill(hiter->second.csize,
std::distance(hit_range.first, hit_range.second));
302 for (
auto hit_iter = hit_range.first; hit_iter != hit_range.second; ++hit_iter)
304 const auto& hit_key = hit_iter->second;
309 fill(hiter->second.csize_phi, phibins.size());
310 fill(hiter->second.csize_z, zbins.size());
323 for (
auto iter = range.first; iter != range.second; ++iter)
326 const auto& hit_key = iter->second;
333 for (
auto truth_iter = g4hit_map.begin(); truth_iter != g4hit_map.end(); ++truth_iter)
336 const auto g4hit_key = truth_iter->second.second;
342 if (g4hit) out.insert(g4hit);