diff --git a/RecoTracker/LSTCore/interface/Constants.h b/RecoTracker/LSTCore/interface/Constants.h index 5934a1af969ee..fe78da2cca67e 100644 --- a/RecoTracker/LSTCore/interface/Constants.h +++ b/RecoTracker/LSTCore/interface/Constants.h @@ -41,6 +41,7 @@ namespace lst { constexpr unsigned int n_max_pixel_md_per_modules = 2 * n_max_pixel_segments_per_module; + constexpr unsigned int n_max_pixel_segments = 15000; constexpr unsigned int n_max_pixel_triplets = 5000; constexpr unsigned int n_max_pixel_quintuplets = 15000; @@ -60,6 +61,9 @@ namespace lst { struct Params_T3 { static constexpr int kLayers = 3, kHits = 6; }; + struct Params_pT2 { + static constexpr int kLayers = 4, kHits = 8; + }; struct Params_pT3 { static constexpr int kLayers = 5, kHits = 10; }; diff --git a/RecoTracker/LSTCore/src/alpaka/Event.dev.cc b/RecoTracker/LSTCore/src/alpaka/Event.dev.cc index 2d15457bc3c38..aaffda08a31ac 100644 --- a/RecoTracker/LSTCore/src/alpaka/Event.dev.cc +++ b/RecoTracker/LSTCore/src/alpaka/Event.dev.cc @@ -12,6 +12,7 @@ void lst::Event::init(bool verbose) { trackCandidatesInGPU = nullptr; pixelTripletsInGPU = nullptr; pixelQuintupletsInGPU = nullptr; + PT2sInGPU = nullptr; rangesInGPU = nullptr; hitsInCPU = nullptr; @@ -24,6 +25,7 @@ void lst::Event::init(bool verbose) { quintupletsInCPU = nullptr; pixelTripletsInCPU = nullptr; pixelQuintupletsInCPU = nullptr; + PT2sInCPU = nullptr; //reset the arrays for (int i = 0; i < 6; i++) { @@ -107,6 +109,11 @@ void lst::Event::resetEvent() { delete pixelQuintupletsBuffers; pixelQuintupletsInGPU = nullptr; } + if (PT2sInGPU) { + delete PT2sInGPU; + delete PT2sBuffers; + PT2sInGPU = nullptr; + } if (hitsInCPU != nullptr) { delete hitsInCPU; @@ -140,6 +147,10 @@ void lst::Event::resetEvent() { delete pixelQuintupletsInCPU; pixelQuintupletsInCPU = nullptr; } + if (PT2sInCPU != nullptr) { + delete PT2sInCPU; + PT2sInCPU = nullptr; + } if (trackCandidatesInCPU != nullptr) { delete trackCandidatesInCPU; trackCandidatesInCPU = nullptr; @@ -730,6 +741,39 @@ void lst::Event::createTrackCandidates(bool no_pls_dupclean, bool tc_pls_ alpaka::enqueue(queue, checkHitspLSTask); } +/* + Vec3D const threadsPerBlock_crossCleanPT2{1, 16, 64}; + Vec3D const blocksPerGrid_crossCleanPT2{1, 4, 20}; + WorkDiv3D const crossClean2_workDiv = + createWorkDiv(blocksPerGrid_crossCleanPT2, threadsPerBlock_crossCleanPT2, elementsPerThread); + + lst::crossCleanPT2 crossCleanPT2_kernel; + auto const crossCleanPT2Task(alpaka::createTaskKernel(crossCleanPT2_workDiv, + crossCleanPT2_kernel, + *modulesBuffers_->data(), + *rangesInGPU, + *PT2sInGPU, + *segmentsInGPU, + *pixelTripletsInGPU, + *pixelQuintupletsInGPU)); + + alpaka::enqueue(queue, crossCleanPT2Task); +*/ + Vec3D const threadsPerBlock_addPT2asTrackCandidatesInGPU{1, 1, 512}; + Vec3D const blocksPerGrid_addPT2asTrackCandidatesInGPU{1, 1, 1}; + WorkDiv3D const addPT2asTrackCandidatesInGPU_workDiv = createWorkDiv( + blocksPerGrid_addPT2asTrackCandidatesInGPU, threadsPerBlock_addPT2asTrackCandidatesInGPU, elementsPerThread); + + lst::addPT2asTrackCandidatesInGPU addPT2asTrackCandidatesInGPU_kernel; + auto const addPT2asTrackCandidatesInGPUTask(alpaka::createTaskKernel(addPT2asTrackCandidatesInGPU_workDiv, + addPT2asTrackCandidatesInGPU_kernel, + nLowerModules_, + *PT2sInGPU, + *trackCandidatesInGPU, + *segmentsInGPU, + *rangesInGPU)); + + alpaka::enqueue(queue, addPT2asTrackCandidatesInGPUTask); Vec3D const threadsPerBlock_crossCleanpLS{1, 16, 32}; Vec3D const blocksPerGrid_crossCleanpLS{1, 4, 20}; @@ -768,19 +812,22 @@ void lst::Event::createTrackCandidates(bool no_pls_dupclean, bool tc_pls_ // Check if either n_max_pixel_track_candidates or n_max_nonpixel_track_candidates was reached auto nTrackCanpT5Host_buf = allocBufWrapper(devHost, 1, queue); auto nTrackCanpT3Host_buf = allocBufWrapper(devHost, 1, queue); + auto nTrackCanPT2Host_buf = allocBufWrapper(devHost, 1, queue); auto nTrackCanpLSHost_buf = allocBufWrapper(devHost, 1, queue); auto nTrackCanT5Host_buf = allocBufWrapper(devHost, 1, queue); alpaka::memcpy(queue, nTrackCanpT5Host_buf, trackCandidatesBuffers->nTrackCandidatespT5_buf); alpaka::memcpy(queue, nTrackCanpT3Host_buf, trackCandidatesBuffers->nTrackCandidatespT3_buf); + alpaka::memcpy(queue, nTrackCanPT2Host_buf, trackCandidatesBuffers->nTrackCandidatesPT2_buf); alpaka::memcpy(queue, nTrackCanpLSHost_buf, trackCandidatesBuffers->nTrackCandidatespLS_buf); alpaka::memcpy(queue, nTrackCanT5Host_buf, trackCandidatesBuffers->nTrackCandidatesT5_buf); alpaka::wait(queue); int nTrackCandidatespT5 = *alpaka::getPtrNative(nTrackCanpT5Host_buf); int nTrackCandidatespT3 = *alpaka::getPtrNative(nTrackCanpT3Host_buf); + int nTrackCandidatesPT2 = *alpaka::getPtrNative(nTrackCanPT2Host_buf); int nTrackCandidatespLS = *alpaka::getPtrNative(nTrackCanpLSHost_buf); int nTrackCandidatesT5 = *alpaka::getPtrNative(nTrackCanT5Host_buf); - if ((nTrackCandidatespT5 + nTrackCandidatespT3 + nTrackCandidatespLS == n_max_pixel_track_candidates) || + if ((nTrackCandidatespT5 + nTrackCandidatespT3 + nTrackCandidatespLS + nTrackCandidatesPT2 == n_max_pixel_track_candidates) || (nTrackCandidatesT5 == n_max_nonpixel_track_candidates)) { printf( "****************************************************************************************************\n" @@ -911,6 +958,125 @@ void lst::Event::createPixelTriplets() { alpaka::wait(queue); } +void lst::Event::createPT2s() { + if (PT2sInGPU == nullptr) { + PT2sInGPU = new lst::PT2s(); + PT2sBuffers = new lst::PT2sBuffer(n_max_pixel_segments, devAcc, queue); + PT2sInGPU->setData(*PT2sBuffers); + } + + unsigned int nInnerSegments; + auto nInnerSegments_src_view = alpaka::createView(devHost, &nInnerSegments, (size_t)1u); + + auto dev_view_nSegments = alpaka::createSubView(segmentsBuffers->nSegments_buf, (Idx)1u, (Idx)nLowerModules_); + + alpaka::memcpy(queue, nInnerSegments_src_view, dev_view_nSegments); + alpaka::wait(queue); + + auto superbins_buf = allocBufWrapper(devHost, n_max_pixel_segments_per_module, queue); + auto pixelTypes_buf = allocBufWrapper(devHost, n_max_pixel_segments_per_module, queue); + + alpaka::memcpy(queue, superbins_buf, segmentsBuffers->superbin_buf); + alpaka::memcpy(queue, pixelTypes_buf, segmentsBuffers->pixelType_buf); + alpaka::wait(queue); + + auto connectedPixelSize_host_buf = allocBufWrapper(devHost, nInnerSegments, queue); + auto connectedPixelIndex_host_buf = allocBufWrapper(devHost, nInnerSegments, queue); + auto connectedPixelSize_dev_buf = allocBufWrapper(devAcc, nInnerSegments, queue); + auto connectedPixelIndex_dev_buf = allocBufWrapper(devAcc, nInnerSegments, queue); + + int* superbins = alpaka::getPtrNative(superbins_buf); + int8_t* pixelTypes = alpaka::getPtrNative(pixelTypes_buf); + unsigned int* connectedPixelSize_host = alpaka::getPtrNative(connectedPixelSize_host_buf); + unsigned int* connectedPixelIndex_host = alpaka::getPtrNative(connectedPixelIndex_host_buf); + alpaka::wait(queue); + + int pixelIndexOffsetPos = + pixelMapping_->connectedPixelsIndex[size_superbins - 1] + pixelMapping_->connectedPixelsSizes[size_superbins - 1]; + int pixelIndexOffsetNeg = pixelMapping_->connectedPixelsIndexPos[size_superbins - 1] + + pixelMapping_->connectedPixelsSizesPos[size_superbins - 1] + pixelIndexOffsetPos; + + // TODO: check if a map/reduction to just eligible pLSs would speed up the kernel + // the current selection still leaves a significant fraction of unmatchable pLSs + for (unsigned int i = 0; i < nInnerSegments; i++) { // loop over # pLS + int8_t pixelType = pixelTypes[i]; // Get pixel type for this pLS + int superbin = superbins[i]; // Get superbin for this pixel + if ((superbin < 0) or (superbin >= (int)size_superbins) or (pixelType > 2) or (pixelType < 0)) { + connectedPixelSize_host[i] = 0; + connectedPixelIndex_host[i] = 0; + continue; + } + + // Used pixel type to select correct size-index arrays + if (pixelType == 0) { + connectedPixelSize_host[i] = + pixelMapping_->connectedPixelsSizes[superbin]; // number of connected modules to this pixel + auto connectedIdxBase = pixelMapping_->connectedPixelsIndex[superbin]; + connectedPixelIndex_host[i] = + connectedIdxBase; // index to get start of connected modules for this superbin in map + } else if (pixelType == 1) { + connectedPixelSize_host[i] = + pixelMapping_->connectedPixelsSizesPos[superbin]; // number of pixel connected modules + auto connectedIdxBase = pixelMapping_->connectedPixelsIndexPos[superbin] + pixelIndexOffsetPos; + connectedPixelIndex_host[i] = connectedIdxBase; // index to get start of connected pixel modules + } else if (pixelType == 2) { + connectedPixelSize_host[i] = + pixelMapping_->connectedPixelsSizesNeg[superbin]; // number of pixel connected modules + auto connectedIdxBase = pixelMapping_->connectedPixelsIndexNeg[superbin] + pixelIndexOffsetNeg; + connectedPixelIndex_host[i] = connectedIdxBase; // index to get start of connected pixel modules + } + } + + alpaka::memcpy(queue, connectedPixelSize_dev_buf, connectedPixelSize_host_buf, nInnerSegments); + alpaka::memcpy(queue, connectedPixelIndex_dev_buf, connectedPixelIndex_host_buf, nInnerSegments); + alpaka::wait(queue); + + Vec3D const threadsPerBlock{1, 4, 32}; + Vec3D const blocksPerGrid{16 /* above median of connected modules*/, 4096, 1}; + WorkDiv3D const createPT2sInGPUFromMapv2_workDiv = + createWorkDiv(blocksPerGrid, threadsPerBlock, elementsPerThread); + + lst::createPT2sInGPUFromMapv2 createPT2sInGPUFromMapv2_kernel; + auto const createPT2sInGPUFromMapv2Task( + alpaka::createTaskKernel(createPT2sInGPUFromMapv2_workDiv, + createPT2sInGPUFromMapv2_kernel, + *modulesBuffers_->data(), + *rangesInGPU, + *mdsInGPU, + *segmentsInGPU, + *PT2sInGPU, + alpaka::getPtrNative(connectedPixelSize_dev_buf), + alpaka::getPtrNative(connectedPixelIndex_dev_buf), + nInnerSegments, + ptCut)); + + alpaka::enqueue(queue, createPT2sInGPUFromMapv2Task); + alpaka::wait(queue); + +#ifdef WARNINGS + auto nPT2s_buf = allocBufWrapper(devHost, 1, queue); + + alpaka::memcpy(queue, nPT2s_buf, pT2sBuffers->nPT2s_buf); + alpaka::wait(queue); + + std::cout << "number of pixel triplets = " << *alpaka::getPtrNative(nPT2s_buf) << std::endl; +#endif + + //pT3s can be cleaned here because they're not used in making pT5s! + Vec3D const threadsPerBlockDupPixTrip{1, 16, 16}; + //seems like more blocks lead to conflicting writes + Vec3D const blocksPerGridDupPixTrip{1, 40, 1}; + WorkDiv3D const removeDupPT2sInGPUFromMap_workDiv = + createWorkDiv(blocksPerGridDupPixTrip, threadsPerBlockDupPixTrip, elementsPerThread); + + lst::removeDupPT2sInGPUFromMap removeDupPT2sInGPUFromMap_kernel; + auto const removeDupPT2sInGPUFromMapTask(alpaka::createTaskKernel( + removeDupPT2sInGPUFromMap_workDiv, removeDupPT2sInGPUFromMap_kernel, *PT2sInGPU)); + + alpaka::enqueue(queue, removeDupPT2sInGPUFromMapTask); + alpaka::wait(queue); +} + void lst::Event::createQuintuplets() { Vec3D const threadsPerBlockCreateQuints{1, 1, 1024}; Vec3D const blocksPerGridCreateQuints{1, 1, 1}; @@ -1382,6 +1548,17 @@ unsigned int lst::Event::getNumberOfTripletsByLayerEndcap(unsigned int la return n_triplets_by_layer_endcap_[layer]; } +int lst::Event::getNumberOfPT2s() { + auto nPT2s_buf = allocBufWrapper(devHost, 1, queue); + + alpaka::memcpy(queue, nPT2s_buf, PT2sBuffers->nPT2s_buf); + alpaka::wait(queue); + + int nPT2s = *alpaka::getPtrNative(nPT2s_buf); + + return nPT2s; +} + int lst::Event::getNumberOfPixelTriplets() { auto nPixelTriplets_buf = allocBufWrapper(devHost, 1, queue); @@ -1464,6 +1641,17 @@ int lst::Event::getNumberOfPT3TrackCandidates() { return nTrackCandidatesPT3; } +int lst::Event::getNumberOfPT2TrackCandidates() { + auto nTrackCandidatesPT2_buf = allocBufWrapper(devHost, 1, queue); + + alpaka::memcpy(queue, nTrackCandidatesPT2_buf, trackCandidatesBuffers->nTrackCandidatesPT2_buf); + alpaka::wait(queue); + + int nTrackCandidatesPT2 = *alpaka::getPtrNative(nTrackCandidatesPT2_buf); + + return nTrackCandidatesPT2; +} + int lst::Event::getNumberOfPLSTrackCandidates() { auto nTrackCandidatesPLS_buf = allocBufWrapper(devHost, 1, queue); @@ -1572,6 +1760,10 @@ lst::MiniDoubletsBuffer* lst::Event::getMiniDoublets() { alpaka::memcpy(queue, mdsInCPU->dphichanges_buf, miniDoubletsBuffers->dphichanges_buf, nMemHost); alpaka::memcpy(queue, mdsInCPU->nMDs_buf, miniDoubletsBuffers->nMDs_buf); alpaka::memcpy(queue, mdsInCPU->totOccupancyMDs_buf, miniDoubletsBuffers->totOccupancyMDs_buf); + alpaka::memcpy(queue, mdsInCPU->anchorX_buf, miniDoubletsBuffers->anchorX_buf); + alpaka::memcpy(queue, mdsInCPU->anchorY_buf, miniDoubletsBuffers->anchorY_buf); + alpaka::memcpy(queue, mdsInCPU->anchorZ_buf, miniDoubletsBuffers->anchorZ_buf); + alpaka::memcpy(queue, mdsInCPU->anchorRt_buf, miniDoubletsBuffers->anchorRt_buf); alpaka::wait(queue); } return mdsInCPU; @@ -1687,6 +1879,36 @@ lst::QuintupletsBuffer* lst::Event::getQuintuplets() { return quintupletsInCPU; } +lst::PT2sBuffer* lst::Event::getPT2s() { + if (PT2sInCPU == nullptr) { + // Get nPT2s parameter to initialize host based PT2sInCPU + auto nPT2s_buf = allocBufWrapper(devHost, 1, queue); + alpaka::memcpy(queue, nPT2s_buf, PT2sBuffers->nPT2s_buf); + alpaka::wait(queue); + + unsigned int nPT2s = *alpaka::getPtrNative(nPT2s_buf); + PT2sInCPU = new lst::PT2sBuffer(nPT2s, devHost, queue); + PT2sInCPU->setData(*PT2sInCPU); + + *alpaka::getPtrNative(PT2sInCPU->nPT2s_buf) = nPT2s; + alpaka::memcpy( + queue, PT2sInCPU->totOccupancyPT2s_buf, PT2sBuffers->totOccupancyPT2s_buf); + alpaka::memcpy(queue, PT2sInCPU->rzChiSquared_buf, PT2sBuffers->rzChiSquared_buf, nPT2s); + alpaka::memcpy( + queue, PT2sInCPU->segmentIndices_buf, PT2sBuffers->segmentIndices_buf, nPT2s); + alpaka::memcpy(queue, + PT2sInCPU->pixelSegmentIndices_buf, + PT2sBuffers->pixelSegmentIndices_buf, + nPT2s); + alpaka::memcpy(queue, PT2sInCPU->pixelRadius_buf, PT2sBuffers->pixelRadius_buf, nPT2s); + alpaka::memcpy(queue, PT2sInCPU->isDup_buf, PT2sBuffers->isDup_buf, nPT2s); + alpaka::memcpy(queue, PT2sInCPU->eta_buf, PT2sBuffers->eta_buf, nPT2s); + alpaka::memcpy(queue, PT2sInCPU->phi_buf, PT2sBuffers->phi_buf, nPT2s); + alpaka::wait(queue); + } + return PT2sInCPU; +} + lst::PixelTripletsBuffer* lst::Event::getPixelTriplets() { if (pixelTripletsInCPU == nullptr) { // Get nPixelTriplets parameter to initialize host based quintupletsInCPU diff --git a/RecoTracker/LSTCore/src/alpaka/Event.h b/RecoTracker/LSTCore/src/alpaka/Event.h index 3df3728cc8854..1e88aac42c89c 100644 --- a/RecoTracker/LSTCore/src/alpaka/Event.h +++ b/RecoTracker/LSTCore/src/alpaka/Event.h @@ -13,6 +13,7 @@ #include "MiniDoublet.h" #include "PixelQuintuplet.h" #include "PixelTriplet.h" +#include "PT2.h" #include "TrackCandidate.h" #include "HeterogeneousCore/AlpakaInterface/interface/host.h" @@ -62,6 +63,8 @@ namespace lst { QuintupletsBuffer* quintupletsBuffers; TrackCandidates* trackCandidatesInGPU; TrackCandidatesBuffer* trackCandidatesBuffers; + PT2s* PT2sInGPU; + PT2sBuffer* PT2sBuffers; PixelTriplets* pixelTripletsInGPU; PixelTripletsBuffer* pixelTripletsBuffers; PixelQuintuplets* pixelQuintupletsInGPU; @@ -76,6 +79,7 @@ namespace lst { TrackCandidatesBuffer* trackCandidatesInCPU; ModulesBuffer* modulesInCPU; QuintupletsBuffer* quintupletsInCPU; + PT2sBuffer* PT2sInCPU; PixelTripletsBuffer* pixelTripletsInCPU; PixelQuintupletsBuffer* pixelQuintupletsInCPU; @@ -155,6 +159,7 @@ namespace lst { void createTrackCandidates(bool no_pls_dupclean, bool tc_pls_triplets); void createExtendedTracks(); void createQuintuplets(); + void createPT2s(); void createPixelTriplets(); void createPixelQuintuplets(); void pixelLineSegmentCleaning(bool no_pls_dupclean); @@ -183,6 +188,7 @@ namespace lst { int getNumberOfPixelTrackCandidates(); int getNumberOfPT5TrackCandidates(); int getNumberOfPT3TrackCandidates(); + int getNumberOfPT2TrackCandidates(); int getNumberOfT5TrackCandidates(); int getNumberOfPLSTrackCandidates(); @@ -191,6 +197,7 @@ namespace lst { unsigned int getNumberOfQuintupletsByLayerBarrel(unsigned int layer); unsigned int getNumberOfQuintupletsByLayerEndcap(unsigned int layer); + int getNumberOfPT2s(); int getNumberOfPixelTriplets(); int getNumberOfPixelQuintuplets(); @@ -203,6 +210,7 @@ namespace lst { QuintupletsBuffer* getQuintuplets(); TrackCandidatesBuffer* getTrackCandidates(); TrackCandidatesBuffer* getTrackCandidatesInCMSSW(); + PT2sBuffer* getPT2s(); PixelTripletsBuffer* getPixelTriplets(); PixelQuintupletsBuffer* getPixelQuintuplets(); ModulesBuffer* getModules(bool isFull = false); diff --git a/RecoTracker/LSTCore/src/alpaka/Kernels.h b/RecoTracker/LSTCore/src/alpaka/Kernels.h index 77b1033349845..e9e36d0048bf1 100644 --- a/RecoTracker/LSTCore/src/alpaka/Kernels.h +++ b/RecoTracker/LSTCore/src/alpaka/Kernels.h @@ -12,6 +12,7 @@ #include "Quintuplet.h" #include "PixelQuintuplet.h" #include "PixelTriplet.h" +#include "PT2.h" namespace lst { ALPAKA_FN_ACC ALPAKA_FN_INLINE void rmQuintupletFromMemory(lst::Quintuplets& quintupletsInGPU, @@ -25,6 +26,11 @@ namespace lst { pixelTripletsInGPU.isDup[pixelTripletIndex] = true; }; + ALPAKA_FN_ACC ALPAKA_FN_INLINE void rmPT2FromMemory(lst::PT2s& PT2sInGPU, + unsigned int PT2Index) { + PT2sInGPU.isDup[PT2Index] = true; + }; + ALPAKA_FN_ACC ALPAKA_FN_INLINE void rmPixelQuintupletFromMemory(lst::PixelQuintuplets& pixelQuintupletsInGPU, unsigned int pixelQuintupletIndex) { pixelQuintupletsInGPU.isDup[pixelQuintupletIndex] = true; @@ -142,6 +148,58 @@ namespace lst { matched[1] = nMatched; }; + ALPAKA_FN_ACC ALPAKA_FN_INLINE void checkHitspT2(unsigned int ix, + unsigned int jx, + lst::PT2s const& PT2sInGPU, + int* matched) { + int phits1[Params_pLS::kHits]; + int phits2[Params_pLS::kHits]; + + for (int i = 0; i < Params_pLS::kHits; i++) { + phits1[i] = PT2sInGPU.hitIndices[Params_pT2::kHits * ix + i]; + phits2[i] = PT2sInGPU.hitIndices[Params_pT2::kHits * jx + i]; + } + + int npMatched = 0; + for (int i = 0; i < Params_pLS::kHits; i++) { + bool pmatched = false; + for (int j = 0; j < Params_pLS::kHits; j++) { + if (phits1[i] == phits2[j]) { + pmatched = true; + break; + } + } + if (pmatched) { + npMatched++; + } + } + + int hits1[Params_LS::kHits]; + int hits2[Params_LS::kHits]; + + for (int i = 0; i < Params_LS::kHits; i++) { + hits1[i] = PT2sInGPU.hitIndices[Params_pT2::kHits * ix + i + 4]; // Omitting the pLS hits + hits2[i] = PT2sInGPU.hitIndices[Params_pT2::kHits * jx + i + 4]; // Omitting the pLS hits + } + + int nMatched = 0; + for (int i = 0; i < Params_LS::kHits; i++) { + bool tmatched = false; + for (int j = 0; j < Params_LS::kHits; j++) { + if (hits1[i] == hits2[j]) { + tmatched = true; + break; + } + } + if (tmatched) { + nMatched++; + } + } + + matched[0] = npMatched; + matched[1] = nMatched; + }; + struct removeDupQuintupletsInGPUAfterBuild { template ALPAKA_FN_ACC void operator()(TAcc const& acc, @@ -274,6 +332,41 @@ namespace lst { } }; + struct removeDupPT2sInGPUFromMap { + template + ALPAKA_FN_ACC void operator()(TAcc const& acc, lst::PT2s PT2sInGPU) const { + auto const globalThreadIdx = alpaka::getIdx(acc); + auto const gridThreadExtent = alpaka::getWorkDiv(acc); + + for (unsigned int ix = globalThreadIdx[1]; ix < *PT2sInGPU.nPT2s; ix += gridThreadExtent[1]) { + for (unsigned int jx = globalThreadIdx[2]; jx < *PT2sInGPU.nPT2s; jx += gridThreadExtent[2]) { + if (ix == jx) + continue; + + int nMatched[2]; + checkHitspT2(ix, jx, PT2sInGPU, nMatched); + const int minNHitsForDup_PT2 = 6; + if ((nMatched[0] + nMatched[1]) >= minNHitsForDup_PT2) { + // Check the layers + if (PT2sInGPU.logicalLayers[Params_pT2::kLayers * jx + 2] < + PT2sInGPU.logicalLayers[Params_pT2::kLayers * ix + 2]) { + rmPT2FromMemory(PT2sInGPU, ix); + break; + } else if (PT2sInGPU.logicalLayers[Params_pT2::kLayers * ix + 2] == + PT2sInGPU.logicalLayers[Params_pT2::kLayers * jx + 2]) { + rmPT2FromMemory(PT2sInGPU, ix); + break; + } else if (PT2sInGPU.logicalLayers[Params_pT2::kLayers * ix + 2] == + PT2sInGPU.logicalLayers[Params_pT2::kLayers * jx + 2] && (ix < jx)) { + rmPT2FromMemory(PT2sInGPU, ix); + break; + } + } + } + } + } + }; + struct removeDupPixelTripletsInGPUFromMap { template ALPAKA_FN_ACC void operator()(TAcc const& acc, lst::PixelTriplets pixelTripletsInGPU) const { diff --git a/RecoTracker/LSTCore/src/alpaka/LST.dev.cc b/RecoTracker/LSTCore/src/alpaka/LST.dev.cc index 063404b963f2a..4a9a091627930 100644 --- a/RecoTracker/LSTCore/src/alpaka/LST.dev.cc +++ b/RecoTracker/LSTCore/src/alpaka/LST.dev.cc @@ -419,6 +419,10 @@ void lst::LST::run(Queue& queue, if (verbose) printf("# of Pixel T3s produced: %d\n", event.getNumberOfPixelTriplets()); + event.createPT2s(); + if (verbose) + printf("# of PT2s produced: %d\n", event.getNumberOfPT2s()); + event.createTrackCandidates(no_pls_dupclean, tc_pls_triplets); if (verbose) { printf("# of TrackCandidates produced: %d\n", event.getNumberOfTrackCandidates()); @@ -427,6 +431,7 @@ void lst::LST::run(Queue& queue, printf(" # of pT3 TrackCandidates produced: %d\n", event.getNumberOfPT3TrackCandidates()); printf(" # of pLS TrackCandidates produced: %d\n", event.getNumberOfPLSTrackCandidates()); printf(" # of T5 TrackCandidates produced: %d\n", event.getNumberOfT5TrackCandidates()); + printf(" # of pT2 TrackCandidates produced: %d\n", event.getNumberOfPT2TrackCandidates()); } getOutput(event); diff --git a/RecoTracker/LSTCore/src/alpaka/PT2.h b/RecoTracker/LSTCore/src/alpaka/PT2.h new file mode 100644 index 0000000000000..803bf671bef98 --- /dev/null +++ b/RecoTracker/LSTCore/src/alpaka/PT2.h @@ -0,0 +1,343 @@ +#ifndef RecoTracker_LSTCore_src_alpaka_PT2_h +#define RecoTracker_LSTCore_src_alpaka_PT2_h + +#include "RecoTracker/LSTCore/interface/alpaka/Constants.h" +#include "RecoTracker/LSTCore/interface/Module.h" + +#include "Segment.h" +#include "MiniDoublet.h" +#include "Hit.h" +#include "ObjectRanges.h" + +namespace lst { + // One pixel segment, one outer tracker triplet! + struct PT2s { + unsigned int* pixelSegmentIndices; + unsigned int* segmentIndices; + unsigned int* nPT2s; + unsigned int* totOccupancyPT2s; + + float* rzChiSquared; + + FPX* pixelRadius; + FPX* pt; + FPX* eta; + FPX* phi; + FPX* eta_pix; + FPX* phi_pix; + bool* isDup; + bool* partOfPT3; + bool* partOfPT5; + + uint8_t* logicalLayers; + unsigned int* hitIndices; + uint16_t* lowerModuleIndices; + FPX* centerX; + FPX* centerY; + + template + void setData(TBuff& buf) { + pixelSegmentIndices = alpaka::getPtrNative(buf.pixelSegmentIndices_buf); + segmentIndices = alpaka::getPtrNative(buf.segmentIndices_buf); + nPT2s = alpaka::getPtrNative(buf.nPT2s_buf); + totOccupancyPT2s = alpaka::getPtrNative(buf.totOccupancyPT2s_buf); + pixelRadius = alpaka::getPtrNative(buf.pixelRadius_buf); + pt = alpaka::getPtrNative(buf.pt_buf); + eta = alpaka::getPtrNative(buf.eta_buf); + phi = alpaka::getPtrNative(buf.phi_buf); + eta_pix = alpaka::getPtrNative(buf.eta_pix_buf); + phi_pix = alpaka::getPtrNative(buf.phi_pix_buf); + isDup = alpaka::getPtrNative(buf.isDup_buf); + partOfPT3 = alpaka::getPtrNative(buf.partOfPT3_buf); + partOfPT5 = alpaka::getPtrNative(buf.partOfPT5_buf); + logicalLayers = alpaka::getPtrNative(buf.logicalLayers_buf); + hitIndices = alpaka::getPtrNative(buf.hitIndices_buf); + lowerModuleIndices = alpaka::getPtrNative(buf.lowerModuleIndices_buf); + centerX = alpaka::getPtrNative(buf.centerX_buf); + centerY = alpaka::getPtrNative(buf.centerY_buf); + rzChiSquared = alpaka::getPtrNative(buf.rzChiSquared_buf); + } + }; + + template + struct PT2sBuffer { + Buf pixelSegmentIndices_buf; + Buf segmentIndices_buf; + Buf nPT2s_buf; + Buf totOccupancyPT2s_buf; + Buf pixelRadius_buf; + Buf pt_buf; + Buf eta_buf; + Buf phi_buf; + Buf eta_pix_buf; + Buf phi_pix_buf; + Buf isDup_buf; + Buf partOfPT3_buf; + Buf partOfPT5_buf; + Buf logicalLayers_buf; + Buf hitIndices_buf; + Buf lowerModuleIndices_buf; + Buf centerX_buf; + Buf centerY_buf; + Buf pixelRadiusError_buf; + Buf rzChiSquared_buf; + + PT2s data_; + + template + PT2sBuffer(unsigned int maxPT2s, TDevAcc const& devAccIn, TQueue& queue) + : pixelSegmentIndices_buf(allocBufWrapper(devAccIn, maxPT2s, queue)), + segmentIndices_buf(allocBufWrapper(devAccIn, maxPT2s, queue)), + nPT2s_buf(allocBufWrapper(devAccIn, 1, queue)), + totOccupancyPT2s_buf(allocBufWrapper(devAccIn, 1, queue)), + pixelRadius_buf(allocBufWrapper(devAccIn, maxPT2s, queue)), + pt_buf(allocBufWrapper(devAccIn, maxPT2s, queue)), + eta_buf(allocBufWrapper(devAccIn, maxPT2s, queue)), + phi_buf(allocBufWrapper(devAccIn, maxPT2s, queue)), + eta_pix_buf(allocBufWrapper(devAccIn, maxPT2s, queue)), + phi_pix_buf(allocBufWrapper(devAccIn, maxPT2s, queue)), + isDup_buf(allocBufWrapper(devAccIn, maxPT2s, queue)), + partOfPT3_buf(allocBufWrapper(devAccIn, maxPT2s, queue)), + partOfPT5_buf(allocBufWrapper(devAccIn, maxPT2s, queue)), + logicalLayers_buf(allocBufWrapper(devAccIn, maxPT2s * Params_pT2::kLayers, queue)), + hitIndices_buf(allocBufWrapper(devAccIn, maxPT2s * Params_pT2::kHits, queue)), + lowerModuleIndices_buf(allocBufWrapper(devAccIn, maxPT2s * Params_pT2::kLayers, queue)), + centerX_buf(allocBufWrapper(devAccIn, maxPT2s, queue)), + centerY_buf(allocBufWrapper(devAccIn, maxPT2s, queue)), + pixelRadiusError_buf(allocBufWrapper(devAccIn, maxPT2s, queue)), + rzChiSquared_buf(allocBufWrapper(devAccIn, maxPT2s, queue)) { + alpaka::memset(queue, nPT2s_buf, 0u); + alpaka::memset(queue, totOccupancyPT2s_buf, 0u); + alpaka::memset(queue, partOfPT3_buf, false); + alpaka::memset(queue, partOfPT5_buf, false); + alpaka::wait(queue); + } + + inline PT2s const* data() const { return &data_; } + inline void setData(PT2sBuffer& buf) { data_.setData(buf); } + }; + + ALPAKA_FN_ACC ALPAKA_FN_INLINE void addPT2ToMemory(lst::MiniDoublets const& mdsInGPU, + lst::Segments const& segmentsInGPU, + lst::PT2s& PT2sInGPU, + unsigned int pixelSegmentIndex, + unsigned int segmentIndex, + float pixelRadius, + float centerX, + float centerY, + float rzChiSquared, + unsigned int PT2Index, + float pt, + float eta, + float phi, + float eta_pix, + float phi_pix) { + PT2sInGPU.pixelSegmentIndices[PT2Index] = pixelSegmentIndex; + PT2sInGPU.segmentIndices[PT2Index] = segmentIndex; + PT2sInGPU.pixelRadius[PT2Index] = __F2H(pixelRadius); + PT2sInGPU.pt[PT2Index] = __F2H(pt); + PT2sInGPU.eta[PT2Index] = __F2H(eta); + PT2sInGPU.phi[PT2Index] = __F2H(phi); + PT2sInGPU.eta_pix[PT2Index] = __F2H(eta_pix); + PT2sInGPU.phi_pix[PT2Index] = __F2H(phi_pix); + PT2sInGPU.isDup[PT2Index] = false; + + PT2sInGPU.centerX[PT2Index] = __F2H(centerX); + PT2sInGPU.centerY[PT2Index] = __F2H(centerY); + PT2sInGPU.logicalLayers[Params_pT2::kLayers * PT2Index] = 0; + PT2sInGPU.logicalLayers[Params_pT2::kLayers * PT2Index + 1] = 0; + PT2sInGPU.logicalLayers[Params_pT2::kLayers * PT2Index + 2] = + segmentsInGPU.logicalLayers[segmentIndex * Params_LS::kLayers]; + PT2sInGPU.logicalLayers[Params_pT2::kLayers * PT2Index + 3] = + segmentsInGPU.logicalLayers[segmentIndex * Params_LS::kLayers + 1]; + + PT2sInGPU.lowerModuleIndices[Params_pT2::kLayers * PT2Index] = + segmentsInGPU.innerLowerModuleIndices[pixelSegmentIndex]; + PT2sInGPU.lowerModuleIndices[Params_pT2::kLayers * PT2Index + 1] = + segmentsInGPU.outerLowerModuleIndices[pixelSegmentIndex]; + PT2sInGPU.lowerModuleIndices[Params_pT2::kLayers * PT2Index + 2] = + segmentsInGPU.innerLowerModuleIndices[segmentIndex]; + PT2sInGPU.lowerModuleIndices[Params_pT2::kLayers * PT2Index + 3] = + segmentsInGPU.outerLowerModuleIndices[segmentIndex]; + + + unsigned int pixelInnerMD = segmentsInGPU.mdIndices[2 * pixelSegmentIndex]; + unsigned int pixelOuterMD = segmentsInGPU.mdIndices[2 * pixelSegmentIndex + 1]; + + PT2sInGPU.hitIndices[Params_pT2::kHits * PT2Index] = mdsInGPU.anchorHitIndices[pixelInnerMD]; + PT2sInGPU.hitIndices[Params_pT2::kHits * PT2Index + 1] = mdsInGPU.outerHitIndices[pixelInnerMD]; + PT2sInGPU.hitIndices[Params_pT2::kHits * PT2Index + 2] = mdsInGPU.anchorHitIndices[pixelOuterMD]; + PT2sInGPU.hitIndices[Params_pT2::kHits * PT2Index + 3] = mdsInGPU.outerHitIndices[pixelOuterMD]; + + PT2sInGPU.hitIndices[Params_pT2::kHits * PT2Index + 4] = + segmentsInGPU.innerMiniDoubletAnchorHitIndices[segmentIndex]; + PT2sInGPU.hitIndices[Params_pT2::kHits * PT2Index + 5] = + segmentsInGPU.innerMiniDoubletOuterHitIndices[segmentIndex]; + PT2sInGPU.hitIndices[Params_pT2::kHits * PT2Index + 6] = + segmentsInGPU.outerMiniDoubletAnchorHitIndices[segmentIndex]; + PT2sInGPU.hitIndices[Params_pT2::kHits * PT2Index + 7] = + segmentsInGPU.outerMiniDoubletOuterHitIndices[segmentIndex]; + PT2sInGPU.rzChiSquared[PT2Index] = rzChiSquared; + }; + + template + ALPAKA_FN_ACC ALPAKA_FN_INLINE bool runPT2DefaultAlgo(TAcc const& acc, + lst::Modules const& modulesInGPU, + lst::ObjectRanges const& rangesInGPU, + lst::MiniDoublets const& mdsInGPU, + lst::Segments const& segmentsInGPU, + unsigned int pixelSegmentIndex, + unsigned int segmentIndex, + float& pixelRadius, + float& centerX, + float& centerY, + float& rzChiSquared, + const float ptCut, + bool runChiSquaredCuts = true) { + //run pT4 compatibility between the pixel segment and inner segment, and between the pixel and outer segment of the triplet +// uint16_t pixelModuleIndex = segmentsInGPU.innerLowerModuleIndices[pixelSegmentIndex]; + +// uint16_t lowerModuleIndex = segmentsInGPU.innerLowerModuleIndices[segmentIndex]; +// uint16_t upperModuleIndex = segmentsInGPU.outerLowerModuleIndices[segmentIndex]; + + rzChiSquared = -1; + centerX = 0; + centerY = 0; + return true; + }; + + struct createPT2sInGPUFromMapv2 { + template + ALPAKA_FN_ACC void operator()(TAcc const& acc, + lst::Modules modulesInGPU, + lst::ObjectRanges rangesInGPU, + lst::MiniDoublets mdsInGPU, + lst::Segments segmentsInGPU, + lst::PT2s PT2sInGPU, + unsigned int* connectedPixelSize, + unsigned int* connectedPixelIndex, + unsigned int nPixelSegments, + const float ptCut) const { + auto const globalBlockIdx = alpaka::getIdx(acc); + auto const globalThreadIdx = alpaka::getIdx(acc); + auto const gridBlockExtent = alpaka::getWorkDiv(acc); + auto const gridThreadExtent = alpaka::getWorkDiv(acc); + + for (unsigned int i_pLS = globalThreadIdx[1]; i_pLS < nPixelSegments; i_pLS += gridThreadExtent[1]) { + auto iLSModule_max = connectedPixelIndex[i_pLS] + connectedPixelSize[i_pLS]; + + for (unsigned int iLSModule = connectedPixelIndex[i_pLS] + globalBlockIdx[0]; iLSModule < iLSModule_max; + iLSModule += gridBlockExtent[0]) { + uint16_t segmentLowerModuleIndex = + modulesInGPU + .connectedPixels[iLSModule]; //connected pixels will have the appropriate lower module index by default! +#ifdef WARNINGS + if (tripletLowerModuleIndex >= *modulesInGPU.nLowerModules) { + printf("tripletLowerModuleIndex %d >= modulesInGPU.nLowerModules %d \n", + tripletLowerModuleIndex, + *modulesInGPU.nLowerModules); + continue; //sanity check + } +#endif + //Removes 2S-2S :FIXME: filter these out in the pixel map + if (modulesInGPU.moduleType[segmentLowerModuleIndex] == lst::TwoS) + continue; + + uint16_t pixelModuleIndex = *modulesInGPU.nLowerModules; + unsigned int nOuterSegments = segmentsInGPU.nSegments[segmentLowerModuleIndex]; + if (nOuterSegments == 0) + continue; + + unsigned int pixelSegmentIndex = rangesInGPU.segmentModuleIndices[pixelModuleIndex] + i_pLS; + + if (segmentsInGPU.isDup[i_pLS]) + continue; + if (segmentsInGPU.partOfPT5[i_pLS]) + continue; //don't make pT3s for those pixels that are part of pT5 + if (segmentsInGPU.partOfPT3[i_pLS]) + continue; + + short layer2_adjustment; + if (modulesInGPU.layers[segmentLowerModuleIndex] == 1) { + layer2_adjustment = 1; + } //get upper segment to be in second layer + else if (modulesInGPU.layers[segmentLowerModuleIndex] == 2) { + layer2_adjustment = 0; + } // get lower segment to be in second layer + else { + continue; + } + + //fetch the segment + for (unsigned int outerSegmentArrayIndex = globalThreadIdx[2]; outerSegmentArrayIndex < nOuterSegments; + outerSegmentArrayIndex += gridThreadExtent[2]) { + unsigned int outerSegmentIndex = + rangesInGPU.segmentModuleIndices[segmentLowerModuleIndex] + outerSegmentArrayIndex; + if (modulesInGPU.moduleType[segmentsInGPU.outerLowerModuleIndices[outerSegmentIndex]] == lst::TwoS) + continue; //REMOVES PS-2S + + if (segmentsInGPU.partOfPT5[outerSegmentIndex]) + continue; //don't create pT2s for T2s accounted in pT5s + + if (segmentsInGPU.partOfPT3[outerSegmentIndex]) + continue; //don't create pT2s for T2s accounted in pT3s + + float pixelRadius, rzChiSquared, centerX, centerY; + bool success = runPT2DefaultAlgo(acc, + modulesInGPU, + rangesInGPU, + mdsInGPU, + segmentsInGPU, + pixelSegmentIndex, + outerSegmentIndex, + pixelRadius, + centerX, + centerY, + rzChiSquared, + ptCut); + + if (success) { + float phi = + mdsInGPU.anchorPhi[segmentsInGPU.mdIndices[outerSegmentIndex + + layer2_adjustment]]; + float eta = + mdsInGPU.anchorEta[segmentsInGPU.mdIndices[outerSegmentIndex + + layer2_adjustment]]; + float eta_pix = segmentsInGPU.eta[i_pLS]; + float phi_pix = segmentsInGPU.phi[i_pLS]; + float pt = segmentsInGPU.ptIn[i_pLS]; + unsigned int totOccupancyPT2s = + alpaka::atomicOp(acc, PT2sInGPU.totOccupancyPT2s, 1u); + if (totOccupancyPT2s >= n_max_pixel_segments) { +#ifdef WARNINGS + printf("Pixel Triplet excess alert!\n"); +#endif + } else { + unsigned int PT2Index = + alpaka::atomicOp(acc, PT2sInGPU.nPT2s, 1u); + addPT2ToMemory(mdsInGPU, + segmentsInGPU, + PT2sInGPU, + pixelSegmentIndex, + outerSegmentIndex, + pixelRadius, + centerX, + centerY, + rzChiSquared, + PT2Index, + pt, + eta, + phi, + eta_pix, + phi_pix); + } + } + } // for outerTripletArrayIndex + } // for iLSModule < iLSModule_max + } // for i_pLS + } + }; + + +} // namespace lst +#endif diff --git a/RecoTracker/LSTCore/src/alpaka/PixelQuintuplet.h b/RecoTracker/LSTCore/src/alpaka/PixelQuintuplet.h index f4bbcf04c69a8..90a71a0a44ef8 100644 --- a/RecoTracker/LSTCore/src/alpaka/PixelQuintuplet.h +++ b/RecoTracker/LSTCore/src/alpaka/PixelQuintuplet.h @@ -872,6 +872,7 @@ namespace lst { tripletsInGPU.partOfPT5[quintupletsInGPU.tripletIndices[2 * quintupletIndex + 1]] = true; segmentsInGPU.partOfPT5[i_pLS] = true; quintupletsInGPU.partOfPT5[quintupletIndex] = true; + segmentsInGPU.partOfPT5[i_pLS] = true; } // tot occupancy } // end success } // end T5 diff --git a/RecoTracker/LSTCore/src/alpaka/PixelTriplet.h b/RecoTracker/LSTCore/src/alpaka/PixelTriplet.h index 6234244a413d6..481e646d28065 100644 --- a/RecoTracker/LSTCore/src/alpaka/PixelTriplet.h +++ b/RecoTracker/LSTCore/src/alpaka/PixelTriplet.h @@ -1017,6 +1017,7 @@ namespace lst { phi_pix, score); tripletsInGPU.partOfPT3[outerTripletIndex] = true; + segmentsInGPU.partOfPT3[i_pLS] = true; } } } // for outerTripletArrayIndex diff --git a/RecoTracker/LSTCore/src/alpaka/Segment.h b/RecoTracker/LSTCore/src/alpaka/Segment.h index 22243a0e776e3..fe3b40261f62a 100644 --- a/RecoTracker/LSTCore/src/alpaka/Segment.h +++ b/RecoTracker/LSTCore/src/alpaka/Segment.h @@ -24,8 +24,11 @@ namespace lst { unsigned int* seedIdx; unsigned int* mdIndices; unsigned int* nMemoryLocations; + uint8_t* logicalLayers; unsigned int* innerMiniDoubletAnchorHitIndices; unsigned int* outerMiniDoubletAnchorHitIndices; + unsigned int* innerMiniDoubletOuterHitIndices; + unsigned int* outerMiniDoubletOuterHitIndices; int* charge; int* superbin; unsigned int* nSegments; //number of segments per inner lower module @@ -35,6 +38,7 @@ namespace lst { char* isQuad; char* isDup; bool* partOfPT5; + bool* partOfPT3; float* ptIn; float* ptErr; float* px; @@ -61,8 +65,11 @@ namespace lst { seedIdx = alpaka::getPtrNative(buf.seedIdx_buf); mdIndices = alpaka::getPtrNative(buf.mdIndices_buf); nMemoryLocations = alpaka::getPtrNative(buf.nMemoryLocations_buf); + logicalLayers = alpaka::getPtrNative(buf.logicalLayers_buf); innerMiniDoubletAnchorHitIndices = alpaka::getPtrNative(buf.innerMiniDoubletAnchorHitIndices_buf); outerMiniDoubletAnchorHitIndices = alpaka::getPtrNative(buf.outerMiniDoubletAnchorHitIndices_buf); + innerMiniDoubletOuterHitIndices = alpaka::getPtrNative(buf.innerMiniDoubletOuterHitIndices_buf); + outerMiniDoubletOuterHitIndices = alpaka::getPtrNative(buf.outerMiniDoubletOuterHitIndices_buf); charge = alpaka::getPtrNative(buf.charge_buf); superbin = alpaka::getPtrNative(buf.superbin_buf); nSegments = alpaka::getPtrNative(buf.nSegments_buf); @@ -72,6 +79,7 @@ namespace lst { isQuad = alpaka::getPtrNative(buf.isQuad_buf); isDup = alpaka::getPtrNative(buf.isDup_buf); partOfPT5 = alpaka::getPtrNative(buf.partOfPT5_buf); + partOfPT3 = alpaka::getPtrNative(buf.partOfPT3_buf); ptIn = alpaka::getPtrNative(buf.ptIn_buf); ptErr = alpaka::getPtrNative(buf.ptErr_buf); px = alpaka::getPtrNative(buf.px_buf); @@ -100,8 +108,11 @@ namespace lst { Buf seedIdx_buf; Buf mdIndices_buf; Buf nMemoryLocations_buf; + Buf logicalLayers_buf; Buf innerMiniDoubletAnchorHitIndices_buf; Buf outerMiniDoubletAnchorHitIndices_buf; + Buf innerMiniDoubletOuterHitIndices_buf; + Buf outerMiniDoubletOuterHitIndices_buf; Buf charge_buf; Buf superbin_buf; Buf nSegments_buf; @@ -111,6 +122,7 @@ namespace lst { Buf isQuad_buf; Buf isDup_buf; Buf partOfPT5_buf; + Buf partOfPT3_buf; Buf ptIn_buf; Buf ptErr_buf; Buf px_buf; @@ -143,8 +155,11 @@ namespace lst { seedIdx_buf(allocBufWrapper(devAccIn, maxPixelSegments, queue)), mdIndices_buf(allocBufWrapper(devAccIn, nMemoryLocationsIn * 2, queue)), nMemoryLocations_buf(allocBufWrapper(devAccIn, 1, queue)), + logicalLayers_buf(allocBufWrapper(devAccIn, nMemoryLocationsIn * Params_LS::kLayers, queue)), innerMiniDoubletAnchorHitIndices_buf(allocBufWrapper(devAccIn, nMemoryLocationsIn, queue)), outerMiniDoubletAnchorHitIndices_buf(allocBufWrapper(devAccIn, nMemoryLocationsIn, queue)), + innerMiniDoubletOuterHitIndices_buf(allocBufWrapper(devAccIn, nMemoryLocationsIn, queue)), + outerMiniDoubletOuterHitIndices_buf(allocBufWrapper(devAccIn, nMemoryLocationsIn, queue)), charge_buf(allocBufWrapper(devAccIn, maxPixelSegments, queue)), superbin_buf(allocBufWrapper(devAccIn, maxPixelSegments, queue)), nSegments_buf(allocBufWrapper(devAccIn, nLowerModules + 1, queue)), @@ -154,6 +169,7 @@ namespace lst { isQuad_buf(allocBufWrapper(devAccIn, maxPixelSegments, queue)), isDup_buf(allocBufWrapper(devAccIn, maxPixelSegments, queue)), partOfPT5_buf(allocBufWrapper(devAccIn, maxPixelSegments, queue)), + partOfPT3_buf(allocBufWrapper(devAccIn, maxPixelSegments, queue)), ptIn_buf(allocBufWrapper(devAccIn, maxPixelSegments, queue)), ptErr_buf(allocBufWrapper(devAccIn, maxPixelSegments, queue)), px_buf(allocBufWrapper(devAccIn, maxPixelSegments, queue)), @@ -169,6 +185,7 @@ namespace lst { alpaka::memset(queue, nSegments_buf, 0u); alpaka::memset(queue, totOccupancySegments_buf, 0u); alpaka::memset(queue, partOfPT5_buf, false); + alpaka::memset(queue, partOfPT3_buf, false); alpaka::memset(queue, pLSHitsIdxs_buf, 0u); alpaka::wait(queue); } @@ -340,12 +357,15 @@ namespace lst { }; ALPAKA_FN_ACC ALPAKA_FN_INLINE void addSegmentToMemory(lst::Segments& segmentsInGPU, + lst::Modules const& modulesInGPU, unsigned int lowerMDIndex, unsigned int upperMDIndex, uint16_t innerLowerModuleIndex, uint16_t outerLowerModuleIndex, unsigned int innerMDAnchorHitIndex, unsigned int outerMDAnchorHitIndex, + unsigned int innerMDOuterHitIndex, + unsigned int outerMDOuterHitIndex, float dPhi, float dPhiMin, float dPhiMax, @@ -359,6 +379,13 @@ namespace lst { segmentsInGPU.outerLowerModuleIndices[idx] = outerLowerModuleIndex; segmentsInGPU.innerMiniDoubletAnchorHitIndices[idx] = innerMDAnchorHitIndex; segmentsInGPU.outerMiniDoubletAnchorHitIndices[idx] = outerMDAnchorHitIndex; + segmentsInGPU.innerMiniDoubletOuterHitIndices[idx] = innerMDOuterHitIndex; + segmentsInGPU.outerMiniDoubletOuterHitIndices[idx] = outerMDOuterHitIndex; + + segmentsInGPU.logicalLayers[idx* Params_LS::kLayers] = + modulesInGPU.layers[innerLowerModuleIndex] + (modulesInGPU.subdets[innerLowerModuleIndex] == 4) * 6; + segmentsInGPU.logicalLayers[idx * Params_LS::kLayers + 1] = + modulesInGPU.layers[outerLowerModuleIndex] + (modulesInGPU.subdets[outerLowerModuleIndex] == 4) * 6; segmentsInGPU.dPhis[idx] = __F2H(dPhi); segmentsInGPU.dPhiMins[idx] = __F2H(dPhiMin); @@ -735,6 +762,8 @@ namespace lst { unsigned int innerMiniDoubletAnchorHitIndex = mdsInGPU.anchorHitIndices[innerMDIndex]; unsigned int outerMiniDoubletAnchorHitIndex = mdsInGPU.anchorHitIndices[outerMDIndex]; + unsigned int innerMiniDoubletOuterHitIndex = mdsInGPU.outerHitIndices[innerMDIndex]; + unsigned int outerMiniDoubletOuterHitIndex = mdsInGPU.outerHitIndices[outerMDIndex]; dPhiMin = 0; dPhiMax = 0; dPhiChangeMin = 0; @@ -767,12 +796,15 @@ namespace lst { unsigned int segmentIdx = rangesInGPU.segmentModuleIndices[innerLowerModuleIndex] + segmentModuleIdx; addSegmentToMemory(segmentsInGPU, + modulesInGPU, innerMDIndex, outerMDIndex, innerLowerModuleIndex, outerLowerModuleIndex, innerMiniDoubletAnchorHitIndex, outerMiniDoubletAnchorHitIndex, + innerMiniDoubletOuterHitIndex, + outerMiniDoubletOuterHitIndex, dPhi, dPhiMin, dPhiMax, diff --git a/RecoTracker/LSTCore/src/alpaka/TrackCandidate.h b/RecoTracker/LSTCore/src/alpaka/TrackCandidate.h index ede4dd9471e8e..17f7e770afaa5 100644 --- a/RecoTracker/LSTCore/src/alpaka/TrackCandidate.h +++ b/RecoTracker/LSTCore/src/alpaka/TrackCandidate.h @@ -9,15 +9,17 @@ #include "MiniDoublet.h" #include "PixelTriplet.h" #include "Quintuplet.h" +#include "PT2.h" #include "Hit.h" #include "ObjectRanges.h" namespace lst { struct TrackCandidates { - short* trackCandidateType; // 4-T5 5-pT3 7-pT5 8-pLS + short* trackCandidateType; // 4-T5 5-pT3 7-pT5 8-pLS 9-T4 10-PT2 unsigned int* directObjectIndices; // Will hold direct indices to each type containers unsigned int* objectIndices; // Will hold tracklet and triplet indices - check the type!! unsigned int* nTrackCandidates; + unsigned int* nTrackCandidatesPT2; unsigned int* nTrackCandidatespT3; unsigned int* nTrackCandidatespT5; unsigned int* nTrackCandidatespLS; @@ -38,6 +40,7 @@ namespace lst { directObjectIndices = alpaka::getPtrNative(buf.directObjectIndices_buf); objectIndices = alpaka::getPtrNative(buf.objectIndices_buf); nTrackCandidates = alpaka::getPtrNative(buf.nTrackCandidates_buf); + nTrackCandidatesPT2 = alpaka::getPtrNative(buf.nTrackCandidatesPT2_buf); nTrackCandidatespT3 = alpaka::getPtrNative(buf.nTrackCandidatespT3_buf); nTrackCandidatespT5 = alpaka::getPtrNative(buf.nTrackCandidatespT5_buf); nTrackCandidatespLS = alpaka::getPtrNative(buf.nTrackCandidatespLS_buf); @@ -60,6 +63,7 @@ namespace lst { Buf directObjectIndices_buf; Buf objectIndices_buf; Buf nTrackCandidates_buf; + Buf nTrackCandidatesPT2_buf; Buf nTrackCandidatespT3_buf; Buf nTrackCandidatespT5_buf; Buf nTrackCandidatespLS_buf; @@ -82,6 +86,7 @@ namespace lst { directObjectIndices_buf(allocBufWrapper(devAccIn, maxTrackCandidates, queue)), objectIndices_buf(allocBufWrapper(devAccIn, 2 * maxTrackCandidates, queue)), nTrackCandidates_buf(allocBufWrapper(devAccIn, 1, queue)), + nTrackCandidatesPT2_buf(allocBufWrapper(devAccIn, 1, queue)), nTrackCandidatespT3_buf(allocBufWrapper(devAccIn, 1, queue)), nTrackCandidatespT5_buf(allocBufWrapper(devAccIn, 1, queue)), nTrackCandidatespLS_buf(allocBufWrapper(devAccIn, 1, queue)), @@ -95,6 +100,7 @@ namespace lst { radius_buf(allocBufWrapper(devAccIn, maxTrackCandidates, queue)) { alpaka::memset(queue, nTrackCandidates_buf, 0u); alpaka::memset(queue, nTrackCandidatesT5_buf, 0u); + alpaka::memset(queue, nTrackCandidatesPT2_buf, 0u); alpaka::memset(queue, nTrackCandidatespT3_buf, 0u); alpaka::memset(queue, nTrackCandidatespT5_buf, 0u); alpaka::memset(queue, nTrackCandidatespLS_buf, 0u); @@ -148,9 +154,10 @@ namespace lst { trackCandidatesInGPU.objectIndices[2 * trackCandidateIndex] = innerTrackletIndex; trackCandidatesInGPU.objectIndices[2 * trackCandidateIndex + 1] = outerTrackletIndex; - size_t limits = trackCandidateType == 7 - ? Params_pT5::kLayers - : Params_pT3::kLayers; // 7 means pT5, Params_pT3::kLayers = Params_T5::kLayers = 5 + size_t limits = trackCandidateType == 7 ? Params_pT5::kLayers : (trackCandidateType == 10 ? Params_pT2::kLayers : Params_pT3::kLayers); +// size_t limits = trackCandidateType == 7 +// ? Params_pT5::kLayers +// : Params_pT3::kLayers; // 7 means pT5, Params_pT3::kLayers = Params_T5::kLayers = 5 //send the starting pointer to the logicalLayer and hitIndices for (size_t i = 0; i < limits; i++) { @@ -434,6 +441,57 @@ namespace lst { } }; + + struct addPT2asTrackCandidatesInGPU { + template + ALPAKA_FN_ACC void operator()(TAcc const& acc, + uint16_t nLowerModules, + lst::PT2s PT2sInGPU, + lst::TrackCandidates trackCandidatesInGPU, + lst::Segments segmentsInGPU, + lst::ObjectRanges rangesInGPU) const { + auto const globalThreadIdx = alpaka::getIdx(acc); + auto const gridThreadExtent = alpaka::getWorkDiv(acc); + + unsigned int nPT2s = *PT2sInGPU.nPT2s; + unsigned int pLS_offset = rangesInGPU.segmentModuleIndices[nLowerModules]; + for (unsigned int PT2Index = globalThreadIdx[2]; PT2Index < nPT2s; + PT2Index += gridThreadExtent[2]) { + if ((PT2sInGPU.isDup[PT2Index])) + continue; + + unsigned int trackCandidateIdx = + alpaka::atomicOp(acc, trackCandidatesInGPU.nTrackCandidates, 1u); +// if (trackCandidateIdx >= n_max_pixel_track_candidates) // This is done before any non-pixel TCs are added +// { +//#ifdef WARNINGS +// printf("Track Candidate excess alert! Type = PT2"); +//#endif +// alpaka::atomicOp(acc, trackCandidatesInGPU.nTrackCandidates, 1u); +// break; + +// } else { + alpaka::atomicOp(acc, trackCandidatesInGPU.nTrackCandidatesPT2, 1u); + + float radius = __H2F(PT2sInGPU.pixelRadius[PT2Index]); + unsigned int PT2PixelIndex = PT2sInGPU.pixelSegmentIndices[PT2Index]; + addTrackCandidateToMemory(trackCandidatesInGPU, + 10 /*track candidate type PT2=10*/, + PT2Index, + PT2Index, + &PT2sInGPU.logicalLayers[Params_pT2::kLayers * PT2Index], + &PT2sInGPU.lowerModuleIndices[Params_pT2::kLayers * PT2Index], + &PT2sInGPU.hitIndices[Params_pT2::kHits * PT2Index], + segmentsInGPU.seedIdx[PT2PixelIndex - pLS_offset], + __H2F(PT2sInGPU.centerX[PT2Index]), + __H2F(PT2sInGPU.centerY[PT2Index]), + radius, + trackCandidateIdx, + PT2Index); +// } + } + } + }; struct addT5asTrackCandidateInGPU { template ALPAKA_FN_ACC void operator()(TAcc const& acc, diff --git a/RecoTracker/LSTCore/standalone/bin/lst.cc b/RecoTracker/LSTCore/standalone/bin/lst.cc index 8cf0681d986c2..2908495e46a0a 100644 --- a/RecoTracker/LSTCore/standalone/bin/lst.cc +++ b/RecoTracker/LSTCore/standalone/bin/lst.cc @@ -414,6 +414,7 @@ void run_lst() { float timing_pLS; float timing_pT5; float timing_pT3; + float timing_PT2; float timing_TC; #pragma omp for // nowait// private(event) @@ -454,6 +455,7 @@ void run_lst() { timing_pLS = runPixelLineSegment(events.at(omp_get_thread_num()), ana.no_pls_dupclean); timing_pT5 = runPixelQuintuplet(events.at(omp_get_thread_num())); timing_pT3 = runpT3(events.at(omp_get_thread_num())); + timing_PT2 = runPT2(events.at(omp_get_thread_num())); timing_TC = runTrackCandidate(events.at(omp_get_thread_num()), ana.no_pls_dupclean, ana.tc_pls_triplets); if (ana.verbose == 4) { @@ -497,6 +499,7 @@ void run_lst() { timing_pLS, timing_pT5, timing_pT3, + timing_PT2, timing_TC, timing_resetEvent}); } diff --git a/RecoTracker/LSTCore/standalone/code/core/AccessHelper.cc b/RecoTracker/LSTCore/standalone/code/core/AccessHelper.cc index 76cfa9760b71a..7af9fbc0f6476 100644 --- a/RecoTracker/LSTCore/standalone/code/core/AccessHelper.cc +++ b/RecoTracker/LSTCore/standalone/code/core/AccessHelper.cc @@ -219,6 +219,89 @@ std::tuple, std::vector> getHitIdxsAndHi return convertHitsToHitIdxsAndHitTypes(event, getHitsFromT5(event, T5)); } +// =============== +// ----* pT2 *---- +// =============== + +//____________________________________________________________________________________________ +unsigned int getPixelLSFrompT2(lst::Event* event, unsigned int PT2) { + lst::PT2s const* PT2s = event->getPT2s()->data(); + lst::ObjectRanges const* rangesEvt = event->getRanges()->data(); + lst::Modules const* modulesEvt = event->getModules()->data(); + const unsigned int pLS_offset = rangesEvt->segmentModuleIndices[*(modulesEvt->nLowerModules)]; + return PT2s->pixelSegmentIndices[PT2] - pLS_offset; +} + +//____________________________________________________________________________________________ +unsigned int getLSsFrompT2(lst::Event* event, unsigned int PT2) { + lst::PT2s const* PT2s = event->getPT2s()->data(); + return PT2s->segmentIndices[PT2]; +} + +//____________________________________________________________________________________________ +std::vector getMDsFrompT2(lst::Event* event, unsigned int PT2) { + unsigned int LS = getLSsFrompT2(event, PT2); + return getMDsFromLS(event, LS); +} + +//____________________________________________________________________________________________ +std::vector getOuterTrackerHitsFrompT2(lst::Event* event, unsigned int PT2) { + unsigned int LS = getLSsFrompT2(event, PT2); + return getHitsFromLS(event, LS); +} + +//____________________________________________________________________________________________ +std::vector getPixelHitsFrompT2(lst::Event* event, unsigned int PT2) { + unsigned int pLS = getPixelLSFrompT2(event, PT2); + return getPixelHitsFrompLS(event, pLS); +} + +//____________________________________________________________________________________________ +std::vector getHitsFrompT2(lst::Event* event, unsigned int PT2) { + unsigned int pLS = getPixelLSFrompT2(event, PT2); + unsigned int LS = getLSsFrompT2(event, PT2); + std::vector pixelHits = getPixelHitsFrompLS(event, pLS); + std::vector outerTrackerHits = getHitsFromLS(event, LS); + pixelHits.insert(pixelHits.end(), outerTrackerHits.begin(), outerTrackerHits.end()); + return pixelHits; +} + +//____________________________________________________________________________________________ +std::vector getHitIdxsFrompT2(lst::Event* event, unsigned int PT2) { + lst::Hits const* hitsEvt = event->getHits()->data(); + std::vector hits = getHitsFrompT2(event, PT2); + std::vector hitidxs; + for (auto& hit : hits) + hitidxs.push_back(hitsEvt->idxs[hit]); + return hitidxs; +} +//____________________________________________________________________________________________ +std::vector getModuleIdxsFrompT2(lst::Event* event, unsigned int PT2) { + std::vector hits = getOuterTrackerHitsFrompT2(event, PT2); + std::vector module_idxs; + lst::Hits const* hitsEvt = event->getHits()->data(); + for (auto& hitIdx : hits) { + module_idxs.push_back(hitsEvt->moduleIndices[hitIdx]); + } + return module_idxs; +} +//____________________________________________________________________________________________ +std::vector getHitTypesFrompT2(lst::Event* event, unsigned int PT2) { + unsigned int pLS = getPixelLSFrompT2(event, PT2); + std::vector pixelHits = getPixelHitsFrompLS(event, pLS); + // pixel Hits list will be either 3 or 4 and depending on it return accordingly + if (pixelHits.size() == 3) + return {0, 0, 0, 4, 4, 4, 4}; + else + return {0, 0, 0, 0, 4, 4, 4, 4}; +} + +//____________________________________________________________________________________________ +std::tuple, std::vector> getHitIdxsAndHitTypesFrompT2(lst::Event* event, + unsigned PT2) { + return convertHitsToHitIdxsAndHitTypes(event, getHitsFrompT2(event, PT2)); +} + // =============== // ----* pT3 *---- // =============== @@ -422,6 +505,12 @@ std::vector getLSsFromTC(lst::Event* event, unsigned int TC case kpT3: return getLSsFrompT3(event, objidx); break; + case kpT2: + { + std::vector PT2SegmentIndexVector = {getLSsFrompT2(event, objidx)}; + return PT2SegmentIndexVector; + break; + } case kT5: return getLSsFromT5(event, objidx); break; @@ -445,6 +534,9 @@ std::tuple, std::vector> getHitIdxsAndHi case kpT3: return getHitIdxsAndHitTypesFrompT3(event, objidx); break; + case kpT2: + return getHitIdxsAndHitTypesFrompT2(event, objidx); + break; case kT5: return getHitIdxsAndHitTypesFromT5(event, objidx); break; diff --git a/RecoTracker/LSTCore/standalone/code/core/AccessHelper.h b/RecoTracker/LSTCore/standalone/code/core/AccessHelper.h index d0924518eeb4d..a23664633d958 100644 --- a/RecoTracker/LSTCore/standalone/code/core/AccessHelper.h +++ b/RecoTracker/LSTCore/standalone/code/core/AccessHelper.h @@ -7,7 +7,7 @@ using LSTEvent = lst::Event; -enum { kpT5 = 7, kpT3 = 5, kT5 = 4, kpLS = 8 }; +enum { kpT5 = 7, kpT3 = 5, kT5 = 4, kpLS = 8, kpT2 = 10 }; // ----* Hit *---- std::tuple, std::vector> convertHitsToHitIdxsAndHitTypes( @@ -78,6 +78,19 @@ std::vector getModuleIdxsFrompT5(LSTEvent* event, unsigned int pT5 std::tuple, std::vector> getHitIdxsAndHitTypesFrompT5(LSTEvent* event, unsigned pT5); +// ----* pT2 *---- +unsigned int getPixelLSFrompT2(LSTEvent* event, unsigned int pT2); +unsigned int getT2FrompT2(LSTEvent* event, unsigned int pT2); +unsigned int getLSsFrompT2(LSTEvent* event, unsigned int pT2); +std::vector getMDsFrompT2(LSTEvent* event, unsigned int pT2); +std::vector getPixelHitsFrompT2(LSTEvent* event, unsigned int pT2); +std::vector getHitsFrompT2(LSTEvent* event, unsigned int pT2); +std::vector getHitIdxsFrompT2(LSTEvent* event, unsigned int pT2); +std::vector getHitTypesFrompT2(LSTEvent* event, unsigned int pT2); +std::vector getModuleIdxsFrompT2(LSTEvent* event, unsigned int pT2); +std::tuple, std::vector> getHitIdxsAndHitTypesFrompT2(LSTEvent* event, + unsigned pT2); + // ----* TC *---- std::vector getLSsFromTC(LSTEvent* event, unsigned int TC); std::vector getHitsFromTC(LSTEvent* event, unsigned int TC); diff --git a/RecoTracker/LSTCore/standalone/code/core/trkCore.cc b/RecoTracker/LSTCore/standalone/code/core/trkCore.cc index bcda760ae7c72..d5189d4a7ed88 100644 --- a/RecoTracker/LSTCore/standalone/code/core/trkCore.cc +++ b/RecoTracker/LSTCore/standalone/code/core/trkCore.cc @@ -242,6 +242,22 @@ float runPixelQuintuplet(lst::Event *event) { return pt5_elapsed; } +//___________________________________________________________________________________________________________________________________________________________________________________________ +float runPT2(lst::Event *event) { + TStopwatch my_timer; + if (ana.verbose >= 2) + std::cout << "Reco PT2 start" << std::endl; + my_timer.Start(); + event->createPT2s(); + float PT2_elapsed = my_timer.RealTime(); + if (ana.verbose >= 2) + std::cout << "Reco PT2 processing time: " << PT2_elapsed << " secs" << std::endl; + if (ana.verbose >= 2) + std::cout << "# of PT2s produced: " << event->getNumberOfPT2s() << std::endl; + + return PT2_elapsed; +} + //___________________________________________________________________________________________________________________________________________________________________________________________ float runTrackCandidate(lst::Event *event, bool no_pls_dupclean, bool tc_pls_triplets) { TStopwatch my_timer; @@ -261,6 +277,8 @@ float runTrackCandidate(lst::Event *event, bool no_pls_dupclean, bool tc_ std::cout << " # of pT5 TrackCandidates produced: " << event->getNumberOfPT5TrackCandidates() << std::endl; if (ana.verbose >= 2) std::cout << " # of pT3 TrackCandidates produced: " << event->getNumberOfPT3TrackCandidates() << std::endl; + if (ana.verbose >= 2) + std::cout << " # of pT2 TrackCandidates produced: " << event->getNumberOfPT2TrackCandidates() << std::endl; if (ana.verbose >= 2) std::cout << " # of pLS TrackCandidates produced: " << event->getNumberOfPLSTrackCandidates() << std::endl; if (ana.verbose >= 2) @@ -956,8 +974,9 @@ void printTimingInformation(std::vector> &timing_information, timing_total_short += timing[4] * 1000; // T5 timing_total_short += timing[6] * 1000; // pT5 timing_total_short += timing[7] * 1000; // pT3 - timing_total_short += timing[8] * 1000; // TC - timing_total_short += timing[9] * 1000; // Reset + timing_total_short += timing[8] * 1000; // pT2 + timing_total_short += timing[9] * 1000; // TC + timing_total_short += timing[10] * 1000; // Reset std::cout << std::setw(6) << ievt; std::cout << " " << std::setw(6) << timing[0] * 1000; // Hits std::cout << " " << std::setw(6) << timing[1] * 1000; // MD @@ -967,8 +986,9 @@ void printTimingInformation(std::vector> &timing_information, std::cout << " " << std::setw(6) << timing[5] * 1000; // pLS std::cout << " " << std::setw(6) << timing[6] * 1000; // pT5 std::cout << " " << std::setw(6) << timing[7] * 1000; // pT3 - std::cout << " " << std::setw(6) << timing[8] * 1000; // TC - std::cout << " " << std::setw(6) << timing[9] * 1000; // Reset + std::cout << " " << std::setw(6) << timing[8] * 1000; // pT2 + std::cout << " " << std::setw(6) << timing[9] * 1000; // TC + std::cout << " " << std::setw(6) << timing[10] * 1000; // Reset std::cout << " " << std::setw(7) << timing_total; // Total time std::cout << " " << std::setw(7) << timing_total_short; // Total time std::cout << std::endl; @@ -980,8 +1000,9 @@ void printTimingInformation(std::vector> &timing_information, timing_sum_information[5] += timing[5] * 1000; // pLS timing_sum_information[6] += timing[6] * 1000; // pT5 timing_sum_information[7] += timing[7] * 1000; // pT3 - timing_sum_information[8] += timing[8] * 1000; // TC - timing_sum_information[9] += timing[9] * 1000; // Reset + timing_sum_information[8] += timing[8] * 1000; // pT2 + timing_sum_information[9] += timing[9] * 1000; // TC + timing_sum_information[10] += timing[10] * 1000; // Reset timing_shortlist.push_back(timing_total_short); // short total timing_list.push_back(timing_total); // short total } @@ -993,8 +1014,9 @@ void printTimingInformation(std::vector> &timing_information, timing_sum_information[5] /= timing_information.size(); // pLS timing_sum_information[6] /= timing_information.size(); // pT5 timing_sum_information[7] /= timing_information.size(); // pT3 - timing_sum_information[8] /= timing_information.size(); // TC - timing_sum_information[9] /= timing_information.size(); // Reset + timing_sum_information[8] /= timing_information.size(); // pT2 + timing_sum_information[9] /= timing_information.size(); // TC + timing_sum_information[10] /= timing_information.size(); // Reset float timing_total_avg = 0.0; timing_total_avg += timing_sum_information[0]; // Hits @@ -1005,8 +1027,9 @@ void printTimingInformation(std::vector> &timing_information, timing_total_avg += timing_sum_information[5]; // pLS timing_total_avg += timing_sum_information[6]; // pT5 timing_total_avg += timing_sum_information[7]; // pT3 - timing_total_avg += timing_sum_information[8]; // TC - timing_total_avg += timing_sum_information[9]; // Reset + timing_total_avg += timing_sum_information[8]; // pT2 + timing_total_avg += timing_sum_information[9]; // TC + timing_total_avg += timing_sum_information[10]; // Reset float timing_totalshort_avg = 0.0; timing_totalshort_avg += timing_sum_information[1]; // MD timing_totalshort_avg += timing_sum_information[2]; // LS @@ -1014,8 +1037,9 @@ void printTimingInformation(std::vector> &timing_information, timing_totalshort_avg += timing_sum_information[4]; // T5 timing_totalshort_avg += timing_sum_information[6]; // pT5 timing_totalshort_avg += timing_sum_information[7]; // pT3 - timing_totalshort_avg += timing_sum_information[8]; // TC - timing_totalshort_avg += timing_sum_information[9]; // Reset + timing_totalshort_avg += timing_sum_information[8]; // pT2 + timing_totalshort_avg += timing_sum_information[9]; // TC + timing_totalshort_avg += timing_sum_information[10]; // Reset float standardDeviation = 0.0; for (auto shorttime : timing_shortlist) { @@ -1033,6 +1057,7 @@ void printTimingInformation(std::vector> &timing_information, std::cout << " " << std::setw(6) << "pLS"; std::cout << " " << std::setw(6) << "pT5"; std::cout << " " << std::setw(6) << "pT3"; + std::cout << " " << std::setw(6) << "pT2"; std::cout << " " << std::setw(6) << "TC"; std::cout << " " << std::setw(6) << "Reset"; std::cout << " " << std::setw(7) << "Total"; @@ -1047,8 +1072,9 @@ void printTimingInformation(std::vector> &timing_information, std::cout << " " << std::setw(6) << timing_sum_information[5]; // pLS std::cout << " " << std::setw(6) << timing_sum_information[6]; // pT5 std::cout << " " << std::setw(6) << timing_sum_information[7]; // pT3 - std::cout << " " << std::setw(6) << timing_sum_information[8]; // TC - std::cout << " " << std::setw(6) << timing_sum_information[9]; // Reset + std::cout << " " << std::setw(6) << timing_sum_information[8]; // pT2 + std::cout << " " << std::setw(6) << timing_sum_information[9]; // TC + std::cout << " " << std::setw(6) << timing_sum_information[10]; // Reset std::cout << " " << std::setw(7) << timing_total_avg; // Average total time std::cout << " " << std::setw(7) << timing_totalshort_avg; // Average total time std::cout << "+/- " << std::setw(4) << stdDev; diff --git a/RecoTracker/LSTCore/standalone/code/core/trkCore.h b/RecoTracker/LSTCore/standalone/code/core/trkCore.h index ad12232eb1678..2c37fe1b9ce6f 100644 --- a/RecoTracker/LSTCore/standalone/code/core/trkCore.h +++ b/RecoTracker/LSTCore/standalone/code/core/trkCore.h @@ -27,6 +27,7 @@ float runQuintuplet(LSTEvent *event); float runPixelQuintuplet(LSTEvent *event); float runPixelLineSegment(LSTEvent *event, bool no_pls_dupclean); float runpT3(LSTEvent *event); +float runPT2(LSTEvent *event); // --------------------- ======================== --------------------- diff --git a/RecoTracker/LSTCore/standalone/code/core/write_lst_ntuple.cc b/RecoTracker/LSTCore/standalone/code/core/write_lst_ntuple.cc index 2108cd0aca26d..84507da889379 100644 --- a/RecoTracker/LSTCore/standalone/code/core/write_lst_ntuple.cc +++ b/RecoTracker/LSTCore/standalone/code/core/write_lst_ntuple.cc @@ -102,6 +102,24 @@ void createOptionalOutputBranches() { ana.tx->createBranch>("pT3_layer_binary"); ana.tx->createBranch>("pT3_moduleType_binary"); + // PT2 branches + ana.tx->createBranch>("sim_pT2_matched"); + ana.tx->createBranch>("pT2_isFake"); + ana.tx->createBranch>("pT2_isDup_reco_test"); + ana.tx->createBranch>("pT2_isDuplicate"); + ana.tx->createBranch>("pT2_pt"); + ana.tx->createBranch>("pT2_eta"); + ana.tx->createBranch>("pT2_phi"); + ana.tx->createBranch>("pT2_layer_binary"); + ana.tx->createBranch>("pT2_moduleType_binary"); + + ana.tx->createBranch>>("pT2_matched_simIdx"); +// ana.tx->createBranch>>("pT2_hitIdxs"); +// ana.tx->createBranch>("pT2_pixelRadius"); +// ana.tx->createBranch>("pT2_pixelRadiusError"); +// ana.tx->createBranch>>("pT2_matched_pt"); +// ana.tx->createBranch>("pT2_rzChiSquared"); + // pLS branches ana.tx->createBranch>("sim_pLS_matched"); ana.tx->createBranch>>("sim_pLS_types"); @@ -150,6 +168,7 @@ void createOptionalOutputBranches() { ana.tx->createBranch>("t3_occupancies"); ana.tx->createBranch("tc_occupancies"); ana.tx->createBranch>("t5_occupancies"); + ana.tx->createBranch("pT2_occupancies"); ana.tx->createBranch("pT3_occupancies"); ana.tx->createBranch("pT5_occupancies"); @@ -333,6 +352,7 @@ void setOptionalOutputBranches(lst::Event* event) { setPixelQuintupletOutputBranches(event); setQuintupletOutputBranches(event); setPixelTripletOutputBranches(event); + setPT2OutputBranches(event); setOccupancyBranches(event); setT5DNNBranches(event); @@ -348,6 +368,7 @@ void setOccupancyBranches(lst::Event* event) { lst::Quintuplets const& quintupletsInGPU = (*event->getQuintuplets()->data()); lst::PixelQuintuplets const& pixelQuintupletsInGPU = (*event->getPixelQuintuplets()->data()); lst::PixelTriplets const& pixelTripletsInGPU = (*event->getPixelTriplets()->data()); + lst::PT2s const& PT2sInGPU = (*event->getPT2s()->data()); lst::TrackCandidates const& trackCandidatesInGPU = (*event->getTrackCandidates()->data()); std::vector moduleLayer; @@ -396,6 +417,7 @@ void setOccupancyBranches(lst::Event* event) { ana.tx->setBranch>("sg_occupancies", segmentOccupancy); ana.tx->setBranch>("t3_occupancies", tripletOccupancy); ana.tx->setBranch("tc_occupancies", *(trackCandidatesInGPU.nTrackCandidates)); + ana.tx->setBranch("pT2_occupancies", *(PT2sInGPU.totOccupancyPT2s)); ana.tx->setBranch("pT3_occupancies", *(pixelTripletsInGPU.totOccupancyPixelTriplets)); ana.tx->setBranch>("t5_occupancies", quintupletOccupancy); ana.tx->setBranch("pT5_occupancies", *(pixelQuintupletsInGPU.totOccupancyPixelQuintuplets)); @@ -617,6 +639,69 @@ void setPixelTripletOutputBranches(lst::Event* event) { ana.tx->setBranch>("pT3_isDuplicate", pT3_isDuplicate); } +//________________________________________________________________________________________________________________________________ +void setPT2OutputBranches(lst::Event* event) { + lst::PT2s const* PT2s = event->getPT2s()->data(); + lst::Modules const* modules = event->getModules()->data(); + lst::Segments const* segments = event->getSegments()->data(); + int n_accepted_simtrk = ana.tx->getBranch>("sim_TC_matched").size(); + + unsigned int nPT2s = *PT2s->nPT2s; + std::vector sim_pT2_matched(n_accepted_simtrk); + std::vector> pT2_matched_simIdx; + + for (unsigned int PT2 = 0; PT2 < nPT2s; PT2++) { + //unsigned int LSIndex = getLSsFrompT2(event, PT2); + unsigned int pLSIndex = getPixelLSFrompT2(event, PT2); + const float pt = segments->ptIn[pLSIndex]; + + float eta = segments->eta[pLSIndex]; + float phi = segments->phi[pLSIndex]; + std::vector hit_idx = getHitIdxsFrompT2(event, PT2); + std::vector hit_type = getHitTypesFrompT2(event, PT2); + + std::vector simidx = matchedSimTrkIdxs(hit_idx, hit_type); + std::vector module_idx = getModuleIdxsFrompT2(event, PT2); + int layer_binary = 1; + int moduleType_binary = 0; + for (size_t i = 0; i < module_idx.size(); i += 2) { + layer_binary |= (1 << (modules->layers[module_idx[i]] + 6 * (modules->subdets[module_idx[i]] == 4))); + moduleType_binary |= (modules->moduleType[module_idx[i]] << i); + } + ana.tx->pushbackToBranch("pT2_isFake", static_cast(simidx.size() == 0)); + ana.tx->pushbackToBranch("pT2_pt", pt); + ana.tx->pushbackToBranch("pT2_eta", eta); + ana.tx->pushbackToBranch("pT2_phi", phi); + ana.tx->pushbackToBranch("pT2_layer_binary", layer_binary); + ana.tx->pushbackToBranch("pT2_moduleType_binary", moduleType_binary); + + pT2_matched_simIdx.push_back(simidx); + + for (auto& idx : simidx) { + if (idx < n_accepted_simtrk) { + sim_pT2_matched.at(idx) += 1; + } + } + } + + std::vector pT2_isDuplicate(pT2_matched_simIdx.size()); + for (unsigned int i = 0; i < pT2_matched_simIdx.size(); i++) { + bool isDuplicate = true; + for (unsigned int isim = 0; isim < pT2_matched_simIdx[i].size(); isim++) { + int simidx = pT2_matched_simIdx[i][isim]; + if (simidx < n_accepted_simtrk) { + if (sim_pT2_matched[simidx] > 1) { + isDuplicate = true; + } + } + } + pT2_isDuplicate[i] = isDuplicate; + } + ana.tx->setBranch>("sim_pT2_matched", sim_pT2_matched); + ana.tx->setBranch>>("pT2_matched_simIdx", pT2_matched_simIdx); + ana.tx->setBranch>("pT2_isDuplicate", pT2_isDuplicate); +} + //________________________________________________________________________________________________________________________________ void fillT5DNNBranches(lst::Event* event, unsigned int iT3) { lst::Hits const* hits = event->getHits()->data(); @@ -936,7 +1021,7 @@ std::tuple> parseTrackCandidate( lst::TrackCandidates const* trackCandidates = event->getTrackCandidates()->data(); short type = trackCandidates->trackCandidateType[idx]; - enum { pT5 = 7, pT3 = 5, T5 = 4, pLS = 8 }; + enum { pT5 = 7, pT3 = 5, T5 = 4, pLS = 8, pT2 = 10 }; // Compute pt eta phi and hit indices that will be used to figure out whether the TC matched float pt, eta, phi; @@ -948,6 +1033,9 @@ std::tuple> parseTrackCandidate( case pT3: std::tie(pt, eta, phi, hit_idx, hit_type) = parsepT3(event, idx); break; + case pT2: + std::tie(pt, eta, phi, hit_idx, hit_type) = parsepT2(event, idx); + break; case T5: std::tie(pt, eta, phi, hit_idx, hit_type) = parseT5(event, idx); break; @@ -1110,6 +1198,38 @@ std::tuple, std::vector, std::vector> parsepT2(lst::Event* event, + unsigned int idx) { + // Get relevant information + lst::TrackCandidates const* trackCandidates = event->getTrackCandidates()->data(); + lst::Segments const* segments = event->getSegments()->data(); + + // + // pictorial representation of a pT3 + // + // inner tracker outer tracker + // ------------- -------------------------- + // pLS 01 23 45 (anchor hit of a minidoublet is always the first of the pair) + // **** oo -- oo -- oo pT3 + unsigned int PT2 = trackCandidates->directObjectIndices[idx]; + unsigned int pLS = getPixelLSFrompT2(event, PT2); + //unsigned int LS = getLSFromPT2(event, PT2); + + // pixel pt + const float pt_pLS = segments->ptIn[pLS]; + const float eta_pLS = segments->eta[pLS]; + const float phi_pLS = segments->phi[pLS]; + // average pt + const float pt = pt_pLS ; + + // Form the hit idx/type std::vector + std::vector hit_idx = getHitIdxsFrompT2(event, PT2); + std::vector hit_type = getHitTypesFrompT2(event, PT2); + + return {pt, eta_pLS, phi_pLS, hit_idx, hit_type}; +} + //________________________________________________________________________________________________________________________________ std::tuple, std::vector> parseT5(lst::Event* event, unsigned int idx) { diff --git a/RecoTracker/LSTCore/standalone/code/core/write_lst_ntuple.h b/RecoTracker/LSTCore/standalone/code/core/write_lst_ntuple.h index 64f166ab997f8..0d5469e183801 100644 --- a/RecoTracker/LSTCore/standalone/code/core/write_lst_ntuple.h +++ b/RecoTracker/LSTCore/standalone/code/core/write_lst_ntuple.h @@ -27,6 +27,7 @@ void setOccupancyBranches(LSTEvent* event); void setPixelQuintupletOutputBranches(LSTEvent* event); void setQuintupletOutputBranches(LSTEvent* event); void setPixelTripletOutputBranches(LSTEvent* event); +void setPT2OutputBranches(LSTEvent* event); void setGnnNtupleBranches(LSTEvent* event); void setGnnNtupleMiniDoublet(LSTEvent* event, unsigned int MD); void fillT5DNNBranches(LSTEvent* event, unsigned int T3); @@ -41,6 +42,8 @@ std::tuple, std::vector, std::vector> parsepLS(LSTEvent* event, unsigned int); +std::tuple, std::vector> parsepT2(LSTEvent* event, + unsigned int); // Print multiplicities void printMiniDoubletMultiplicities(LSTEvent* event); @@ -54,6 +57,7 @@ void printLSs(LSTEvent* event); void printpLSs(LSTEvent* event); void printT3s(LSTEvent* event); void printT4s(LSTEvent* event); +void printPT2s(LSTEvent* event); void printTCs(LSTEvent* event); // Print anomalous multiplicities diff --git a/RecoTracker/LSTCore/standalone/efficiency/python/lst_plot_performance.py b/RecoTracker/LSTCore/standalone/efficiency/python/lst_plot_performance.py index f504ac48826ac..6b5442aef9f37 100755 --- a/RecoTracker/LSTCore/standalone/efficiency/python/lst_plot_performance.py +++ b/RecoTracker/LSTCore/standalone/efficiency/python/lst_plot_performance.py @@ -10,7 +10,7 @@ sel_choices = ["base", "loweta", "xtr", "vtr", "none"] metric_choices = ["eff", "fakerate", "duplrate"] variable_choices = ["pt", "ptmtv", "ptlow", "eta", "phi", "dxy", "dz", "vxy"] -objecttype_choices = ["TC", "pT5", "T5", "pT3", "pLS", "pT5_lower", "pT3_lower", "T5_lower"] +objecttype_choices = ["TC", "pT5", "T5", "pT3", "pLS", "pT2", "pT5_lower", "pT3_lower", "T5_lower", "pT2_lower"] #lowerObjectType = ["pT5_lower", "pT3_lower", "T5_lower"] r.gROOT.SetBatch(True) @@ -118,7 +118,7 @@ def plot(args): numer = [] numer.append(params["input_file"].Get(params["numer"]).Clone()) - breakdown_hist_types = ["pT5", "pT3", "T5", "pLS"] + breakdown_hist_types = ["pT5", "pT3", "T5", "pLS", "pT2"] print("breakdown = ", params["breakdown"]) if params["breakdown"]: for breakdown_hist_type in breakdown_hist_types: @@ -136,7 +136,7 @@ def plot(args): if params["breakdown"]: - params["legend_labels"] = ["TC" ,"pT5" ,"pT3" ,"T5" ,"pLS"] + params["legend_labels"] = ["TC" ,"pT5" ,"pT3" ,"T5" ,"pLS", "pT2"] else: params["legend_labels"] = [args.objecttype] @@ -377,6 +377,8 @@ def parse_plot_name(output_name): rtnstr.append("Quadruplet w/ gap") elif "pT3_" in output_name: rtnstr.append("Pixel Triplet") + elif "pT2_" in output_name: + rtnstr.append("Pixel T2") elif "pT5_" in output_name: rtnstr.append("Pixel Quintuplet") elif "T3_" in output_name: @@ -537,8 +539,8 @@ def draw_plot(effs, nums, dens, params): effs[0].SetTitle(parse_plot_name(output_name)) # Draw the efficiency graphs - colors = [1, 2, 3, 4, 6] - markerstyles = [20, 26, 28, 24, 27] + colors = [1, 2, 3, 4, 6, 8] + markerstyles = [20, 26, 28, 24, 27, 22] markersize = 1.2 linewidth = 2 for i, eff in enumerate(effs): @@ -672,9 +674,11 @@ def plot_standard_performance_plots(args): "pT3": [False], "T5": [False], "pLS": [False], + "pT2": [False], "pT5_lower":[False], "pT3_lower":[False], "T5_lower":[False], + "pT2_lower":[False], }, "fakerate":{ "TC": [True, False], @@ -682,9 +686,11 @@ def plot_standard_performance_plots(args): "pT3": [False], "T5": [False], "pLS": [False], + "pT2": [False], "pT5_lower":[False], "pT3_lower":[False], "T5_lower":[False], + "pT2_lower":[False], }, "duplrate":{ "TC": [True, False], @@ -692,9 +698,11 @@ def plot_standard_performance_plots(args): "pT3": [False], "T5": [False], "pLS": [False], + "pT2": [False], "pT5_lower":[False], "pT3_lower":[False], "T5_lower":[False], + "pT2_lower":[False], }, } pdgids = { @@ -743,16 +751,16 @@ def plot_standard_performance_plots(args): if args.individual: # Only eff / TC matters here - breakdowns = {"eff":{"TC":[False], "pT5_lower":[False], "pT3_lower":[False], "T5_lower":[False]}, - "fakerate": {"TC":[False], "pT5_lower":[False], "pT3_lower":[False], "T5_lower":[False]}, - "duplrate": {"TC":[False], "pT5_lower":[False], "pT3_lower":[False], "T5_lower":[False]}} + breakdowns = {"eff":{"TC":[False], "pT5_lower":[False], "pT3_lower":[False], "T5_lower":[False], "pT2_lower":[False]}, + "fakerate": {"TC":[False], "pT5_lower":[False], "pT3_lower":[False], "T5_lower":[False], "pT2_lower":[False]}, + "duplrate": {"TC":[False], "pT5_lower":[False], "pT3_lower":[False], "T5_lower":[False], "pT2_lower":[False]}} else: # Only eff / TC matters here - breakdowns = {"eff":{"TC":[True], "pT5_lower":[False], "pT3_lower":[False], "T5_lower":[False]}, - "fakerate": {"TC":[True], "pT5_lower":[False], "pT3_lower":[False], "T5_lower":[False]}, - "duplrate": {"TC":[True], "pT5_lower":[False], "pT3_lower":[False], "T5_lower":[False]}} + breakdowns = {"eff":{"TC":[True], "pT5_lower":[False], "pT3_lower":[False], "T5_lower":[False], "pT2_lower":[False]}, + "fakerate": {"TC":[True], "pT5_lower":[False], "pT3_lower":[False], "T5_lower":[False], "pT2_lower":[False]}, + "duplrate": {"TC":[True], "pT5_lower":[False], "pT3_lower":[False], "T5_lower":[False], "pT2_lower":[False]}} if args.yzoom: args.yzooms = [args.yzoom] diff --git a/RecoTracker/LSTCore/standalone/efficiency/src/LSTEff.cc b/RecoTracker/LSTCore/standalone/efficiency/src/LSTEff.cc index 6474f9a17ae04..03c696da4086b 100644 --- a/RecoTracker/LSTCore/standalone/efficiency/src/LSTEff.cc +++ b/RecoTracker/LSTCore/standalone/efficiency/src/LSTEff.cc @@ -122,6 +122,13 @@ void LSTEff::Init(TTree *tree) { pT3_isDuplicate_branch->SetAddress(&pT3_isDuplicate_); } } + pT2_isDuplicate_branch = 0; + if (tree->GetBranch("pT2_isDuplicate") != 0) { + pT2_isDuplicate_branch = tree->GetBranch("pT2_isDuplicate"); + if (pT2_isDuplicate_branch) { + pT2_isDuplicate_branch->SetAddress(&pT2_isDuplicate_); + } + } tc_isDuplicate_branch = 0; if (tree->GetBranch("tc_isDuplicate") != 0) { tc_isDuplicate_branch = tree->GetBranch("tc_isDuplicate"); @@ -136,6 +143,13 @@ void LSTEff::Init(TTree *tree) { pT3_eta_2_branch->SetAddress(&pT3_eta_2_); } } + pT2_eta_2_branch = 0; + if (tree->GetBranch("pT2_eta_2") != 0) { + pT2_eta_2_branch = tree->GetBranch("pT2_eta_2"); + if (pT2_eta_2_branch) { + pT2_eta_2_branch->SetAddress(&pT2_eta_2_); + } + } sim_pT3_matched_branch = 0; if (tree->GetBranch("sim_pT3_matched") != 0) { sim_pT3_matched_branch = tree->GetBranch("sim_pT3_matched"); @@ -143,6 +157,13 @@ void LSTEff::Init(TTree *tree) { sim_pT3_matched_branch->SetAddress(&sim_pT3_matched_); } } + sim_pT2_matched_branch = 0; + if (tree->GetBranch("sim_pT2_matched") != 0) { + sim_pT2_matched_branch = tree->GetBranch("sim_pT2_matched"); + if (sim_pT2_matched_branch) { + sim_pT2_matched_branch->SetAddress(&sim_pT2_matched_); + } + } pureTCE_rzChiSquared_branch = 0; if (tree->GetBranch("pureTCE_rzChiSquared") != 0) { pureTCE_rzChiSquared_branch = tree->GetBranch("pureTCE_rzChiSquared"); @@ -423,6 +444,13 @@ void LSTEff::Init(TTree *tree) { pT3_isFake_branch->SetAddress(&pT3_isFake_); } } + pT2_isFake_branch = 0; + if (tree->GetBranch("pT2_isFake") != 0) { + pT2_isFake_branch = tree->GetBranch("pT2_isFake"); + if (pT2_isFake_branch) { + pT2_isFake_branch->SetAddress(&pT2_isFake_); + } + } tce_nLayerOverlaps_branch = 0; if (tree->GetBranch("tce_nLayerOverlaps") != 0) { tce_nLayerOverlaps_branch = tree->GetBranch("tce_nLayerOverlaps"); @@ -549,6 +577,13 @@ void LSTEff::Init(TTree *tree) { pT3_pt_branch->SetAddress(&pT3_pt_); } } + pT2_pt_branch = 0; + if (tree->GetBranch("pT2_pt") != 0) { + pT2_pt_branch = tree->GetBranch("pT2_pt"); + if (pT2_pt_branch) { + pT2_pt_branch->SetAddress(&pT2_pt_); + } + } tc_pt_branch = 0; if (tree->GetBranch("tc_pt") != 0) { tc_pt_branch = tree->GetBranch("tc_pt"); @@ -563,6 +598,13 @@ void LSTEff::Init(TTree *tree) { pT3_phi_2_branch->SetAddress(&pT3_phi_2_); } } + pT2_phi_2_branch = 0; + if (tree->GetBranch("pT2_phi_2") != 0) { + pT2_phi_2_branch = tree->GetBranch("pT2_phi_2"); + if (pT2_phi_2_branch) { + pT2_phi_2_branch->SetAddress(&pT2_phi_2_); + } + } pT5_pt_branch = 0; if (tree->GetBranch("pT5_pt") != 0) { pT5_pt_branch = tree->GetBranch("pT5_pt"); @@ -668,6 +710,13 @@ void LSTEff::Init(TTree *tree) { pT3_phi_branch->SetAddress(&pT3_phi_); } } + pT2_phi_branch = 0; + if (tree->GetBranch("pT2_phi") != 0) { + pT2_phi_branch = tree->GetBranch("pT2_phi"); + if (pT2_phi_branch) { + pT2_phi_branch->SetAddress(&pT2_phi_); + } + } pT5_eta_branch = 0; if (tree->GetBranch("pT5_eta") != 0) { pT5_eta_branch = tree->GetBranch("pT5_eta"); @@ -724,6 +773,13 @@ void LSTEff::Init(TTree *tree) { pT3_eta_branch->SetAddress(&pT3_eta_); } } + pT2_eta_branch = 0; + if (tree->GetBranch("pT2_eta") != 0) { + pT2_eta_branch = tree->GetBranch("pT2_eta"); + if (pT2_eta_branch) { + pT2_eta_branch->SetAddress(&pT2_eta_); + } + } sim_parentVtxIdx_branch = 0; if (tree->GetBranch("sim_parentVtxIdx") != 0) { sim_parentVtxIdx_branch = tree->GetBranch("sim_parentVtxIdx"); @@ -1046,9 +1102,12 @@ void LSTEff::GetEntry(unsigned int idx) { pT5_isDuplicate_isLoaded = false; sim_tce_matched_isLoaded = false; pT3_isDuplicate_isLoaded = false; + pT2_isDuplicate_isLoaded = false; tc_isDuplicate_isLoaded = false; pT3_eta_2_isLoaded = false; + pT2_eta_2_isLoaded = false; sim_pT3_matched_isLoaded = false; + sim_pT2_matched_isLoaded = false; pureTCE_rzChiSquared_isLoaded = false; t4_isDuplicate_isLoaded = false; pureTCE_eta_isLoaded = false; @@ -1089,6 +1148,7 @@ void LSTEff::GetEntry(unsigned int idx) { tce_pt_isLoaded = false; tc_isFake_isLoaded = false; pT3_isFake_isLoaded = false; + pT2_isFake_isLoaded = false; tce_nLayerOverlaps_isLoaded = false; tc_sim_isLoaded = false; sim_pLS_types_isLoaded = false; @@ -1107,8 +1167,10 @@ void LSTEff::GetEntry(unsigned int idx) { sim_T4_matched_isLoaded = false; sim_isGood_isLoaded = false; pT3_pt_isLoaded = false; + pT2_pt_isLoaded = false; tc_pt_isLoaded = false; pT3_phi_2_isLoaded = false; + pT2_phi_2_isLoaded = false; pT5_pt_isLoaded = false; pureTCE_rPhiChiSquared_isLoaded = false; pT5_score_isLoaded = false; @@ -1124,6 +1186,7 @@ void LSTEff::GetEntry(unsigned int idx) { sim_T3_matched_isLoaded = false; pLS_score_isLoaded = false; pT3_phi_isLoaded = false; + pT2_phi_isLoaded = false; pT5_eta_isLoaded = false; tc_phi_isLoaded = false; t4_eta_isLoaded = false; @@ -1132,6 +1195,7 @@ void LSTEff::GetEntry(unsigned int idx) { sim_bunchCrossing_isLoaded = false; tc_partOfExtension_isLoaded = false; pT3_eta_isLoaded = false; + pT2_eta_isLoaded = false; sim_parentVtxIdx_isLoaded = false; pureTCE_layer_binary_isLoaded = false; sim_pT4_matched_isLoaded = false; @@ -1211,12 +1275,18 @@ void LSTEff::LoadAllBranches() { sim_tce_matched(); if (pT3_isDuplicate_branch != 0) pT3_isDuplicate(); + if (pT2_isDuplicate_branch != 0) + pT2_isDuplicate(); if (tc_isDuplicate_branch != 0) tc_isDuplicate(); if (pT3_eta_2_branch != 0) pT3_eta_2(); + if (pT2_eta_2_branch != 0) + pT2_eta_2(); if (sim_pT3_matched_branch != 0) sim_pT3_matched(); + if (sim_pT2_matched_branch != 0) + sim_pT2_matched(); if (pureTCE_rzChiSquared_branch != 0) pureTCE_rzChiSquared(); if (t4_isDuplicate_branch != 0) @@ -1297,6 +1367,8 @@ void LSTEff::LoadAllBranches() { tc_isFake(); if (pT3_isFake_branch != 0) pT3_isFake(); + if (pT2_isFake_branch != 0) + pT2_isFake(); if (tce_nLayerOverlaps_branch != 0) tce_nLayerOverlaps(); if (tc_sim_branch != 0) @@ -1333,10 +1405,14 @@ void LSTEff::LoadAllBranches() { sim_isGood(); if (pT3_pt_branch != 0) pT3_pt(); + if (pT2_pt_branch != 0) + pT2_pt(); if (tc_pt_branch != 0) tc_pt(); if (pT3_phi_2_branch != 0) pT3_phi_2(); + if (pT2_phi_2_branch != 0) + pT2_phi_2(); if (pT5_pt_branch != 0) pT5_pt(); if (pureTCE_rPhiChiSquared_branch != 0) @@ -1367,6 +1443,8 @@ void LSTEff::LoadAllBranches() { pLS_score(); if (pT3_phi_branch != 0) pT3_phi(); + if (pT2_phi_branch != 0) + pT2_phi(); if (pT5_eta_branch != 0) pT5_eta(); if (tc_phi_branch != 0) @@ -1383,6 +1461,8 @@ void LSTEff::LoadAllBranches() { tc_partOfExtension(); if (pT3_eta_branch != 0) pT3_eta(); + if (pT2_eta_branch != 0) + pT2_eta(); if (sim_parentVtxIdx_branch != 0) sim_parentVtxIdx(); if (pureTCE_layer_binary_branch != 0) @@ -1674,6 +1754,18 @@ const std::vector &LSTEff::pT3_isDuplicate() { } return *pT3_isDuplicate_; } +const std::vector &LSTEff::pT2_isDuplicate() { + if (not pT2_isDuplicate_isLoaded) { + if (pT2_isDuplicate_branch != 0) { + pT2_isDuplicate_branch->GetEntry(index); + } else { + printf("branch pT2_isDuplicate_branch does not exist!\n"); + exit(1); + } + pT2_isDuplicate_isLoaded = true; + } + return *pT2_isDuplicate_; +} const std::vector &LSTEff::tc_isDuplicate() { if (not tc_isDuplicate_isLoaded) { if (tc_isDuplicate_branch != 0) { @@ -1698,6 +1790,18 @@ const std::vector &LSTEff::pT3_eta_2() { } return *pT3_eta_2_; } +const std::vector &LSTEff::pT2_eta_2() { + if (not pT2_eta_2_isLoaded) { + if (pT2_eta_2_branch != 0) { + pT2_eta_2_branch->GetEntry(index); + } else { + printf("branch pT2_eta_2_branch does not exist!\n"); + exit(1); + } + pT2_eta_2_isLoaded = true; + } + return *pT2_eta_2_; +} const std::vector &LSTEff::sim_pT3_matched() { if (not sim_pT3_matched_isLoaded) { if (sim_pT3_matched_branch != 0) { @@ -1710,6 +1814,20 @@ const std::vector &LSTEff::sim_pT3_matched() { } return *sim_pT3_matched_; } + +const std::vector &LSTEff::sim_pT2_matched() { + if (not sim_pT2_matched_isLoaded) { + if (sim_pT2_matched_branch != 0) { + sim_pT2_matched_branch->GetEntry(index); + } else { + printf("branch sim_pT2_matched_branch does not exist!\n"); + exit(1); + } + sim_pT2_matched_isLoaded = true; + } + return *sim_pT2_matched_; +} + const std::vector &LSTEff::pureTCE_rzChiSquared() { if (not pureTCE_rzChiSquared_isLoaded) { if (pureTCE_rzChiSquared_branch != 0) { @@ -2190,6 +2308,18 @@ const std::vector &LSTEff::pT3_isFake() { } return *pT3_isFake_; } +const std::vector &LSTEff::pT2_isFake() { + if (not pT2_isFake_isLoaded) { + if (pT2_isFake_branch != 0) { + pT2_isFake_branch->GetEntry(index); + } else { + printf("branch pT2_isFake_branch does not exist!\n"); + exit(1); + } + pT2_isFake_isLoaded = true; + } + return *pT2_isFake_; +} const std::vector > &LSTEff::tce_nLayerOverlaps() { if (not tce_nLayerOverlaps_isLoaded) { if (tce_nLayerOverlaps_branch != 0) { @@ -2406,6 +2536,18 @@ const std::vector &LSTEff::pT3_pt() { } return *pT3_pt_; } +const std::vector &LSTEff::pT2_pt() { + if (not pT2_pt_isLoaded) { + if (pT2_pt_branch != 0) { + pT2_pt_branch->GetEntry(index); + } else { + printf("branch pT2_pt_branch does not exist!\n"); + exit(1); + } + pT2_pt_isLoaded = true; + } + return *pT2_pt_; +} const std::vector &LSTEff::tc_pt() { if (not tc_pt_isLoaded) { if (tc_pt_branch != 0) { @@ -2430,6 +2572,18 @@ const std::vector &LSTEff::pT3_phi_2() { } return *pT3_phi_2_; } +const std::vector &LSTEff::pT2_phi_2() { + if (not pT2_phi_2_isLoaded) { + if (pT2_phi_2_branch != 0) { + pT2_phi_2_branch->GetEntry(index); + } else { + printf("branch pT2_phi_2_branch does not exist!\n"); + exit(1); + } + pT2_phi_2_isLoaded = true; + } + return *pT2_phi_2_; +} const std::vector &LSTEff::pT5_pt() { if (not pT5_pt_isLoaded) { if (pT5_pt_branch != 0) { @@ -2610,6 +2764,18 @@ const std::vector &LSTEff::pT3_phi() { } return *pT3_phi_; } +const std::vector &LSTEff::pT2_phi() { + if (not pT2_phi_isLoaded) { + if (pT2_phi_branch != 0) { + pT2_phi_branch->GetEntry(index); + } else { + printf("branch pT2_phi_branch does not exist!\n"); + exit(1); + } + pT2_phi_isLoaded = true; + } + return *pT2_phi_; +} const std::vector &LSTEff::pT5_eta() { if (not pT5_eta_isLoaded) { if (pT5_eta_branch != 0) { @@ -2706,6 +2872,18 @@ const std::vector &LSTEff::pT3_eta() { } return *pT3_eta_; } +const std::vector &LSTEff::pT2_eta() { + if (not pT2_eta_isLoaded) { + if (pT2_eta_branch != 0) { + pT2_eta_branch->GetEntry(index); + } else { + printf("branch pT2_eta_branch does not exist!\n"); + exit(1); + } + pT2_eta_isLoaded = true; + } + return *pT2_eta_; +} const std::vector &LSTEff::sim_parentVtxIdx() { if (not sim_parentVtxIdx_isLoaded) { if (sim_parentVtxIdx_branch != 0) { @@ -3261,9 +3439,12 @@ namespace tas { const std::vector &pT5_isDuplicate() { return lstEff.pT5_isDuplicate(); } const std::vector &sim_tce_matched() { return lstEff.sim_tce_matched(); } const std::vector &pT3_isDuplicate() { return lstEff.pT3_isDuplicate(); } + const std::vector &pT2_isDuplicate() { return lstEff.pT2_isDuplicate(); } const std::vector &tc_isDuplicate() { return lstEff.tc_isDuplicate(); } const std::vector &pT3_eta_2() { return lstEff.pT3_eta_2(); } + const std::vector &pT2_eta_2() { return lstEff.pT2_eta_2(); } const std::vector &sim_pT3_matched() { return lstEff.sim_pT3_matched(); } + const std::vector &sim_pT2_matched() { return lstEff.sim_pT2_matched(); } const std::vector &pureTCE_rzChiSquared() { return lstEff.pureTCE_rzChiSquared(); } const std::vector &t4_isDuplicate() { return lstEff.t4_isDuplicate(); } const std::vector &pureTCE_eta() { return lstEff.pureTCE_eta(); } @@ -3304,6 +3485,7 @@ namespace tas { const std::vector &tce_pt() { return lstEff.tce_pt(); } const std::vector &tc_isFake() { return lstEff.tc_isFake(); } const std::vector &pT3_isFake() { return lstEff.pT3_isFake(); } + const std::vector &pT2_isFake() { return lstEff.pT2_isFake(); } const std::vector > &tce_nLayerOverlaps() { return lstEff.tce_nLayerOverlaps(); } const std::vector &tc_sim() { return lstEff.tc_sim(); } const std::vector > &sim_pLS_types() { return lstEff.sim_pLS_types(); } @@ -3322,8 +3504,10 @@ namespace tas { const std::vector &sim_T4_matched() { return lstEff.sim_T4_matched(); } const std::vector &sim_isGood() { return lstEff.sim_isGood(); } const std::vector &pT3_pt() { return lstEff.pT3_pt(); } + const std::vector &pT2_pt() { return lstEff.pT2_pt(); } const std::vector &tc_pt() { return lstEff.tc_pt(); } const std::vector &pT3_phi_2() { return lstEff.pT3_phi_2(); } + const std::vector &pT2_phi_2() { return lstEff.pT2_phi_2(); } const std::vector &pT5_pt() { return lstEff.pT5_pt(); } const std::vector &pureTCE_rPhiChiSquared() { return lstEff.pureTCE_rPhiChiSquared(); } const std::vector &pT5_score() { return lstEff.pT5_score(); } @@ -3339,6 +3523,7 @@ namespace tas { const std::vector &sim_T3_matched() { return lstEff.sim_T3_matched(); } const std::vector &pLS_score() { return lstEff.pLS_score(); } const std::vector &pT3_phi() { return lstEff.pT3_phi(); } + const std::vector &pT2_phi() { return lstEff.pT2_phi(); } const std::vector &pT5_eta() { return lstEff.pT5_eta(); } const std::vector &tc_phi() { return lstEff.tc_phi(); } const std::vector &t4_eta() { return lstEff.t4_eta(); } @@ -3347,6 +3532,7 @@ namespace tas { const std::vector &sim_bunchCrossing() { return lstEff.sim_bunchCrossing(); } const std::vector &tc_partOfExtension() { return lstEff.tc_partOfExtension(); } const std::vector &pT3_eta() { return lstEff.pT3_eta(); } + const std::vector &pT2_eta() { return lstEff.pT2_eta(); } const std::vector &sim_parentVtxIdx() { return lstEff.sim_parentVtxIdx(); } const std::vector &pureTCE_layer_binary() { return lstEff.pureTCE_layer_binary(); } const std::vector &sim_pT4_matched() { return lstEff.sim_pT4_matched(); } diff --git a/RecoTracker/LSTCore/standalone/efficiency/src/LSTEff.h b/RecoTracker/LSTCore/standalone/efficiency/src/LSTEff.h index 9c53fb92c7f62..ad704858bc8fa 100644 --- a/RecoTracker/LSTCore/standalone/efficiency/src/LSTEff.h +++ b/RecoTracker/LSTCore/standalone/efficiency/src/LSTEff.h @@ -68,15 +68,24 @@ class LSTEff { std::vector *pT3_isDuplicate_; TBranch *pT3_isDuplicate_branch; bool pT3_isDuplicate_isLoaded; + std::vector *pT2_isDuplicate_; + TBranch *pT2_isDuplicate_branch; + bool pT2_isDuplicate_isLoaded; std::vector *tc_isDuplicate_; TBranch *tc_isDuplicate_branch; bool tc_isDuplicate_isLoaded; std::vector *pT3_eta_2_; TBranch *pT3_eta_2_branch; bool pT3_eta_2_isLoaded; + std::vector *pT2_eta_2_; + TBranch *pT2_eta_2_branch; + bool pT2_eta_2_isLoaded; std::vector *sim_pT3_matched_; TBranch *sim_pT3_matched_branch; bool sim_pT3_matched_isLoaded; + std::vector *sim_pT2_matched_; + TBranch *sim_pT2_matched_branch; + bool sim_pT2_matched_isLoaded; std::vector *pureTCE_rzChiSquared_; TBranch *pureTCE_rzChiSquared_branch; bool pureTCE_rzChiSquared_isLoaded; @@ -197,6 +206,9 @@ class LSTEff { std::vector *pT3_isFake_; TBranch *pT3_isFake_branch; bool pT3_isFake_isLoaded; + std::vector *pT2_isFake_; + TBranch *pT2_isFake_branch; + bool pT2_isFake_isLoaded; std::vector > *tce_nLayerOverlaps_; TBranch *tce_nLayerOverlaps_branch; bool tce_nLayerOverlaps_isLoaded; @@ -251,12 +263,18 @@ class LSTEff { std::vector *pT3_pt_; TBranch *pT3_pt_branch; bool pT3_pt_isLoaded; + std::vector *pT2_pt_; + TBranch *pT2_pt_branch; + bool pT2_pt_isLoaded; std::vector *tc_pt_; TBranch *tc_pt_branch; bool tc_pt_isLoaded; std::vector *pT3_phi_2_; TBranch *pT3_phi_2_branch; bool pT3_phi_2_isLoaded; + std::vector *pT2_phi_2_; + TBranch *pT2_phi_2_branch; + bool pT2_phi_2_isLoaded; std::vector *pT5_pt_; TBranch *pT5_pt_branch; bool pT5_pt_isLoaded; @@ -302,6 +320,9 @@ class LSTEff { std::vector *pT3_phi_; TBranch *pT3_phi_branch; bool pT3_phi_isLoaded; + std::vector *pT2_phi_; + TBranch *pT2_phi_branch; + bool pT2_phi_isLoaded; std::vector *pT5_eta_; TBranch *pT5_eta_branch; bool pT5_eta_isLoaded; @@ -326,6 +347,9 @@ class LSTEff { std::vector *pT3_eta_; TBranch *pT3_eta_branch; bool pT3_eta_isLoaded; + std::vector *pT2_eta_; + TBranch *pT2_eta_branch; + bool pT2_eta_isLoaded; std::vector *sim_parentVtxIdx_; TBranch *sim_parentVtxIdx_branch; bool sim_parentVtxIdx_isLoaded; @@ -477,9 +501,12 @@ class LSTEff { const std::vector &pT5_isDuplicate(); const std::vector &sim_tce_matched(); const std::vector &pT3_isDuplicate(); + const std::vector &pT2_isDuplicate(); const std::vector &tc_isDuplicate(); const std::vector &pT3_eta_2(); + const std::vector &pT2_eta_2(); const std::vector &sim_pT3_matched(); + const std::vector &sim_pT2_matched(); const std::vector &pureTCE_rzChiSquared(); const std::vector &t4_isDuplicate(); const std::vector &pureTCE_eta(); @@ -520,6 +547,7 @@ class LSTEff { const std::vector &tce_pt(); const std::vector &tc_isFake(); const std::vector &pT3_isFake(); + const std::vector &pT2_isFake(); const std::vector > &tce_nLayerOverlaps(); const std::vector &tc_sim(); const std::vector > &sim_pLS_types(); @@ -538,8 +566,10 @@ class LSTEff { const std::vector &sim_T4_matched(); const std::vector &sim_isGood(); const std::vector &pT3_pt(); + const std::vector &pT2_pt(); const std::vector &tc_pt(); const std::vector &pT3_phi_2(); + const std::vector &pT2_phi_2(); const std::vector &pT5_pt(); const std::vector &pureTCE_rPhiChiSquared(); const std::vector &pT5_score(); @@ -555,6 +585,7 @@ class LSTEff { const std::vector &sim_T3_matched(); const std::vector &pLS_score(); const std::vector &pT3_phi(); + const std::vector &pT2_phi(); const std::vector &pT5_eta(); const std::vector &tc_phi(); const std::vector &t4_eta(); @@ -563,6 +594,7 @@ class LSTEff { const std::vector &sim_bunchCrossing(); const std::vector &tc_partOfExtension(); const std::vector &pT3_eta(); + const std::vector &pT2_eta(); const std::vector &sim_parentVtxIdx(); const std::vector &pureTCE_layer_binary(); const std::vector &sim_pT4_matched(); @@ -632,9 +664,12 @@ namespace tas { const std::vector &pT5_isDuplicate(); const std::vector &sim_tce_matched(); const std::vector &pT3_isDuplicate(); + const std::vector &pT2_isDuplicate(); const std::vector &tc_isDuplicate(); const std::vector &pT3_eta_2(); + const std::vector &pT2_eta_2(); const std::vector &sim_pT3_matched(); + const std::vector &sim_pT2_matched(); const std::vector &pureTCE_rzChiSquared(); const std::vector &t4_isDuplicate(); const std::vector &pureTCE_eta(); @@ -675,6 +710,7 @@ namespace tas { const std::vector &tce_pt(); const std::vector &tc_isFake(); const std::vector &pT3_isFake(); + const std::vector &pT2_isFake(); const std::vector > &tce_nLayerOverlaps(); const std::vector &tc_sim(); const std::vector > &sim_pLS_types(); @@ -693,8 +729,10 @@ namespace tas { const std::vector &sim_T4_matched(); const std::vector &sim_isGood(); const std::vector &pT3_pt(); + const std::vector &pT2_pt(); const std::vector &tc_pt(); const std::vector &pT3_phi_2(); + const std::vector &pT2_phi_2(); const std::vector &pT5_pt(); const std::vector &pureTCE_rPhiChiSquared(); const std::vector &pT5_score(); @@ -710,6 +748,7 @@ namespace tas { const std::vector &sim_T3_matched(); const std::vector &pLS_score(); const std::vector &pT3_phi(); + const std::vector &pT2_phi(); const std::vector &pT5_eta(); const std::vector &tc_phi(); const std::vector &t4_eta(); @@ -718,6 +757,7 @@ namespace tas { const std::vector &sim_bunchCrossing(); const std::vector &tc_partOfExtension(); const std::vector &pT3_eta(); + const std::vector &pT2_eta(); const std::vector &sim_parentVtxIdx(); const std::vector &pureTCE_layer_binary(); const std::vector &sim_pT4_matched(); diff --git a/RecoTracker/LSTCore/standalone/efficiency/src/performance.cc b/RecoTracker/LSTCore/standalone/efficiency/src/performance.cc index fa7bbc0edf8b0..638687946478e 100644 --- a/RecoTracker/LSTCore/standalone/efficiency/src/performance.cc +++ b/RecoTracker/LSTCore/standalone/efficiency/src/performance.cc @@ -1,6 +1,6 @@ #include "performance.h" -enum { pT5 = 7, pT3 = 5, T5 = 4, pLS = 8 }; +enum { pT5 = 7, pT3 = 5, T5 = 4, pLS = 8, pT2 = 10}; //__________________________________________________________________________________________________________________________________________________________________________ int main(int argc, char** argv) { @@ -74,7 +74,13 @@ int main(int argc, char** argv) { /* q */ charge, /* pass */ [&](unsigned int isim) { return lstEff.sim_TC_matched_mask().at(isim) & (1 << pLS); }, /* sel */ sels[isel])); - + list_effSetDef.push_back(SimTrackSetDefinition( + /* name */ + TString("pT2_") + selnames[isel], + /* pdgid */ pdgid, + /* q */ charge, + /* pass */ [&](unsigned int isim) { return lstEff.sim_TC_matched_mask().at(isim) & (1 << pT2); }, + /* sel */ sels[isel])); if (ana.do_lower_level) { //lower objects - name will have pT5_lower_, T5_lower_, pT3_lower_ list_effSetDef.push_back(SimTrackSetDefinition( @@ -98,6 +104,14 @@ int main(int argc, char** argv) { /* q */ charge, /* pass */ [&](unsigned int isim) { return lstEff.sim_pT3_matched().at(isim) > 0; }, /* sel */ sels[isel])); + + list_effSetDef.push_back(SimTrackSetDefinition( + /* name */ + TString("pT2_lower_") + selnames[isel], + /* pdgid */ pdgid, + /* q */ charge, + /* pass */ [&](unsigned int isim) { return lstEff.sim_pT2_matched().at(isim) > 0; }, + /* sel */ sels[isel])); } } } @@ -152,6 +166,15 @@ int main(int argc, char** argv) { /* eta */ tas::tc_eta, /* phi */ tas::tc_phi, /* type */ tas::tc_type)); + list_FRSetDef.push_back( + RecoTrackSetDefinition(/* name */ + "pT2", + /* pass */ [&](unsigned int itc) { return lstEff.tc_isFake().at(itc) > 0; }, + /* sel */ [&](unsigned int itc) { return lstEff.tc_type().at(itc) == pT2; }, + /* pt */ tas::tc_pt, + /* eta */ tas::tc_eta, + /* phi */ tas::tc_phi, + /* type */ tas::tc_type)); if (ana.do_lower_level) { list_FRSetDef.push_back(RecoTrackSetDefinition( @@ -181,6 +204,16 @@ int main(int argc, char** argv) { /* eta */ tas::pT3_eta, /* phi */ tas::pT3_phi, /* type */ [&]() { return static_cast>(std::vector(tas::pT3_pt().size(), 1)); })); + list_FRSetDef.push_back(RecoTrackSetDefinition( + /* name */ + "pT2_lower", + /* pass */ [&](unsigned int ipT2) { return lstEff.pT2_isFake().at(ipT2) > 0; }, + /* sel */ [&](unsigned int ipT2) { return 1; }, + /* pt */ tas::pT2_pt, + /* eta */ tas::pT2_eta, + /* phi */ tas::pT2_phi, + /* type */ [&]() { return static_cast>(std::vector(tas::pT2_pt().size(), 1)); })); + } bookFakeRateSets(list_FRSetDef); @@ -232,6 +265,15 @@ int main(int argc, char** argv) { /* eta */ tas::tc_eta, /* phi */ tas::tc_phi, /* type */ tas::tc_type)); + list_DRSetDef.push_back( + RecoTrackSetDefinition(/* name */ + "pT2", + /* pass */ [&](unsigned int itc) { return lstEff.tc_isDuplicate().at(itc) > 0; }, + /* sel */ [&](unsigned int itc) { return lstEff.tc_type().at(itc) == pT2; }, + /* pt */ tas::tc_pt, + /* eta */ tas::tc_eta, + /* phi */ tas::tc_phi, + /* type */ tas::tc_type)); if (ana.do_lower_level) { list_DRSetDef.push_back(RecoTrackSetDefinition( @@ -261,6 +303,16 @@ int main(int argc, char** argv) { /* eta */ tas::pT3_eta, /* phi */ tas::pT3_phi, /* type */ [&]() { return static_cast>(std::vector(tas::pT3_pt().size(), 1)); })); + list_DRSetDef.push_back(RecoTrackSetDefinition( + /* name */ + "pT2_lower", + /* pass */ [&](unsigned int ipT2) { return lstEff.pT2_isDuplicate().at(ipT2) > 0; }, + /* sel */ [&](unsigned int ipT2) { return 1; }, + /* pt */ tas::pT2_pt, + /* eta */ tas::pT2_eta, + /* phi */ tas::pT2_phi, + /* type */ [&]() { return static_cast>(std::vector(tas::pT2_pt().size(), 1)); })); + } bookDuplicateRateSets(list_DRSetDef);