Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

LSTOD LS/MD Changes #412

Draft
wants to merge 5 commits into
base: lstod
Choose a base branch
from
Draft
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions SDL/MiniDoublet.h
Original file line number Diff line number Diff line change
@@ -1,5 +1,6 @@
#ifndef MiniDoublet_cuh
#define MiniDoublet_cuh
#define OUTPUT_MD_CUTS

#include "Constants.h"
#include "EndcapGeometry.h"
Expand Down
139 changes: 127 additions & 12 deletions SDL/Segment.h
Original file line number Diff line number Diff line change
@@ -1,5 +1,6 @@
#ifndef Segment_cuh
#define Segment_cuh
#define OUTPUT_LS_CUTS

#include "Constants.h"
#include "EndcapGeometry.h"
Expand All @@ -15,6 +16,15 @@ namespace SDL {
FPX* dPhiChanges;
FPX* dPhiChangeMins;
FPX* dPhiChangeMaxs;
#ifdef OUTPUT_LS_CUTS
FPX* zHis;
FPX* zLos;
FPX* rtHis;
FPX* rtLos;
FPX* dAlphaInners;
FPX* dAlphaOuters;
FPX* dAlphaInnerOuters;
#endif
uint16_t* innerLowerModuleIndices;
uint16_t* outerLowerModuleIndices;
unsigned int* seedIdx;
Expand Down Expand Up @@ -52,6 +62,17 @@ namespace SDL {
dPhiChanges = alpaka::getPtrNative(segmentsbuf.dPhiChanges_buf);
dPhiChangeMins = alpaka::getPtrNative(segmentsbuf.dPhiChangeMins_buf);
dPhiChangeMaxs = alpaka::getPtrNative(segmentsbuf.dPhiChangeMaxs_buf);

#ifdef OUTPUT_LS_CUTS
zHis = alpaka::getPtrNative(segmentsbuf.zHis_buf);
zLos = alpaka::getPtrNative(segmentsbuf.zLos_buf);
rtHis = alpaka::getPtrNative(segmentsbuf.rtHis_buf);
rtLos = alpaka::getPtrNative(segmentsbuf.rtLos_buf);
dAlphaInners = alpaka::getPtrNative(segmentsbuf.dAlphaInners_buf);
dAlphaOuters = alpaka::getPtrNative(segmentsbuf.dAlphaOuters_buf);
dAlphaInnerOuters = alpaka::getPtrNative(segmentsbuf.dAlphaInnerOuters_buf);

#endif
innerLowerModuleIndices = alpaka::getPtrNative(segmentsbuf.innerLowerModuleIndices_buf);
outerLowerModuleIndices = alpaka::getPtrNative(segmentsbuf.outerLowerModuleIndices_buf);
seedIdx = alpaka::getPtrNative(segmentsbuf.seedIdx_buf);
Expand Down Expand Up @@ -91,6 +112,17 @@ namespace SDL {
Buf<TDev, FPX> dPhiChanges_buf;
Buf<TDev, FPX> dPhiChangeMins_buf;
Buf<TDev, FPX> dPhiChangeMaxs_buf;

#ifdef OUTPUT_LS_CUTS
Buf<TDev, FPX> zHis_buf;
Buf<TDev, FPX> zLos_buf;
Buf<TDev, FPX> rtHis_buf;
Buf<TDev, FPX> rtLos_buf;
Buf<TDev, FPX> dAlphaInners_buf;
Buf<TDev, FPX> dAlphaOuters_buf;
Buf<TDev, FPX> dAlphaInnerOuters_buf;
#endif

Buf<TDev, uint16_t> innerLowerModuleIndices_buf;
Buf<TDev, uint16_t> outerLowerModuleIndices_buf;
Buf<TDev, unsigned int> seedIdx_buf;
Expand Down Expand Up @@ -132,6 +164,16 @@ namespace SDL {
dPhiChanges_buf(allocBufWrapper<FPX>(devAccIn, nMemoryLocationsIn, queue)),
dPhiChangeMins_buf(allocBufWrapper<FPX>(devAccIn, nMemoryLocationsIn, queue)),
dPhiChangeMaxs_buf(allocBufWrapper<FPX>(devAccIn, nMemoryLocationsIn, queue)),
#ifdef OUTPUT_LS_CUTS
zHis_buf(allocBufWrapper<FPX>(devAccIn, nMemoryLocationsIn, queue)),
zLos_buf(allocBufWrapper<FPX>(devAccIn, nMemoryLocationsIn, queue)),
rtHis_buf(allocBufWrapper<FPX>(devAccIn, nMemoryLocationsIn, queue)),
rtLos_buf(allocBufWrapper<FPX>(devAccIn, nMemoryLocationsIn, queue)),
dAlphaInners_buf(allocBufWrapper<FPX>(devAccIn, nMemoryLocationsIn,queue)),
dAlphaOuters_buf(allocBufWrapper<FPX>(devAccIn, nMemoryLocationsIn, queue)),
dAlphaInnerOuters_buf(allocBufWrapper<FPX>(devAccIn, nMemoryLocationsIn, queue)),
#endif

innerLowerModuleIndices_buf(allocBufWrapper<uint16_t>(devAccIn, nMemoryLocationsIn, queue)),
outerLowerModuleIndices_buf(allocBufWrapper<uint16_t>(devAccIn, nMemoryLocationsIn, queue)),
seedIdx_buf(allocBufWrapper<unsigned int>(devAccIn, maxPixelSegments, queue)),
Expand Down Expand Up @@ -331,11 +373,8 @@ namespace SDL {
modulesInGPU.sides[innerLowerModuleIndex] == SDL::Center) {
dAlphaThresholdValues[0] = dAlpha_Bfield + alpaka::math::sqrt(acc, dAlpha_res * dAlpha_res + sdMuls * sdMuls);
} else {
dAlphaThresholdValues[0] =
dAlpha_Bfield +
alpaka::math::sqrt(acc, dAlpha_res * dAlpha_res + sdMuls * sdMuls + sdLumForInnerMini * sdLumForInnerMini);
dAlphaThresholdValues[0] = dAlpha_Bfield + alpaka::math::sqrt(acc, dAlpha_res * dAlpha_res + sdMuls *sdMuls + sdLumForInnerMini * sdLumForOuterMini);
}

if (modulesInGPU.subdets[outerLowerModuleIndex] == SDL::Barrel and
modulesInGPU.sides[outerLowerModuleIndex] == SDL::Center) {
dAlphaThresholdValues[1] = dAlpha_Bfield + alpaka::math::sqrt(acc, dAlpha_res * dAlpha_res + sdMuls * sdMuls);
Expand All @@ -348,7 +387,7 @@ namespace SDL {
//Inner to outer
dAlphaThresholdValues[2] = dAlpha_Bfield + alpaka::math::sqrt(acc, dAlpha_res * dAlpha_res + sdMuls * sdMuls);
};

#ifdef OUTPUT_LS_CUTS
ALPAKA_FN_ACC ALPAKA_FN_INLINE void addSegmentToMemory(struct SDL::segments& segmentsInGPU,
unsigned int lowerMDIndex,
unsigned int upperMDIndex,
Expand All @@ -362,6 +401,13 @@ namespace SDL {
float& dPhiChange,
float& dPhiChangeMin,
float& dPhiChangeMax,
float& zHi,
float& zLo,
float& rtHi,
float& rtLo,
float& dAlphaInner,
float& dAlphaOuter,
float& dAlphaInnerOuter,
unsigned int idx) {
//idx will be computed in the kernel, which is the index into which the
//segment will be written
Expand All @@ -380,7 +426,52 @@ namespace SDL {
segmentsInGPU.dPhiChanges[idx] = __F2H(dPhiChange);
segmentsInGPU.dPhiChangeMins[idx] = __F2H(dPhiChangeMin);
segmentsInGPU.dPhiChangeMaxs[idx] = __F2H(dPhiChangeMax);
}

//hubert
segmentsInGPU.zHis[idx] = __F2H(zHi);
segmentsInGPU.zLos[idx] = __F2H(zLo);
segmentsInGPU.rtHis[idx] = __F2H(rtHi);
segmentsInGPU.rtLos[idx] = __F2H(rtLo);
segmentsInGPU.dAlphaInners[idx] = __F2H(dAlphaInner);
segmentsInGPU.dAlphaOuters[idx] = __F2H(dAlphaOuter);
segmentsInGPU.dAlphaInnerOuters[idx] = __F2H(dAlphaInnerOuter);

};
#else

ALPAKA_FN_ACC ALPAKA_FN_INLINE void addSegmentToMemory(struct SDL::segments& segmentsInGPU,
unsigned int lowerMDIndex,
unsigned int upperMDIndex,
uint16_t innerLowerModuleIndex,
uint16_t outerLowerModuleIndex,
unsigned int innerMDAnchorHitIndex,
unsigned int outerMDAnchorHitIndex,
float& dPhi,
float& dPhiMin,
float& dPhiMax,
float& dPhiChange,
float& dPhiChangeMin,
float& dPhiChangeMax,
unsigned int idx) {
//idx will be computed in the kernel, which is the index into which the
//segment will be written
//nSegments will be incremented in the kernel
//printf("seg: %u %u %u %u\n",lowerMDIndex, upperMDIndex,innerLowerModuleIndex,outerLowerModuleIndex);
segmentsInGPU.mdIndices[idx * 2] = lowerMDIndex;
segmentsInGPU.mdIndices[idx * 2 + 1] = upperMDIndex;
segmentsInGPU.innerLowerModuleIndices[idx] = innerLowerModuleIndex;
segmentsInGPU.outerLowerModuleIndices[idx] = outerLowerModuleIndex;
segmentsInGPU.innerMiniDoubletAnchorHitIndices[idx] = innerMDAnchorHitIndex;
segmentsInGPU.outerMiniDoubletAnchorHitIndices[idx] = outerMDAnchorHitIndex;

segmentsInGPU.dPhis[idx] = __F2H(dPhi);
segmentsInGPU.dPhiMins[idx] = __F2H(dPhiMin);
segmentsInGPU.dPhiMaxs[idx] = __F2H(dPhiMax);
segmentsInGPU.dPhiChanges[idx] = __F2H(dPhiChange);
segmentsInGPU.dPhiChangeMins[idx] = __F2H(dPhiChangeMin);
segmentsInGPU.dPhiChangeMaxs[idx] = __F2H(dPhiChangeMax);
};
#endif

template <typename TAcc>
ALPAKA_FN_ACC ALPAKA_FN_INLINE void addPixelSegmentToMemory(TAcc const& acc,
Expand Down Expand Up @@ -883,8 +974,8 @@ namespace SDL {
} else {
unsigned int segmentModuleIdx =
alpaka::atomicOp<alpaka::AtomicAdd>(acc, &segmentsInGPU.nSegments[innerLowerModuleIndex], 1u);
unsigned int segmentIdx = rangesInGPU.segmentModuleIndices[innerLowerModuleIndex] + segmentModuleIdx;

unsigned int segmentIdx = rangesInGPU.segmentModuleIndices[innerLowerModuleIndex] + segmentModuleIdx;
#ifdef OUTPUT_LS_CUTS
addSegmentToMemory(segmentsInGPU,
innerMDIndex,
outerMDIndex,
Expand All @@ -898,7 +989,32 @@ namespace SDL {
dPhiChange,
dPhiChangeMin,
dPhiChangeMax,
zHi,
zLo,
rtHi,
rtLo,
dAlphaInnerMDSegment,
dAlphaOuterMDSegment,
dAlphaInnerMDOuterMD,
segmentIdx);

#else
addSegmentToMemory(segmentsInGPU,
innerMDIndex,
outerMDIndex,
innerLowerModuleIndex,
outerLowerModuleIndex,
innerMiniDoubletAnchorHitIndex,
outerMiniDoubletAnchorHitIndex,
dPhi,
dPhiMin,
dPhiMax,
dPhiChange,
dPhiChangeMin,
dPhiChangeMax,
segmentIdx);

#endif
}
}
}
Expand Down Expand Up @@ -991,11 +1107,10 @@ namespace SDL {
else {
occupancy = 0;
#ifdef Warnings
printf("Unhandled case in createSegmentArrayRanges! Module index = %i\n", i);
printf("Unhandled case in createSegmentArrayRanges! Module index = %i\n",i);
#endif
}

int nTotSegs = alpaka::atomicOp<alpaka::AtomicAdd>(acc, &nTotalSegments, occupancy);
int nTotSegs = alpaka::atomicOp<alpaka::AtomicAdd>(acc, &nTotalSegments, occupancy);
rangesInGPU.segmentModuleIndices[i] = nTotSegs;
rangesInGPU.segmentModuleOccupancy[i] = occupancy;
}
Expand Down Expand Up @@ -1127,4 +1242,4 @@ namespace SDL {
};
} // namespace SDL

#endif
#endif
88 changes: 85 additions & 3 deletions code/core/write_sdl_ntuple.cc
Original file line number Diff line number Diff line change
@@ -1,5 +1,4 @@
#include "write_sdl_ntuple.h"

//________________________________________________________________________________________________________________________________
void createOutputBranches()
{
Expand Down Expand Up @@ -159,6 +158,12 @@ void createGnnNtupleBranches()
ana.tx->createBranch<vector<float>>("MD_eta");
ana.tx->createBranch<vector<float>>("MD_phi");
ana.tx->createBranch<vector<float>>("MD_dphichange");

//hubert
ana.tx->createBranch<vector<float>>("MD_dphi");
ana.tx->createBranch<vector<float>>("MD_dz");


ana.tx->createBranch<vector<int>>("MD_isFake");
ana.tx->createBranch<vector<int>>("MD_tpType");
ana.tx->createBranch<vector<int>>("MD_detId");
Expand Down Expand Up @@ -756,6 +761,10 @@ void setGnnNtupleMiniDoublet(SDL::Event* event, unsigned int MD)

// Obtaining dPhiChange
float dphichange = miniDoubletsInGPU.dphichanges[MD];

//Huberts addition:
float dphi = miniDoubletsInGPU.dphis[MD];
float dz = miniDoubletsInGPU.dzs[MD];

// Computing pt
const float kRinv1GeVf = (2.99792458e-3 * 3.8);
Expand Down Expand Up @@ -784,6 +793,11 @@ void setGnnNtupleMiniDoublet(SDL::Event* event, unsigned int MD)
ana.tx->pushbackToBranch<float>("MD_1_x", hit1_x);
ana.tx->pushbackToBranch<float>("MD_1_y", hit1_y);
ana.tx->pushbackToBranch<float>("MD_1_z", hit1_z);

// hubert
ana.tx->pushbackToBranch<float>("MD_dz", dz);
ana.tx->pushbackToBranch<float>("MD_dphi", dphi);

// ana.tx->pushbackToBranch<int>("MD_sim_idx", simidxs.size() > 0 ? simidxs[0] : -999);
}

Expand Down Expand Up @@ -1448,7 +1462,14 @@ void createOutputBranches_v2()
//
ana.tx->createBranch<vector<float>> ("md_pt"); // pt (computed based on delta phi change)
ana.tx->createBranch<vector<float>> ("md_eta"); // eta (computed based on anchor hit's eta)
ana.tx->createBranch<vector<float>> ("md_phi"); // phi (computed based on anchor hit's phi)
ana.tx->createBranch<vector<float>> ("md_phi"); // phi (computed based on anchor hit's phi)

#ifdef OUTPUT_MD_CUTS
ana.tx->createBranch<vector<float>> ("md_dphi");
ana.tx->createBranch<vector<float>> ("md_dphichange");
ana.tx->createBranch<vector<float>> ("md_dz");
#endif

ana.tx->createBranch<vector<float>> ("md_anchor_x"); // anchor hit x
ana.tx->createBranch<vector<float>> ("md_anchor_y"); // anchor hit y
ana.tx->createBranch<vector<float>> ("md_anchor_z"); // anchor hit z
Expand All @@ -1475,6 +1496,23 @@ void createOutputBranches_v2()
ana.tx->createBranch<vector<int>> ("ls_mdIdx1"); // index to the second MD
ana.tx->createBranch<vector<int>> ("ls_isFake"); // 1 if md is fake 0 other if not
ana.tx->createBranch<vector<int>> ("ls_simIdx"); // idx of best matched (highest nhit and > 75%) simulated track

#ifdef OUTPUT_LS_CUTS
ana.tx->createBranch<vector<float>> ("ls_zLos");
ana.tx->createBranch<vector<float>> ("ls_zHis");
ana.tx->createBranch<vector<float>> ("ls_rtLos");
ana.tx->createBranch<vector<float>> ("ls_rtHis");
ana.tx->createBranch<vector<float>> ("ls_dPhis");
ana.tx->createBranch<vector<float>> ("ls_dPhiMins");
ana.tx->createBranch<vector<float>> ("ls_dPhiMaxs");
ana.tx->createBranch<vector<float>> ("ls_dPhiChanges");
ana.tx->createBranch<vector<float>> ("ls_dPhiChangeMins");
ana.tx->createBranch<vector<float>> ("ls_dPhiChangeMaxs");
ana.tx->createBranch<vector<float>> ("ls_dAlphaInners");
ana.tx->createBranch<vector<float>> ("ls_dAlphaOuters");
ana.tx->createBranch<vector<float>> ("ls_dAlphaInnerOuters");
#endif

ana.tx->createBranch<vector<vector<int>>> ("ls_simIdxAll"); // list of idx of all matched (> 0%) simulated track
ana.tx->createBranch<vector<vector<float>>>("ls_simIdxAllFrac"); // list of idx of all matched (> 0%) simulated track

Expand Down Expand Up @@ -1859,6 +1897,8 @@ void fillOutputBranches_v2(SDL::Event* event)

// Pt is computed via dphichange and the eta and phi are computed based on anchor hit
float dphichange = miniDoubletsInGPU.dphichanges[mdIdx];
float dphi = miniDoubletsInGPU.dphis[mdIdx];
float dz = miniDoubletsInGPU.dzs[mdIdx];
float k2Rinv1GeVf = (2.99792458e-3 * 3.8) / 2;
float pt = anchor_hit.rt() * k2Rinv1GeVf / sin(dphichange);
float eta = anchor_hit.eta();
Expand All @@ -1876,6 +1916,12 @@ void fillOutputBranches_v2(SDL::Event* event)
ana.tx->pushbackToBranch<float>("md_pt", pt);
ana.tx->pushbackToBranch<float>("md_eta", eta);
ana.tx->pushbackToBranch<float>("md_phi", phi);

#ifdef OUTPUT_MD_CUTS
ana.tx->pushbackToBranch<float>("md_dphichange", dphichange);
ana.tx->pushbackToBranch<float>("md_dphi", dphi);
ana.tx->pushbackToBranch<float>("md_dz", dz);
#endif
ana.tx->pushbackToBranch<float>("md_anchor_x", anchor_x);
ana.tx->pushbackToBranch<float>("md_anchor_y", anchor_y);
ana.tx->pushbackToBranch<float>("md_anchor_z", anchor_z);
Expand Down Expand Up @@ -2016,11 +2062,47 @@ void fillOutputBranches_v2(SDL::Event* event)
float pt = SDL::CPU::MathUtil::ptEstimateFromRadius(center.rt());
float eta = hitC.eta();
float phi = hitB.phi();


#ifdef OUTPUT_LS_CUTS
float zHi = segmentsInGPU.zHis[lsIdx];
float zLo = segmentsInGPU.zLos[lsIdx];
float rtHi = segmentsInGPU.rtHis[lsIdx];
float rtLo = segmentsInGPU.rtLos[lsIdx];
float dAlphaInner = segmentsInGPU.dAlphaInners[lsIdx];
float dAlphaOuter = segmentsInGPU.dAlphaOuters[lsIdx];
float dAlphaInnerOuter = segmentsInGPU.dAlphaInnerOuters[lsIdx];
#endif
float dPhi = segmentsInGPU.dPhis[lsIdx];
float dPhiMin = segmentsInGPU.dPhiMins[lsIdx];
float dPhiMax = segmentsInGPU.dPhiMaxs[lsIdx];
float dPhiChange = segmentsInGPU.dPhiChanges[lsIdx];
float dPhiChangeMin = segmentsInGPU.dPhiChangeMins[lsIdx];
float dPhiChangeMax = segmentsInGPU.dPhiChangeMaxs[lsIdx];

// Write out the ntuple
ana.tx->pushbackToBranch<float>("ls_pt", pt);
ana.tx->pushbackToBranch<float>("ls_eta", eta);
ana.tx->pushbackToBranch<float>("ls_phi", phi);

#ifdef OUTPUT_LS_CUTS
ana.tx->pushbackToBranch<float>("ls_zHis",zHi);
ana.tx->pushbackToBranch<float>("ls_zLos",zLo);
ana.tx->pushbackToBranch<float>("ls_rtHis",rtHi);
ana.tx->pushbackToBranch<float>("ls_rtLos",rtLo);
ana.tx->pushbackToBranch<float>("ls_dPhis",dPhi);
ana.tx->pushbackToBranch<float>("ls_dPhiMins",dPhiMin);
ana.tx->pushbackToBranch<float>("ls_dPhiMaxs",dPhiMax);
ana.tx->pushbackToBranch<float>("ls_dPhiChanges",dPhiChange);
ana.tx->pushbackToBranch<float>("ls_dPhiChangeMins",dPhiChangeMin);
ana.tx->pushbackToBranch<float>("ls_dPhiChangeMaxs",dPhiChangeMax);
ana.tx->pushbackToBranch<float>("ls_dAlphaInners",dAlphaInner);
ana.tx->pushbackToBranch<float>("ls_dAlphaOuters",dAlphaOuter);
ana.tx->pushbackToBranch<float>("ls_dAlphaInnerOuters",dAlphaInnerOuter);

#endif


ana.tx->pushbackToBranch<int>("ls_mdIdx0", md_idx_map[mdIdxs[0]]);
ana.tx->pushbackToBranch<int>("ls_mdIdx1", md_idx_map[mdIdxs[1]]);

Expand Down Expand Up @@ -2935,4 +3017,4 @@ void fillOutputBranches_v2(SDL::Event* event)
// Then clear the branches to default values (e.g. -999, or clear the vectors to empty vectors)
ana.tx->clear();

}
}
Loading