Skip to content

Commit dbfcd3f

Browse files
sarahherrmannsarahherrmann
committed
bug corrected and new way to access the collisions
1 parent 3bfee44 commit dbfcd3f

1 file changed

Lines changed: 130 additions & 74 deletions

File tree

‎PWGMM/Mult/Tasks/vertexing-fwd.cxx‎

Lines changed: 130 additions & 74 deletions
Original file line numberDiff line numberDiff line change
@@ -55,7 +55,7 @@ struct vertexingfwd {
5555
{"TracksDCAX", "; DCA_{x} (cm); counts", {HistType::kTH1F, {{100, -10, 10}}}},
5656
{"TracksDCAY", "; DCA_{y} (cm); counts", {HistType::kTH1F, {{100, -10, 10}}}},
5757
{"AmbiguousTracksStatus", "; Status; counts", {HistType::kTH1F, {{6, -0.5, 5.5}}}},
58-
{"CollisionsSize", "; Collisions size", {HistType::kTH1F, {{10, -0.5, 9.5}}}},
58+
{"NbCollComp", "; NbCollComp", {HistType::kTH1F, {{10, -0.5, 9.5}}}},
5959
{"NumberOfContributors", "; N_{tr} for vertexing; counts", {HistType::kTH1F, {{100, 0, 100}}}},
6060
{"CollisionsMatchIndicesMC", "; Rec. minDCA ambitrack coll.ID; Gen. ambitrack coll.ID", {HistType::kTH2F, {{401, -0.5, 1000.5}, {401, -0.5, 1000.5}}}},
6161
{"TracksDCAXYBest", "; DCA_{xy}^{best} (cm); counts", {HistType::kTH1F, {{100, -1, 10}}}},
@@ -89,14 +89,24 @@ struct vertexingfwd {
8989
int getIndexBestCollision(std::vector<double> vecOfDCA, int method = 0)
9090
{
9191
int indice = 0;
92+
if (vecOfDCA.size() == 0) {
93+
return -1;
94+
}
9295
if (method == 0) {
9396
indice = std::distance(vecOfDCA.begin(), std::min_element(vecOfDCA.begin(), vecOfDCA.end()));
9497
}
9598
return indice;
9699
}
97100

98-
void processNew(aod::AmbiguousMFTTracks const& ambitracks, FullBCs const& bcs, soa::Join<o2::aod::MFTTracks, aod::McMFTTrackLabels> const& tracks, FullCollision const& collisions, aod::McParticles const& mcParticles, aod::McCollisions const& mcCollisions)
101+
void processNew(aod::AmbiguousMFTTracks const& ambitracks, aod::BCs const& bcs, soa::Join<o2::aod::MFTTracks, aod::McMFTTrackLabels> const& tracks, FullCollision const& collisions, aod::McParticles const& mcParticles, aod::McCollisions const& mcCollisions)
99102
{
103+
104+
int ntracks = tracks.size();
105+
int nambitracks = ambitracks.size();
106+
107+
registry.fill(HIST("AmbiguousTracksStatus"), 0.0, ntracks);
108+
registry.fill(HIST("AmbiguousTracksStatus"), 1.0, nambitracks);
109+
100110
for (auto& ambitrack : ambitracks) {
101111
vecCollForAmb.clear();
102112
vecDCACollForAmb.clear();
@@ -109,8 +119,6 @@ struct vertexingfwd {
109119

110120
auto track = ambitrack.mfttrack_as<MFTTracksLabeled>(); // Obtain the MFT ambiguous track with the MC labels
111121
// auto extAmbiTrackid = ambitrack.mfttrackId(); // Global index of the MFT ambiguous track
112-
int ntracks = tracks.size();
113-
int nambitracks = ambitracks.size();
114122

115123
if (!track.has_mcParticle()) {
116124
LOGF(warning, "No MC particle for ambiguous track, skip...");
@@ -120,8 +128,6 @@ struct vertexingfwd {
120128
mcCollAmbiID = particle.mcCollisionId();
121129
zVtxMCAmbi = particle.mcCollision().posZ();
122130

123-
registry.fill(HIST("AmbiguousTracksStatus"), 0.0, ntracks);
124-
registry.fill(HIST("AmbiguousTracksStatus"), 1.0, nambitracks);
125131
if (particle.isPhysicalPrimary()) {
126132
registry.fill(HIST("AmbiguousTracksStatus"), 2.0);
127133
} else {
@@ -137,51 +143,78 @@ struct vertexingfwd {
137143
vecAmbTrack.push_back(track.z());
138144
vecAmbTrack.push_back(track.chi2());
139145

140-
auto bcambis = ambitrack.bc_as<soa::Join<aod::BCs, aod::MatchedBCCollisionsSparse>>();
146+
auto bcambis = ambitrack.bc();
141147

142-
for (auto& bcambi : bcambis) {
148+
int collCounter = 0;
149+
for (auto& collision : collisions) {
150+
uint64_t mostProbableBC = collision.bc().globalBC();
151+
//uint64_t meanBC = mostProbableBC - std::lround(collision.collisionTime() / (o2::constants::lhc::LHCBunchSpacingNS / 1000));
152+
//int deltaBC = std::ceil(collision.collisionTimeRes() / (o2::constants::lhc::LHCBunchSpacingNS / 1000) * 4);
143153

144-
if (!bcambi.has_collision()) {
145-
continue;
146-
}
154+
for (auto& bcambi : bcambis) {
147155

148-
auto collision = bcambi.collision_as<soa::Join<aod::Collisions, aod::McCollisionLabels>>();
156+
if (bcambi.globalBC() != mostProbableBC) {
157+
continue;
158+
}
159+
collCounter++;
149160

150-
SMatrix5 tpars(vecAmbTrack[0], vecAmbTrack[1], vecAmbTrack[2], vecAmbTrack[3], vecAmbTrack[4]);
151-
// std::vector<double> v1{extAmbiTrack.cXX(), extAmbiTrack.cXY(), extAmbiTrack.cYY(), extAmbiTrack.cPhiX(), extAmbiTrack.cPhiY(),
152-
// extAmbiTrack.cPhiPhi(), extAmbiTrack.cTglX(), extAmbiTrack.cTglY(), extAmbiTrack.cTglPhi(), extAmbiTrack.cTglTgl(),
153-
// extAmbiTrack.c1PtX(), extAmbiTrack.c1PtY(), extAmbiTrack.c1PtPhi(), extAmbiTrack.c1PtTgl(), extAmbiTrack.c1Pt21Pt2()};
161+
SMatrix5 tpars(vecAmbTrack[0], vecAmbTrack[1], vecAmbTrack[2], vecAmbTrack[3], vecAmbTrack[4]);
162+
// std::vector<double> v1{extAmbiTrack.cXX(), extAmbiTrack.cXY(), extAmbiTrack.cYY(), extAmbiTrack.cPhiX(), extAmbiTrack.cPhiY(),
163+
// extAmbiTrack.cPhiPhi(), extAmbiTrack.cTglX(), extAmbiTrack.cTglY(), extAmbiTrack.cTglPhi(), extAmbiTrack.cTglTgl(),
164+
// extAmbiTrack.c1PtX(), extAmbiTrack.c1PtY(), extAmbiTrack.c1PtPhi(), extAmbiTrack.c1PtTgl(), extAmbiTrack.c1Pt21Pt2()};
154165

155-
std::vector<double> v1; // Temporary null vector for the computation of the covariance matrix
156-
SMatrix55 tcovs(v1.begin(), v1.end());
157-
o2::track::TrackParCovFwd pars1{vecAmbTrack[5], tpars, tcovs, vecAmbTrack[6]};
166+
std::vector<double> v1; // Temporary null vector for the computation of the covariance matrix
167+
SMatrix55 tcovs(v1.begin(), v1.end());
168+
o2::track::TrackParCovFwd pars1{vecAmbTrack[5], tpars, tcovs, vecAmbTrack[6]};
158169

159-
// o2::track::TrackParCovFwd pars1{extAmbiTrack.z(), tpars, tcovs, chi2};
160-
pars1.propagateToZlinear(collision.posZ()); // track parameters propagation to the position of the z vertex
170+
// o2::track::TrackParCovFwd pars1{extAmbiTrack.z(), tpars, tcovs, chi2};
171+
pars1.propagateToZlinear(collision.posZ()); // track parameters propagation to the position of the z vertex
161172

162-
const auto dcaX(pars1.getX() - collision.posX());
163-
const auto dcaY(pars1.getY() - collision.posY());
164-
auto dcaXY = std::sqrt(dcaX * dcaX + dcaY * dcaY);
173+
const auto dcaX(pars1.getX() - collision.posX());
174+
const auto dcaY(pars1.getY() - collision.posY());
175+
auto dcaXY = std::sqrt(dcaX * dcaX + dcaY * dcaY);
165176

166-
registry.fill(HIST("TracksDCAXY"), dcaXY);
167-
registry.fill(HIST("TracksDCAX"), dcaX);
168-
registry.fill(HIST("TracksDCAY"), dcaY);
169-
registry.fill(HIST("NumberOfContributors"), collision.numContrib());
177+
registry.fill(HIST("TracksDCAXY"), dcaXY);
178+
registry.fill(HIST("TracksDCAX"), dcaX);
179+
registry.fill(HIST("TracksDCAY"), dcaY);
180+
registry.fill(HIST("NumberOfContributors"), collision.numContrib());
170181

171-
int mcCollindex = collision.mcCollision().globalIndex();
172-
vecCollForAmb.push_back(mcCollindex);
173-
vecDCACollForAmb.push_back(dcaXY);
174-
vecZposCollForAmb.push_back(collision.mcCollision().posZ());
182+
vecDCACollForAmb.push_back(dcaXY);
175183

176-
registry.fill(HIST("DeltaZvtx"), collision.mcCollision().posZ() - zVtxMCAmbi);
177-
} // BCs loop
184+
if (!collision.has_mcCollision()) {
185+
continue;
186+
}
187+
188+
int mcCollindex = collision.mcCollision().globalIndex();
189+
vecCollForAmb.push_back(mcCollindex);
190+
191+
vecZposCollForAmb.push_back(collision.mcCollision().posZ());
192+
193+
registry.fill(HIST("DeltaZvtx"), collision.mcCollision().posZ() - zVtxMCAmbi);
194+
break;
195+
}
196+
}
197+
198+
registry.fill(HIST("NbCollComp"), collCounter);
199+
registry.fill(HIST("CorrectMatch"), 3.0); // counting for ambiguous track with N collisions >=0
178200

179201
int indexMinDCA = getIndexBestCollision(vecDCACollForAmb, 0); // obtain min value in the stored vector of DCAs
180-
int indexMCcoll = vecCollForAmb[indexMinDCA];
202+
int indexMCcoll = -1;
203+
if (indexMinDCA == -1) {
204+
continue;
205+
}
206+
indexMCcoll = vecCollForAmb[indexMinDCA];
207+
181208
registry.fill(HIST("CollisionsMatchIndicesMC"), mcCollAmbiID, indexMCcoll);
182-
registry.fill(HIST("CollisionsSize"), vecCollForAmb.size());
183-
registry.fill(HIST("CorrectMatch"), 3.0); // counting for amibuous track with N collisions >=0
184-
if (vecCollForAmb.size() == 0) { // do not use the vector with no collisions
209+
210+
if (collCounter == 1) {
211+
printf("strange ambiguous track of mfttrackId %d and trackGI %lld\n", ambitrack.mfttrackId(), track.globalIndex());
212+
if (mcCollAmbiID == indexMCcoll) {
213+
printf("and this is a correct match for the ambiguous track of mfttrackid %d\n", ambitrack.mfttrackId());
214+
}
215+
}
216+
217+
if (vecDCACollForAmb.size() == 0) { // do not use the vector with no collisions
185218
if (!particle.isPhysicalPrimary()) {
186219
registry.fill(HIST("AmbiguousTracksStatus"), 5.0);
187220
}
@@ -205,6 +238,7 @@ struct vertexingfwd {
205238
registry.fill(HIST("CorrectMatch"), value);
206239
registry.fill(HIST("EfficiencyZvtx"), zVtxMCAmbi, value);
207240
registry.fill(HIST("CorrectMatch"), 2.0); // Counting for amibuous track with N collisions > 0
241+
208242
if (value == 0.0) {
209243
registry.fill(HIST("TracksDCAXYBestFalse"), vecDCACollForAmb[indexMinDCA]); // Incorrect association with min DCA
210244
}
@@ -214,8 +248,15 @@ struct vertexingfwd {
214248

215249
PROCESS_SWITCH(vertexingfwd, processNew, "Process ambiguous track DCA", true);
216250

217-
void processOld(aod::AmbiguousTracksMFT const& ambitracks, FullBCs const& bcs, soa::Join<o2::aod::MFTTracks, aod::McMFTTrackLabels> const& tracks, FullCollision const& collisions, aod::McParticles const& mcParticles, aod::McCollisions const& mcCollisions)
251+
void processOld(aod::AmbiguousTracksMFT const& ambitracks, aod::BCs const& bcs, soa::Join<o2::aod::MFTTracks, aod::McMFTTrackLabels> const& tracks, FullCollision const& collisions, aod::McParticles const& mcParticles, aod::McCollisions const& mcCollisions)
218252
{
253+
254+
int ntracks = tracks.size();
255+
int nambitracks = ambitracks.size();
256+
257+
registry.fill(HIST("AmbiguousTracksStatus"), 0.0, ntracks);
258+
registry.fill(HIST("AmbiguousTracksStatus"), 1.0, nambitracks);
259+
219260
for (auto& ambitrack : ambitracks) {
220261
vecCollForAmb.clear();
221262
vecDCACollForAmb.clear();
@@ -228,8 +269,6 @@ struct vertexingfwd {
228269

229270
auto track = ambitrack.mfttrack_as<MFTTracksLabeled>(); // Obtain the MFT ambiguous track with the MC labels
230271
// auto extAmbiTrackid = ambitrack.mfttrackId(); // Global index of the MFT ambiguous track
231-
int ntracks = tracks.size();
232-
int nambitracks = ambitracks.size();
233272

234273
if (!track.has_mcParticle()) {
235274
LOGF(warning, "No MC particle for ambiguous track, skip...");
@@ -239,8 +278,6 @@ struct vertexingfwd {
239278
mcCollAmbiID = particle.mcCollisionId();
240279
zVtxMCAmbi = particle.mcCollision().posZ();
241280

242-
registry.fill(HIST("AmbiguousTracksStatus"), 0.0, ntracks);
243-
registry.fill(HIST("AmbiguousTracksStatus"), 1.0, nambitracks);
244281
if (particle.isPhysicalPrimary()) {
245282
registry.fill(HIST("AmbiguousTracksStatus"), 2.0);
246283
} else {
@@ -256,51 +293,70 @@ struct vertexingfwd {
256293
vecAmbTrack.push_back(track.z());
257294
vecAmbTrack.push_back(track.chi2());
258295

259-
auto bcambis = ambitrack.bc_as<soa::Join<aod::BCs, aod::MatchedBCCollisionsSparse>>();
296+
auto bcambis = ambitrack.bc();
260297

261-
for (auto& bcambi : bcambis) {
298+
int collCounter = 0;
299+
for (auto& collision : collisions) {
300+
uint64_t mostProbableBC = collision.bc().globalBC();
301+
//uint64_t meanBC = mostProbableBC - std::lround(collision.collisionTime() / (o2::constants::lhc::LHCBunchSpacingNS / 1000));
302+
//int deltaBC = std::ceil(collision.collisionTimeRes() / (o2::constants::lhc::LHCBunchSpacingNS / 1000) * 4);
262303

263-
if (!bcambi.has_collision()) {
264-
continue;
265-
}
304+
for (auto& bcambi : bcambis) {
305+
if (bcambi.globalBC() != mostProbableBC) {
306+
continue;
307+
}
308+
collCounter++;
309+
310+
SMatrix5 tpars(vecAmbTrack[0], vecAmbTrack[1], vecAmbTrack[2], vecAmbTrack[3], vecAmbTrack[4]);
311+
// std::vector<double> v1{extAmbiTrack.cXX(), extAmbiTrack.cXY(), extAmbiTrack.cYY(), extAmbiTrack.cPhiX(), extAmbiTrack.cPhiY(),
312+
// extAmbiTrack.cPhiPhi(), extAmbiTrack.cTglX(), extAmbiTrack.cTglY(), extAmbiTrack.cTglPhi(), extAmbiTrack.cTglTgl(),
313+
// extAmbiTrack.c1PtX(), extAmbiTrack.c1PtY(), extAmbiTrack.c1PtPhi(), extAmbiTrack.c1PtTgl(), extAmbiTrack.c1Pt21Pt2()};
266314

267-
auto collision = bcambi.collision_as<soa::Join<aod::Collisions, aod::McCollisionLabels>>();
315+
std::vector<double> v1; // Temporary null vector for the computation of the covariance matrix
316+
SMatrix55 tcovs(v1.begin(), v1.end());
317+
o2::track::TrackParCovFwd pars1{vecAmbTrack[5], tpars, tcovs, vecAmbTrack[6]};
268318

269-
SMatrix5 tpars(vecAmbTrack[0], vecAmbTrack[1], vecAmbTrack[2], vecAmbTrack[3], vecAmbTrack[4]);
270-
// std::vector<double> v1{extAmbiTrack.cXX(), extAmbiTrack.cXY(), extAmbiTrack.cYY(), extAmbiTrack.cPhiX(), extAmbiTrack.cPhiY(),
271-
// extAmbiTrack.cPhiPhi(), extAmbiTrack.cTglX(), extAmbiTrack.cTglY(), extAmbiTrack.cTglPhi(), extAmbiTrack.cTglTgl(),
272-
// extAmbiTrack.c1PtX(), extAmbiTrack.c1PtY(), extAmbiTrack.c1PtPhi(), extAmbiTrack.c1PtTgl(), extAmbiTrack.c1Pt21Pt2()};
319+
// o2::track::TrackParCovFwd pars1{extAmbiTrack.z(), tpars, tcovs, chi2};
320+
pars1.propagateToZlinear(collision.posZ()); // track parameters propagation to the position of the z vertex
273321

274-
std::vector<double> v1; // Temporary null vector for the computation of the covariance matrix
275-
SMatrix55 tcovs(v1.begin(), v1.end());
276-
o2::track::TrackParCovFwd pars1{vecAmbTrack[5], tpars, tcovs, vecAmbTrack[6]};
322+
const auto dcaX(pars1.getX() - collision.posX());
323+
const auto dcaY(pars1.getY() - collision.posY());
324+
auto dcaXY = std::sqrt(dcaX * dcaX + dcaY * dcaY);
277325

278-
// o2::track::TrackParCovFwd pars1{extAmbiTrack.z(), tpars, tcovs, chi2};
279-
pars1.propagateToZlinear(collision.posZ()); // track parameters propagation to the position of the z vertex
326+
registry.fill(HIST("TracksDCAXY"), dcaXY);
327+
registry.fill(HIST("TracksDCAX"), dcaX);
328+
registry.fill(HIST("TracksDCAY"), dcaY);
329+
registry.fill(HIST("NumberOfContributors"), collision.numContrib());
280330

281-
const auto dcaX(pars1.getX() - collision.posX());
282-
const auto dcaY(pars1.getY() - collision.posY());
283-
auto dcaXY = std::sqrt(dcaX * dcaX + dcaY * dcaY);
331+
vecDCACollForAmb.push_back(dcaXY);
284332

285-
registry.fill(HIST("TracksDCAXY"), dcaXY);
286-
registry.fill(HIST("TracksDCAX"), dcaX);
287-
registry.fill(HIST("TracksDCAY"), dcaY);
288-
registry.fill(HIST("NumberOfContributors"), collision.numContrib());
333+
if (!collision.has_mcCollision()) {
334+
continue;
335+
}
289336

290-
int mcCollindex = collision.mcCollision().globalIndex();
291-
vecCollForAmb.push_back(mcCollindex);
292-
vecDCACollForAmb.push_back(dcaXY);
293-
vecZposCollForAmb.push_back(collision.mcCollision().posZ());
337+
int mcCollindex = collision.mcCollision().globalIndex();
338+
vecCollForAmb.push_back(mcCollindex);
294339

295-
registry.fill(HIST("DeltaZvtx"), collision.mcCollision().posZ() - zVtxMCAmbi);
296-
} // BCs loop
340+
vecZposCollForAmb.push_back(collision.mcCollision().posZ());
341+
342+
registry.fill(HIST("DeltaZvtx"), collision.mcCollision().posZ() - zVtxMCAmbi);
343+
break;
344+
}
345+
}
346+
347+
registry.fill(HIST("NbCollComp"), collCounter);
297348

298349
int indexMinDCA = getIndexBestCollision(vecDCACollForAmb, 0); // obtain min value in the stored vector of DCAs
299-
int indexMCcoll = vecCollForAmb[indexMinDCA];
350+
int indexMCcoll = -1;
351+
if (indexMinDCA == -1) {
352+
continue;
353+
}
354+
indexMCcoll = vecCollForAmb[indexMinDCA];
355+
300356
registry.fill(HIST("CollisionsMatchIndicesMC"), mcCollAmbiID, indexMCcoll);
301-
registry.fill(HIST("CollisionsSize"), vecCollForAmb.size());
357+
302358
registry.fill(HIST("CorrectMatch"), 3.0); // counting for amibuous track with N collisions >=0
303-
if (vecCollForAmb.size() == 0) { // do not use the vector with no collisions
359+
if (vecDCACollForAmb.size() == 0) { // do not use the vector with no collisions
304360
if (!particle.isPhysicalPrimary()) {
305361
registry.fill(HIST("AmbiguousTracksStatus"), 5.0);
306362
}

0 commit comments

Comments
 (0)