2424#include < cmath>
2525#include < array>
2626#include < cstdlib>
27+ #include < chrono>
2728
2829#include " Framework/runDataProcessing.h"
2930#include " Framework/AnalysisTask.h"
@@ -73,6 +74,7 @@ struct sigma0builder {
7374 // For manual sliceBy
7475 PresliceUnsorted<V0DerivedMCDatas> perCollisionMCDerived = o2::aod::v0data::straCollisionId;
7576 PresliceUnsorted<V0StandardDerivedDatas> perCollisionSTDDerived = o2::aod::v0data::straCollisionId;
77+ Preslice<V0StandardDerivedDatas> perCollisionSTDSorted = o2::aod::v0data::straCollisionId;
7678 PresliceUnsorted<soa::Join<aod::StraCollisions, aod::StraCents, aod::StraEvSels, aod::StraCollLabels>> perMcCollision = aod::v0data::straMCCollisionId;
7779
7880 // pack track quality but separte also afterburner
@@ -90,6 +92,7 @@ struct sigma0builder {
9092 Configurable<bool > fillBkgQAhistos{" fillBkgQAhistos" , false , " if true, fill MC QA histograms for Bkg study" };
9193 Configurable<bool > doPi0QA{" doPi0QA" , true , " Flag to fill QA histos for pi0 rejection study." };
9294 Configurable<bool > doAssocStudy{" doAssocStudy" , false , " Do v0 to collision association study." };
95+ Configurable<bool > fverbose{" fverbose" , false , " QA printout." };
9396
9497 // Event level
9598 Configurable<bool > doPPAnalysis{" doPPAnalysis" , true , " if in pp, set to true" };
@@ -1317,13 +1320,20 @@ struct sigma0builder {
13171320
13181321 void processRealData (soa::Join<aod::StraCollisions, aod::StraCents, aod::StraEvSels, aod::StraStamps> const & collisions, V0StandardDerivedDatas const & fullV0s, dauTracks const &)
13191322 {
1323+ uint64_t CollIDBuffer = 0 ;
1324+ auto start = std::chrono::high_resolution_clock::now ();
13201325 for (const auto & coll : collisions) {
13211326
13221327 if (!IsEventAccepted (coll, true ))
13231328 continue ;
13241329
13251330 // Do analysis with collision-grouped V0s, retain full collision information
13261331 const uint64_t collIdx = coll.globalIndex ();
1332+ if (collIdx < CollIDBuffer)
1333+ LOGF (fatal, " Collision table unsorted! Previous index: %i, current index: %i" , CollIDBuffer, collIdx);
1334+
1335+ CollIDBuffer = collIdx;
1336+
13271337 auto V0s = fullV0s.sliceBy (perCollisionSTDDerived, collIdx);
13281338 float centrality = doPPAnalysis ? coll.centFT0M () : coll.centFT0C ();
13291339
@@ -1383,8 +1393,240 @@ struct sigma0builder {
13831393 if (nSigmaCandidates % 10000 == 0 )
13841394 LOG (info) << " Sigma0 Candidates built: " << nSigmaCandidates;
13851395 }
1396+ }
1397+ }
1398+ auto end = std::chrono::high_resolution_clock::now ();
1399+ std::chrono::duration<double > elapsed = end - start;
1400+
1401+ if (fverbose) LOGF (info, " [Process function call, PreSliceUnsorted] N. Collisions: %i, N. V0s: %i, Processing time (s): %lf" , collisions.size (), fullV0s.size (), elapsed.count ());
1402+ }
1403+
1404+ void processRealDataSorted (soa::Join<aod::StraCollisions, aod::StraCents, aod::StraEvSels, aod::StraStamps> const & collisions, V0StandardDerivedDatas const & fullV0s, dauTracks const &)
1405+ {
1406+ uint64_t CollIDBuffer = 0 ;
1407+ auto start = std::chrono::high_resolution_clock::now ();
1408+ for (const auto & coll : collisions) {
1409+
1410+ if (!IsEventAccepted (coll, true ))
1411+ continue ;
1412+
1413+ // Do analysis with collision-grouped V0s, retain full collision information
1414+ const uint64_t collIdx = coll.globalIndex ();
1415+ if (collIdx < CollIDBuffer)
1416+ LOGF (fatal, " Collision table unsorted! Previous index: %i, current index: %i" , CollIDBuffer, collIdx);
1417+
1418+ CollIDBuffer = collIdx;
1419+
1420+ auto V0s = fullV0s.sliceBy (perCollisionSTDSorted, collIdx);
1421+ float centrality = doPPAnalysis ? coll.centFT0M () : coll.centFT0C ();
1422+
1423+ // _______________________________________________
1424+ // Retrieving IR info
1425+ float interactionRate = -1 ;
1426+ if (fGetIR ) {
1427+ interactionRate = rateFetcher.fetch (ccdb.service , coll.timestamp (), coll.runNumber (), irSource, fIRCrashOnNull ) * 1 .e -3 ;
1428+ if (interactionRate < 0 )
1429+ histos.get <TH1>(HIST (" GeneralQA/hRunNumberNegativeIR" ))->Fill (Form (" %d" , coll.runNumber ()), 1 );
1430+
1431+ histos.fill (HIST (" GeneralQA/hInteractionRate" ), interactionRate);
1432+ histos.fill (HIST (" GeneralQA/hCentralityVsInteractionRate" ), centrality, interactionRate);
1433+ }
1434+
1435+ std::vector<int > bestGammasArray;
1436+ std::vector<int > bestLambdasArray;
1437+
1438+ // _______________________________________________
1439+ // V0s loop
1440+ for (auto & v0 : V0s) {
1441+ if (processPhotonCandidate (v0, coll)) // selecting photons
1442+ bestGammasArray.push_back (v0.globalIndex ()); // Save indices of best gamma candidates
1443+
1444+ if (processLambdaCandidate (v0, coll)) // selecting lambdas
1445+ bestLambdasArray.push_back (v0.globalIndex ()); // Save indices of best lambda candidates
1446+ }
1447+
1448+ // _______________________________________________
1449+ // Pi0 optional loop
1450+ if (doPi0QA) {
1451+ for (size_t i = 0 ; i < bestGammasArray.size (); ++i) {
1452+ auto gamma1 = fullV0s.rawIteratorAt (bestGammasArray[i]);
1453+ for (size_t j = i + 1 ; j < bestGammasArray.size (); ++j) {
1454+ auto gamma2 = fullV0s.rawIteratorAt (bestGammasArray[j]);
1455+ runPi0QA (gamma1, gamma2, coll);
1456+ }
1457+ }
13861458 }
1459+
1460+ // _______________________________________________
1461+ // Sigma0 nested loop
1462+ for (size_t i = 0 ; i < bestGammasArray.size (); ++i) {
1463+ auto gamma = fullV0s.rawIteratorAt (bestGammasArray[i]);
1464+
1465+ for (size_t j = 0 ; j < bestLambdasArray.size (); ++j) {
1466+ auto lambda = fullV0s.rawIteratorAt (bestLambdasArray[j]);
1467+
1468+ // Building sigma0 candidate
1469+ if (!buildSigma0 (lambda, gamma, coll))
1470+ continue ;
1471+
1472+ // Filling tables with accepted candidates
1473+ fillTables (lambda, gamma, coll);
1474+
1475+ nSigmaCandidates++;
1476+ if (nSigmaCandidates % 10000 == 0 )
1477+ LOG (info) << " Sigma0 Candidates built: " << nSigmaCandidates;
1478+ }
1479+ }
13871480 }
1481+ auto end = std::chrono::high_resolution_clock::now ();
1482+ std::chrono::duration<double > elapsed = end - start;
1483+
1484+ if (fverbose) LOGF (info, " [Process function call, PreSliceSorted] N. Collisions: %i, N. V0s: %i, Processing time (s): %lf" , collisions.size (), fullV0s.size (), elapsed.count ());
1485+ }
1486+
1487+ void processRealDataIterator (soa::Join<aod::StraCollisions, aod::StraCents, aod::StraEvSels, aod::StraStamps>::iterator const & coll, V0StandardDerivedDatas const & fullV0s, dauTracks const &)
1488+ {
1489+ auto start = std::chrono::high_resolution_clock::now ();
1490+
1491+ if (!IsEventAccepted (coll, true ))
1492+ continue ;
1493+
1494+ float centrality = doPPAnalysis ? coll.centFT0M () : coll.centFT0C ();
1495+
1496+ // _______________________________________________
1497+ // Retrieving IR info
1498+ float interactionRate = -1 ;
1499+ if (fGetIR ) {
1500+ interactionRate = rateFetcher.fetch (ccdb.service , coll.timestamp (), coll.runNumber (), irSource, fIRCrashOnNull ) * 1 .e -3 ;
1501+ if (interactionRate < 0 )
1502+ histos.get <TH1>(HIST (" GeneralQA/hRunNumberNegativeIR" ))->Fill (Form (" %d" , coll.runNumber ()), 1 );
1503+
1504+ histos.fill (HIST (" GeneralQA/hInteractionRate" ), interactionRate);
1505+ histos.fill (HIST (" GeneralQA/hCentralityVsInteractionRate" ), centrality, interactionRate);
1506+ }
1507+
1508+ std::vector<int > bestGammasArray;
1509+ std::vector<int > bestLambdasArray;
1510+
1511+ // _______________________________________________
1512+ // V0s loop
1513+ for (auto & v0 : fullV0s) {
1514+ if (processPhotonCandidate (v0, coll)) // selecting photons
1515+ bestGammasArray.push_back (v0.globalIndex () - fullV0s.offset ()); // Save indices of best gamma candidates
1516+
1517+ if (processLambdaCandidate (v0, coll)) // selecting lambdas
1518+ bestLambdasArray.push_back (v0.globalIndex () - fullV0s.offset ()); // Save indices of best lambda candidates
1519+ }
1520+
1521+ // _______________________________________________
1522+ // Pi0 optional loop
1523+ if (doPi0QA) {
1524+ for (size_t i = 0 ; i < bestGammasArray.size (); ++i) {
1525+ auto gamma1 = fullV0s.rawIteratorAt (bestGammasArray[i]);
1526+ for (size_t j = i + 1 ; j < bestGammasArray.size (); ++j) {
1527+ auto gamma2 = fullV0s.rawIteratorAt (bestGammasArray[j]);
1528+ runPi0QA (gamma1, gamma2, coll);
1529+ }
1530+ }
1531+ }
1532+
1533+ // _______________________________________________
1534+ // Sigma0 nested loop
1535+ for (size_t i = 0 ; i < bestGammasArray.size (); ++i) {
1536+ auto gamma = fullV0s.rawIteratorAt (bestGammasArray[i]);
1537+
1538+ for (size_t j = 0 ; j < bestLambdasArray.size (); ++j) {
1539+ auto lambda = fullV0s.rawIteratorAt (bestLambdasArray[j]);
1540+
1541+ // Building sigma0 candidate
1542+ if (!buildSigma0 (lambda, gamma, coll))
1543+ continue ;
1544+
1545+ // Filling tables with accepted candidates
1546+ fillTables (lambda, gamma, coll);
1547+
1548+ nSigmaCandidates++;
1549+ if (nSigmaCandidates % 10000 == 0 )
1550+ LOG (info) << " Sigma0 Candidates built: " << nSigmaCandidates;
1551+ }
1552+ }
1553+
1554+ auto end = std::chrono::high_resolution_clock::now ();
1555+ std::chrono::duration<double > elapsed = end - start;
1556+
1557+ if (fverbose) LOGF (info, " [Process function call, Iterator] N. V0s per collision: %i, Processing time (s): %lf" , fullV0s.size (), elapsed.count ());
1558+ }
1559+
1560+ void processRealDataDavid (soa::Join<aod::StraCollisions, aod::StraCents, aod::StraEvSels, aod::StraStamps> const & collisions, V0StandardDerivedDatas const & fullV0s, dauTracks const &)
1561+ {
1562+ auto start = std::chrono::high_resolution_clock::now ();
1563+
1564+ // brute force grouped index construction
1565+ std::vector<std::vector<int >> v0grouped (collisions.size ());
1566+
1567+ for (const auto & v0 : fullV0s) {
1568+ v0grouped[v0.straCollisionId ()].push_back (v0.globalIndex ());
1569+ }
1570+
1571+ for (const auto & coll : collisions) {
1572+
1573+ if (!IsEventAccepted (coll, true ))
1574+ continue ;
1575+
1576+ float centrality = doPPAnalysis ? coll.centFT0M () : coll.centFT0C ();
1577+
1578+ std::vector<int > bestGammasArray;
1579+ std::vector<int > bestLambdasArray;
1580+
1581+ // _______________________________________________
1582+ // V0s loop
1583+ for (size_t i; i < v0grouped[coll.globalIndex ()].size (); i++) {
1584+ auto v0 = fullV0s.rawIteratorAt (v0grouped[coll.globalIndex ()][i]);
1585+ if (processPhotonCandidate (v0, coll)) // selecting photons
1586+ bestGammasArray.push_back (v0.globalIndex ()); // Save indices of best gamma candidates
1587+
1588+ if (processLambdaCandidate (v0, coll)) // selecting lambdas
1589+ bestLambdasArray.push_back (v0.globalIndex ()); // Save indices of best lambda candidates
1590+ }
1591+
1592+ // _______________________________________________
1593+ // Pi0 optional loop
1594+ if (doPi0QA) {
1595+ for (size_t i = 0 ; i < bestGammasArray.size (); ++i) {
1596+ auto gamma1 = fullV0s.rawIteratorAt (bestGammasArray[i]);
1597+ for (size_t j = i + 1 ; j < bestGammasArray.size (); ++j) {
1598+ auto gamma2 = fullV0s.rawIteratorAt (bestGammasArray[j]);
1599+ runPi0QA (gamma1, gamma2, coll);
1600+ }
1601+ }
1602+ }
1603+
1604+ // _______________________________________________
1605+ // Sigma0 nested loop
1606+ for (size_t i = 0 ; i < bestGammasArray.size (); ++i) {
1607+ auto gamma = fullV0s.rawIteratorAt (bestGammasArray[i]);
1608+
1609+ for (size_t j = 0 ; j < bestLambdasArray.size (); ++j) {
1610+ auto lambda = fullV0s.rawIteratorAt (bestLambdasArray[j]);
1611+
1612+ // Building sigma0 candidate
1613+ if (!buildSigma0 (lambda, gamma, coll))
1614+ continue ;
1615+
1616+ // Filling tables with accepted candidates
1617+ fillTables (lambda, gamma, coll);
1618+
1619+ nSigmaCandidates++;
1620+ if (nSigmaCandidates % 10000 == 0 )
1621+ LOG (info) << " Sigma0 Candidates built: " << nSigmaCandidates;
1622+ }
1623+ }
1624+ }
1625+ auto end = std::chrono::high_resolution_clock::now ();
1626+ std::chrono::duration<double > elapsed = end - start;
1627+
1628+ if (fverbose) LOGF (info, " [David's process function call] N. Collisions: %i, N. V0s: %i, Processing time (s): %lf" , collisions.size (), fullV0s.size (), elapsed.count ());
1629+
13881630 }
13891631
13901632 // Simulated processing in Run 3 (subscribes to MC information too)
@@ -1395,6 +1637,9 @@ struct sigma0builder {
13951637
13961638 PROCESS_SWITCH (sigma0builder, processMonteCarlo, " process as if MC data" , false );
13971639 PROCESS_SWITCH (sigma0builder, processRealData, " process as if real data" , true );
1640+ PROCESS_SWITCH (sigma0builder, processRealDataSorted, " process as if real data. QA only." , true );
1641+ PROCESS_SWITCH (sigma0builder, processRealDataIterator, " process as if real data. QA only." , true );
1642+ PROCESS_SWITCH (sigma0builder, processRealDataDavid, " process as if real data. QA only." , true );
13981643 PROCESS_SWITCH (sigma0builder, processGeneratedRun3, " process generated MC collisions" , false );
13991644};
14001645
0 commit comments