From 646f938c93a404e5748fc00beb4ac0053b2c5656 Mon Sep 17 00:00:00 2001 From: YonsiG Date: Tue, 25 Oct 2022 12:41:34 -0400 Subject: [PATCH 01/25] add print out rz chi2 for T5 inn debug mode --- SDL/Event.cu | 3 +-- SDL/Quintuplet.cu | 23 ++++++++++++----------- SDL/Quintuplet.cuh | 4 ++-- 3 files changed, 15 insertions(+), 15 deletions(-) diff --git a/SDL/Event.cu b/SDL/Event.cu index 49730279..9a3075b1 100644 --- a/SDL/Event.cu +++ b/SDL/Event.cu @@ -1959,10 +1959,9 @@ SDL::quintuplets* SDL::Event::getQuintuplets() quintupletsInCPU->eta = new FPX[nMemoryLocations]; quintupletsInCPU->phi = new FPX[nMemoryLocations]; - quintupletsInCPU->rzChiSquared = new float[nMemoryLocations]; quintupletsInCPU->chiSquared = new float[nMemoryLocations]; quintupletsInCPU->nonAnchorChiSquared = new float[nMemoryLocations]; - + quintupletsInCPU->rzChiSquared = new float[nMemoryLocations]; cudaMemcpyAsync(quintupletsInCPU->nQuintuplets, quintupletsInGPU->nQuintuplets, nLowerModules * sizeof(unsigned int), cudaMemcpyDeviceToHost,stream); cudaMemcpyAsync(quintupletsInCPU->totOccupancyQuintuplets, quintupletsInGPU->totOccupancyQuintuplets, nLowerModules * sizeof(unsigned int), cudaMemcpyDeviceToHost,stream); diff --git a/SDL/Quintuplet.cu b/SDL/Quintuplet.cu index 6af50370..714f5532 100644 --- a/SDL/Quintuplet.cu +++ b/SDL/Quintuplet.cu @@ -283,7 +283,8 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, pass = pass and runQuintupletDefaultAlgo(modulesInGPU, mdsInGPU, segmentsInGPU, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex4, lowerModuleIndex5, firstSegmentIndex, fourthSegmentIndex, firstMDIndex, secondMDIndex, fourthMDIndex, fifthMDIndex, zOut, rtOut, deltaPhiPos, deltaPhi, betaIn, betaOut, pt_beta, zLo, zHi, rtLo, rtHi, zLoPointed, zHiPointed, sdlCut, betaInCut, betaOutCut, deltaBetaCut, kZ); if(not pass) return pass; - pass = pass and passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5); + pass = pass and passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared); + if(not pass) return pass; float x1 = mdsInGPU.anchorX[firstMDIndex]; @@ -517,13 +518,13 @@ __device__ bool SDL::passChiSquaredConstraint(struct SDL::modules& modulesInGPU, } else if(layer1 == 1 and layer2 == 2 and layer3 == 3 and layer4 == 4) { - if(layer5 == 12) + if(layer5 == 5) { - return chiSquared < 0.09461f; + return chiSquared < 0.04725f; } - else if(layer5 == 5) + else if(layer5 == 12) { - return chiSquared < 0.04725f; + return chiSquared < 0.09461f; } } else if(layer1 == 2 and layer2 == 7 and layer3 == 8) @@ -554,17 +555,17 @@ __device__ bool SDL::passChiSquaredConstraint(struct SDL::modules& modulesInGPU, } else if(layer1 == 2 and layer2 == 3 and layer3 == 4) { - if(layer4 == 12 and layer5 == 13) + if(layer4 == 5 and layer5 == 6) { - return chiSquared < 0.10870f; + return chiSquared < 0.08234f; } else if(layer4 == 5 and layer5 == 12) { return chiSquared < 0.10870f; } - else if(layer4 == 5 and layer5 == 6) + else if(layer4 == 12 and layer5 == 13) { - return chiSquared < 0.08234f; + return chiSquared < 0.10870f; } } else if(layer1 == 3 and layer2 == 7 and layer3 == 8 and layer4 == 14 and layer5 == 15) @@ -580,7 +581,7 @@ __device__ bool SDL::passChiSquaredConstraint(struct SDL::modules& modulesInGPU, } //bounds can be found at http://uaf-10.t2.ucsd.edu/~bsathian/SDL/T5_RZFix/t5_rz_thresholds.txt -__device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex, unsigned int fifthMDIndex, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5) +__device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex, unsigned int fifthMDIndex, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, float& rzChiSquared) { const float& rt1 = mdsInGPU.anchorRt[firstMDIndex]; const float& rt2 = mdsInGPU.anchorRt[secondMDIndex]; @@ -626,7 +627,7 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc residual5 = (moduleLayer5 == 0) ? residual5/2.4f : residual5/5.0f; const float RMSE = sqrtf(0.5 * (residual4 * residual4 + residual5 * residual5)); - + rzChiSquared = RMSE; //categories! if(layer1 == 1 and layer2 == 2 and layer3 == 3) { diff --git a/SDL/Quintuplet.cuh b/SDL/Quintuplet.cuh index 6ab1474f..e397fae7 100644 --- a/SDL/Quintuplet.cuh +++ b/SDL/Quintuplet.cuh @@ -70,7 +70,7 @@ namespace SDL CUDA_DEV bool runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, struct SDL::segments& segmentsInGPU, struct SDL::triplets& tripletsInGPU, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, unsigned int& innerTripletIndex, unsigned int& outerTripletIndex, float& innerRadius, float& outerRadius, float& bridgeRadius, float& regressionG, float& regressionF, float& regressionRadius, float& rzChiSquared, float& chiSquared, float& nonAnchorChiSquared); - CUDA_DEV bool passT5RZConstraint(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex, unsigned int fifthMDIndex, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5); + CUDA_DEV bool passT5RZConstraint(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex, unsigned int fifthMDIndex, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, float& rzChiSquared); CUDA_DEV bool T5HasCommonMiniDoublet(struct triplets& tripletsInGPU, struct segments& segmentsInGPU, unsigned int innerTripletIndex, unsigned int outerTripletIndex); @@ -88,7 +88,7 @@ namespace SDL CUDA_DEV bool matchRadiiBBBEE23478(const float& innerRadius, const float& bridgeRadius, const float& outerRadius, const float& innerRadiusMin2S, const float& innerRadiusMax2S, const float& bridgeRadiusMin2S, const float& bridgeRadiusMax2S, const float& outerRadiusMin2S, const float& outerRadiusMax2S,float& innerRadiusMin, float& innerRadiusMax, float& bridgeRadiusMin, float& bridgeRadiusMax, float& outerRadiusMin, float& outerRadiusMax); -CUDA_DEV bool matchRadiiBBBEE34578(const float& innerRadius, const float& bridgeRadius, const float& outerRadius, const float& innerRadiusMin2S, const float& innerRadiusMax2S, const float& bridgeRadiusMin2S, const float& bridgeRadiusMax2S, const float& outerRadiusMin2S, const float& outerRadiusMax2S,float& innerRadiusMin, float& innerRadiusMax, float& bridgeRadiusMin, float& bridgeRadiusMax, float& outerRadiusMin, float& outerRadiusMax); + CUDA_DEV bool matchRadiiBBBEE34578(const float& innerRadius, const float& bridgeRadius, const float& outerRadius, const float& innerRadiusMin2S, const float& innerRadiusMax2S, const float& bridgeRadiusMin2S, const float& bridgeRadiusMax2S, const float& outerRadiusMin2S, const float& outerRadiusMax2S,float& innerRadiusMin, float& innerRadiusMax, float& bridgeRadiusMin, float& bridgeRadiusMax, float& outerRadiusMin, float& outerRadiusMax); CUDA_DEV bool matchRadiiBBEEE(const float& innerRadius, const float& bridgeRadius, const float& outerRadius, const float& innerRadiusMin2S, const float& innerRadiusMax2S, const float& bridgeRadiusMin2S, const float& bridgeRadiusMax2S, const float& outerRadiusMin2S, const float& outerRadiusMax2S,float& innerRadiusMin, float& innerRadiusMax, float& bridgeRadiusMin, float& bridgeRadiusMax, float& outerRadiusMin, float& outerRadiusMax); From 40aa5c0c88806963547f1b6fbad5699893e6fa1d Mon Sep 17 00:00:00 2001 From: YonsiG Date: Fri, 28 Oct 2022 14:10:55 -0400 Subject: [PATCH 02/25] correct the rz chi2 in T5 building --- SDL/Quintuplet.cu | 38 ++++++++++++++++++++++++++++++-------- 1 file changed, 30 insertions(+), 8 deletions(-) diff --git a/SDL/Quintuplet.cu b/SDL/Quintuplet.cu index 714f5532..3266cd03 100644 --- a/SDL/Quintuplet.cu +++ b/SDL/Quintuplet.cu @@ -283,7 +283,14 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, pass = pass and runQuintupletDefaultAlgo(modulesInGPU, mdsInGPU, segmentsInGPU, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex4, lowerModuleIndex5, firstSegmentIndex, fourthSegmentIndex, firstMDIndex, secondMDIndex, fourthMDIndex, fifthMDIndex, zOut, rtOut, deltaPhiPos, deltaPhi, betaIn, betaOut, pt_beta, zLo, zHi, rtLo, rtHi, zLoPointed, zHiPointed, sdlCut, betaInCut, betaOutCut, deltaBetaCut, kZ); if(not pass) return pass; - pass = pass and passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared); + pass = pass and runQuintupletDefaultAlgo(modulesInGPU, mdsInGPU, segmentsInGPU, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex3, lowerModuleIndex4, secondSegmentIndex, thirdSegmentIndex, secondMDIndex, thirdMDIndex, thirdMDIndex, fourthMDIndex, zOut, rtOut, deltaPhiPos, deltaPhi, betaIn, betaOut, pt_beta, zLo, zHi, rtLo, rtHi, zLoPointed, zHiPointed, sdlCut, betaInCut, betaOutCut, deltaBetaCut, kZ); + if(not pass) return pass; + + pass = pass and runQuintupletDefaultAlgo(modulesInGPU, mdsInGPU, segmentsInGPU, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, secondSegmentIndex, fourthSegmentIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, zOut, rtOut, deltaPhiPos, deltaPhi, betaIn, betaOut, pt_beta, zLo, zHi, rtLo, rtHi, zLoPointed, zHiPointed, sdlCut, betaInCut, betaOutCut, deltaBetaCut, kZ); + if(not pass) return pass; + +// pass = pass and passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared); + passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared); if(not pass) return pass; @@ -419,7 +426,7 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, } //compute regression radius right here - this computation is expensive!!! - pass = pass and tempPass; +// pass = pass and tempPass; if(not pass) return pass; float xVec[] = {x1, x2, x3, x4, x5}; @@ -435,7 +442,7 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, //extra chi squared cuts! if(regressionRadius < 5.0f/(2.f * k2Rinv1GeVf)) { - pass = pass and passChiSquaredConstraint(modulesInGPU, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, chiSquared); +// pass = pass and passChiSquaredConstraint(modulesInGPU, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, chiSquared); if(not pass) return pass; } @@ -618,15 +625,26 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc { slope = (z3 - z1)/(rt3 - rt1); } - float residual4 = (layer4 <= 6)? ((z4 - z1) - slope * (rt4 - rt1)) : ((rt4 - rt1) - (z4 - z1)/slope); - float residual5 = (layer4 <= 6) ? ((z5 - z1) - slope * (rt5 - rt1)) : ((rt5 - rt1) - (z5 - z1)/slope); + + //distinguish tilted modules + float side4 = modulesInGPU.sides[lowerModuleIndex4]; + float side5 = modulesInGPU.sides[lowerModuleIndex5]; + float drdz4 = modulesInGPU.drdzs[lowerModuleIndex4]; + float drdz5 = modulesInGPU.drdzs[lowerModuleIndex5]; + + float residual4 = (layer4 <= 6 && ((side4 == SDL::Center) or (drdz4 < 1))) ? ((z4 - z1) - slope * (rt4 - rt1)) : ((rt4 - rt1) - (z4 - z1)/slope); + float residual5 = (layer4 <= 6 && ((side5 == SDL::Center) or (drdz5 < 1))) ? ((z5 - z1) - slope * (rt5 - rt1)) : ((rt5 - rt1) - (z5 - z1)/slope); // creating a chi squared type quantity // 0-> PS, 1->2S - residual4 = (moduleLayer4 == 0) ? residual4/2.4f : residual4/5.0f; - residual5 = (moduleLayer5 == 0) ? residual5/2.4f : residual5/5.0f; + // ps tilted modules correction + float projection4 = ((subdet4 == SDL::Endcap) or (side4 == SDL::Center)) ? 1.f : fmaxf(1.f, drdz4)/sqrtf(1+drdz4*drdz4); + float projection5 = ((subdet5 == SDL::Endcap) or (side5 == SDL::Center)) ? 1.f : fmaxf(1.f, drdz5)/sqrtf(1+drdz5*drdz5); + + residual4 = (moduleType4 == SDL::PS) ? residual4 / (0.15f*projection4) : residual4/5.f; + residual5 = (moduleType5 == SDL::PS) ? residual5 / (0.15f*projection5) : residual5/5.f; - const float RMSE = sqrtf(0.5 * (residual4 * residual4 + residual5 * residual5)); + const float RMSE = sqrtf(residual4 * residual4 + residual5 * residual5); rzChiSquared = RMSE; //categories! if(layer1 == 1 and layer2 == 2 and layer3 == 3) @@ -697,6 +715,10 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc { return RMSE < 0.525f; } + else if(layer4 == 9 and layer5 == 10) + { + printf("here"); + } } else if(layer1 == 2 and layer2 == 7 and layer3 == 13 and layer4 == 14 and layer5 == 15) { From b55ed9a8cd55ae2216fbc66ff85f36648b95c171 Mon Sep 17 00:00:00 2001 From: YonsiG Date: Wed, 2 Nov 2022 09:17:38 -0400 Subject: [PATCH 03/25] add missing point in rz chi2 --- SDL/Quintuplet.cu | 95 ++++++++++++++++++++++++++++------------------- 1 file changed, 56 insertions(+), 39 deletions(-) diff --git a/SDL/Quintuplet.cu b/SDL/Quintuplet.cu index 3266cd03..a33a6688 100644 --- a/SDL/Quintuplet.cu +++ b/SDL/Quintuplet.cu @@ -616,104 +616,121 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc const int moduleLayer4 = modulesInGPU.moduleType[lowerModuleIndex4]; const int moduleLayer5 = modulesInGPU.moduleType[lowerModuleIndex5]; + int missing_points = 0; float slope; + short side_missing; + float drdz_missing; + short subdets_missing; + int moduleType_missing; if(moduleLayer1 == 0 and moduleLayer2 == 0 and moduleLayer3 == 1) //PSPS2S { slope = (z2 -z1)/(rt2 - rt1); + missing_points=3; + side_missing=modulesInGPU.sides[lowerModuleIndex3]; + drdz_missing=modulesInGPU.drdzs[lowerModuleIndex3]; + subdets_missing=modulesInGPU.subdets[lowerModuleIndex3]; + moduleType_missing=modulesInGPU.moduleType[lowerModuleIndex3]; } else { slope = (z3 - z1)/(rt3 - rt1); + missing_points=2; + side_missing=modulesInGPU.sides[lowerModuleIndex2]; + drdz_missing=modulesInGPU.drdzs[lowerModuleIndex2]; + subdets_missing=modulesInGPU.subdets[lowerModuleIndex2]; + moduleType_missing=modulesInGPU.moduleType[lowerModuleIndex2]; } - //distinguish tilted modules - float side4 = modulesInGPU.sides[lowerModuleIndex4]; - float side5 = modulesInGPU.sides[lowerModuleIndex5]; - float drdz4 = modulesInGPU.drdzs[lowerModuleIndex4]; - float drdz5 = modulesInGPU.drdzs[lowerModuleIndex5]; + //correction of tilted modules + float residual_missing; + if (missing_points==2) + { + residual_missing = (layer2 <= 6 && ((side_missing == SDL::Center) or (drdz_missing < 1))) ? ((z2 - z1) - slope * (rt2 - rt1)) : ((rt2 - rt1) - (z2 - z1)/slope); + } + if (missing_points==3) + { + residual_missing = (layer3 <= 6 && ((side_missing == SDL::Center) or (drdz_missing < 1))) ? ((z3 - z1) - slope * (rt3 - rt1)) : ((rt3 - rt1) - (z3 - z1)/slope); + } - float residual4 = (layer4 <= 6 && ((side4 == SDL::Center) or (drdz4 < 1))) ? ((z4 - z1) - slope * (rt4 - rt1)) : ((rt4 - rt1) - (z4 - z1)/slope); - float residual5 = (layer4 <= 6 && ((side5 == SDL::Center) or (drdz5 < 1))) ? ((z5 - z1) - slope * (rt5 - rt1)) : ((rt5 - rt1) - (z5 - z1)/slope); + float projection_missing = ((subdets_missing == SDL::Endcap) or (side_missing == SDL::Center)) ? 1.f : fmaxf(1.f, drdz_missing)/sqrtf(1+drdz_missing*drdz_missing); + residual_missing = (moduleType_missing == SDL::PS) ? residual_missing / (0.15f*projection_missing) : residual_missing/5.f; - // creating a chi squared type quantity - // 0-> PS, 1->2S - // ps tilted modules correction - float projection4 = ((subdet4 == SDL::Endcap) or (side4 == SDL::Center)) ? 1.f : fmaxf(1.f, drdz4)/sqrtf(1+drdz4*drdz4); - float projection5 = ((subdet5 == SDL::Endcap) or (side5 == SDL::Center)) ? 1.f : fmaxf(1.f, drdz5)/sqrtf(1+drdz5*drdz5); + float residual4 = (layer4 <= 6) ? ((z4 - z1) - slope * (rt4 - rt1)) : ((rt4 - rt1) - (z4 - z1)/slope); + float residual5 = (layer4 <= 6) ? ((z5 - z1) - slope * (rt5 - rt1)) : ((rt5 - rt1) - (z5 - z1)/slope); - residual4 = (moduleType4 == SDL::PS) ? residual4 / (0.15f*projection4) : residual4/5.f; - residual5 = (moduleType5 == SDL::PS) ? residual5 / (0.15f*projection5) : residual5/5.f; + residual4 = (moduleLayer4 == 0) ? residual4/0.15f : residual4/5.0f; + residual5 = (moduleLayer5 == 0) ? residual5/0.15f : residual5/5.0f; - const float RMSE = sqrtf(residual4 * residual4 + residual5 * residual5); - rzChiSquared = RMSE; + rzChiSquared = 0.5*12*(residual4 * residual4 + residual5 * residual5 + residual_missing*residual_missing); // 12 is the factor for uniform random variable + //categories! if(layer1 == 1 and layer2 == 2 and layer3 == 3) { if(layer4 == 4 and layer5 == 5) { - return RMSE < 0.545f; + return rzChiSquared < 0.545f; } else if(layer4 == 4 and layer5 == 12) { - return RMSE < 1.105f; + return rzChiSquared < 1.105f; } else if(layer4 == 7 and layer5 == 13) { - return RMSE < 0.775f; + return rzChiSquared < 0.775f; } else if(layer4 == 12 and layer5 == 13) { - return RMSE < 0.625f; + return rzChiSquared < 0.625f; } } else if(layer1 == 1 and layer2 == 2 and layer3 == 7) { if(layer4 == 8 and layer5 == 14) { - return RMSE < 0.835f; + return rzChiSquared < 0.835f; } else if(layer4 == 13 and layer5 == 14) { - return RMSE < 0.575f; + return rzChiSquared < 0.575f; } } else if(layer1 == 1 and layer2 == 7 and layer3 == 8 and layer4 == 9 and layer5 == 15) { - return RMSE < 0.825f; + return rzChiSquared < 0.825f; } else if(layer1 == 2 and layer2 == 3 and layer3 == 4) { if(layer4 == 5 and layer5 == 6) { - return RMSE < 0.845f; + return rzChiSquared < 0.845f; } else if(layer4 == 5 and layer5 == 12) { - return RMSE < 1.365f; + return rzChiSquared < 1.365f; } else if(layer4 == 12 and layer5 == 13) { - return RMSE < 0.675f; + return rzChiSquared < 0.675f; } } else if(layer1 == 2 and layer2 == 3 and layer3 == 7 and layer4 == 13 and layer5 == 14) { - return RMSE < 0.495f; + return rzChiSquared < 0.495f; } else if(layer1 == 2 and layer2 == 3 and layer3 == 12 and layer4 == 13 and layer5 == 14) { - return RMSE < 0.695f; + return rzChiSquared < 0.695f; } else if(layer1 == 2 and layer2 == 7 and layer3 == 8) { if(layer4 == 9 and layer5 == 15) { - return RMSE < 0.735f; + return rzChiSquared < 0.735f; } else if(layer4 == 14 and layer5 == 15) { - return RMSE < 0.525f; + return rzChiSquared < 0.525f; } else if(layer4 == 9 and layer5 == 10) { @@ -722,39 +739,39 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc } else if(layer1 == 2 and layer2 == 7 and layer3 == 13 and layer4 == 14 and layer5 == 15) { - return RMSE < 0.665f; + return rzChiSquared < 0.665f; } else if(layer1 == 3 and layer2 == 4 and layer3 == 5 and layer4 == 12 and layer5 == 13) { - return RMSE < 0.995f; + return rzChiSquared < 0.995f; } else if(layer1 == 3 and layer2 == 4 and layer3 == 12 and layer4 == 13 and layer5 == 14) { - return RMSE < 0.525f; + return rzChiSquared < 0.525f; } else if(layer1 == 3 and layer2 == 7 and layer3 == 8 and layer4 == 14 and layer5 == 15) { - return RMSE < 0.525f; + return rzChiSquared < 0.525f; } else if(layer1 == 3 and layer2 == 7 and layer3 == 13 and layer4 == 14 and layer5 == 15) { - return RMSE < 0.745f; + return rzChiSquared < 0.745f; } else if(layer1 == 3 and layer2 == 12 and layer3 == 13 and layer4 == 14 and layer5 == 15) { - return RMSE < 0.555f; + return rzChiSquared < 0.555f; } else if(layer1 == 7 and layer2 == 8 and layer3 == 9 and layer4 == 15 and layer5 == 16) { - return RMSE < 0.525f; + return rzChiSquared < 0.525f; } else if(layer1 == 7 and layer2 == 8 and layer3 == 14 and layer4 == 15 and layer5 == 16) { - return RMSE < 0.885f; + return rzChiSquared < 0.885f; } else if(layer1 == 7 and layer2 == 13 and layer3 == 14 and layer4 == 15 and layer5 == 16) { - return RMSE < 0.845f; + return rzChiSquared < 0.845f; } return true; From 7fdebd0a8c6f5d5ef7e460b051cc006b97086057 Mon Sep 17 00:00:00 2001 From: YonsiG Date: Wed, 2 Nov 2022 10:46:32 -0400 Subject: [PATCH 04/25] add residual for cut value debug --- SDL/Event.cu | 7 +++++++ SDL/Quintuplet.cu | 47 +++++++++++++++++++++++++++++++++++----------- SDL/Quintuplet.cuh | 10 +++++++--- 3 files changed, 50 insertions(+), 14 deletions(-) diff --git a/SDL/Event.cu b/SDL/Event.cu index 9a3075b1..fb8250c8 100644 --- a/SDL/Event.cu +++ b/SDL/Event.cu @@ -173,6 +173,10 @@ SDL::Event::~Event() delete[] quintupletsInCPU->chiSquared; delete[] quintupletsInCPU->rzChiSquared; delete[] quintupletsInCPU->nonAnchorChiSquared; + delete[] quintupletsInCPU->residual_missing; + delete[] quintupletsInCPU->residual4; + delete[] quintupletsInCPU->residual5; +#endif delete quintupletsInCPU; } @@ -1962,6 +1966,9 @@ SDL::quintuplets* SDL::Event::getQuintuplets() quintupletsInCPU->chiSquared = new float[nMemoryLocations]; quintupletsInCPU->nonAnchorChiSquared = new float[nMemoryLocations]; quintupletsInCPU->rzChiSquared = new float[nMemoryLocations]; + quintupletsInCPU->residual_missing = new float[nMemoryLocations]; + quintupletsInCPU->residual4 = new float[nMemoryLocations]; + quintupletsInCPU->residual5 = new float[nMemoryLocations]; cudaMemcpyAsync(quintupletsInCPU->nQuintuplets, quintupletsInGPU->nQuintuplets, nLowerModules * sizeof(unsigned int), cudaMemcpyDeviceToHost,stream); cudaMemcpyAsync(quintupletsInCPU->totOccupancyQuintuplets, quintupletsInGPU->totOccupancyQuintuplets, nLowerModules * sizeof(unsigned int), cudaMemcpyDeviceToHost,stream); diff --git a/SDL/Quintuplet.cu b/SDL/Quintuplet.cu index a33a6688..026e020b 100644 --- a/SDL/Quintuplet.cu +++ b/SDL/Quintuplet.cu @@ -26,6 +26,10 @@ SDL::quintuplets::quintuplets() chiSquared = nullptr; rzChiSquared = nullptr; nonAnchorChiSquared = nullptr; + residual_missing = nullptr; + residual4 = nullptr; + residual5 = nullptr; +#endif } SDL::quintuplets::~quintuplets() @@ -56,6 +60,10 @@ void SDL::quintuplets::freeMemoryCache() cms::cuda::free_device(dev, rzChiSquared); cms::cuda::free_device(dev, chiSquared); cms::cuda::free_device(dev, nonAnchorChiSquared); + cms::cuda::free_device(dev, residual_missing); + cms::cuda::free_device(dev, residual4); + cms::cuda::free_device(dev, residual5); +#endif } void SDL::quintuplets::freeMemory(cudaStream_t stream) @@ -80,7 +88,11 @@ void SDL::quintuplets::freeMemory(cudaStream_t stream) cudaFree(rzChiSquared); cudaFree(chiSquared); cudaFree(nonAnchorChiSquared); - cudaStreamSynchronize(stream); + cudaFree(residual_missing); + cudaFree(residual4); + cudaFree(residual5); +#endif +cudaStreamSynchronize(stream); } //TODO:Reuse the track candidate one instead of this! __global__ void SDL::createEligibleModulesListForQuintupletsGPU(struct modules& modulesInGPU,struct triplets& tripletsInGPU, struct objectRanges& rangesInGPU) @@ -173,6 +185,10 @@ void SDL::createQuintupletsInExplicitMemory(struct SDL::quintuplets& quintuplets quintupletsInGPU.rzChiSquared = (float*)cms::cuda::allocate_device(dev, nTotalQuintuplets * sizeof(float), stream); quintupletsInGPU.chiSquared = (float*)cms::cuda::allocate_device(dev, nTotalQuintuplets * sizeof(float), stream); quintupletsInGPU.nonAnchorChiSquared = (float*)cms::cuda::allocate_device(dev, nTotalQuintuplets * sizeof(float), stream); + quintupletsInGPU.residual_missing = (float*)cms::cuda::allocate_device(dev, nTotalQuintuplets * sizeof(float), stream); + quintupletsInGPU.residual4 = (float*)cms::cuda::allocate_device(dev, nTotalQuintuplets * sizeof(float), stream); + quintupletsInGPU.residual5 = (float*)cms::cuda::allocate_device(dev, nTotalQuintuplets * sizeof(float), stream); +#endif #else cudaMalloc(&quintupletsInGPU.tripletIndices, 2 * nTotalQuintuplets * sizeof(unsigned int)); cudaMalloc(&quintupletsInGPU.lowerModuleIndices, 5 * nTotalQuintuplets * sizeof(uint16_t)); @@ -194,6 +210,9 @@ void SDL::createQuintupletsInExplicitMemory(struct SDL::quintuplets& quintuplets cudaMalloc(&quintupletsInGPU.rzChiSquared, nTotalQuintuplets * sizeof(float)); cudaMalloc(&quintupletsInGPU.chiSquared, nTotalQuintuplets * sizeof(float)); cudaMalloc(&quintupletsInGPU.nonAnchorChiSquared, nTotalQuintuplets * sizeof(float)); + cudaMalloc(&quintupletsInGPU.residual_missing, nTotalQuintuplets * sizeof(float)); + cudaMalloc(&quintupletsInGPU.residual4, nTotalQuintuplets * sizeof(float)); + cudaMalloc(&quintupletsInGPU.residual5, nTotalQuintuplets * sizeof(float)); cudaMalloc(&quintupletsInGPU.nMemoryLocations, sizeof(unsigned int)); #endif cudaMemsetAsync(quintupletsInGPU.nQuintuplets,0,nLowerModules * sizeof(unsigned int),stream); @@ -208,7 +227,7 @@ void SDL::createQuintupletsInExplicitMemory(struct SDL::quintuplets& quintuplets __device__ void SDL::addQuintupletToMemory(struct SDL::triplets& tripletsInGPU, struct SDL::quintuplets& quintupletsInGPU, unsigned int innerTripletIndex, unsigned int outerTripletIndex, uint16_t& lowerModule1, uint16_t& lowerModule2, uint16_t& lowerModule3, uint16_t& lowerModule4, uint16_t& lowerModule5, float& innerRadius, float& bridgeRadius, float& outerRadius, float& regressionG, float& regressionF, float& regressionRadius, float& rzChiSquared, float& rPhiChiSquared, float& - nonAnchorChiSquared, float pt, float eta, float phi, float scores, uint8_t layer, unsigned int quintupletIndex) + nonAnchorChiSquared, float residual_missing, float residual4, float residual5, float pt, float eta, float phi, float scores, uint8_t layer, unsigned int quintupletIndex) { quintupletsInGPU.tripletIndices[2 * quintupletIndex] = innerTripletIndex; @@ -250,10 +269,12 @@ __device__ void SDL::addQuintupletToMemory(struct SDL::triplets& tripletsInGPU, quintupletsInGPU.rzChiSquared[quintupletIndex] = rzChiSquared; quintupletsInGPU.chiSquared[quintupletIndex] = rPhiChiSquared; quintupletsInGPU.nonAnchorChiSquared[quintupletIndex] = nonAnchorChiSquared; - + quintupletsInGPU.residual_missing[quintupletIndex] = residual_missing; + quintupletsInGPU.residual4[quintupletIndex] = residual4; + quintupletsInGPU.residual5[quintupletIndex] = residual5; } -__device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, struct SDL::segments& segmentsInGPU, struct SDL::triplets& tripletsInGPU, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, unsigned int& innerTripletIndex, unsigned int& outerTripletIndex, float& innerRadius, float& outerRadius, float& bridgeRadius, float& regressionG, float& regressionF, float& regressionRadius, float& rzChiSquared, float& chiSquared, float& nonAnchorChiSquared) +__device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, struct SDL::segments& segmentsInGPU, struct SDL::triplets& tripletsInGPU, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, unsigned int& innerTripletIndex, unsigned int& outerTripletIndex, float& innerRadius, float& outerRadius, float& bridgeRadius, float& regressionG, float& regressionF, float& regressionRadius, float& rzChiSquared, float& chiSquared, float& nonAnchorChiSquared, float& residual_missing, float& residual4, float& residual5) { bool pass = true; unsigned int firstSegmentIndex = tripletsInGPU.segmentIndices[2 * innerTripletIndex]; @@ -290,7 +311,7 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, if(not pass) return pass; // pass = pass and passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared); - passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared); + passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5); if(not pass) return pass; @@ -588,7 +609,7 @@ __device__ bool SDL::passChiSquaredConstraint(struct SDL::modules& modulesInGPU, } //bounds can be found at http://uaf-10.t2.ucsd.edu/~bsathian/SDL/T5_RZFix/t5_rz_thresholds.txt -__device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex, unsigned int fifthMDIndex, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, float& rzChiSquared) +__device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex, unsigned int fifthMDIndex, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, float& rzChiSquared, float& residual_missing, float& residual4, float& residual5) { const float& rt1 = mdsInGPU.anchorRt[firstMDIndex]; const float& rt2 = mdsInGPU.anchorRt[secondMDIndex]; @@ -642,7 +663,7 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc } //correction of tilted modules - float residual_missing; +// float residual_missing; if (missing_points==2) { residual_missing = (layer2 <= 6 && ((side_missing == SDL::Center) or (drdz_missing < 1))) ? ((z2 - z1) - slope * (rt2 - rt1)) : ((rt2 - rt1) - (z2 - z1)/slope); @@ -654,9 +675,13 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc float projection_missing = ((subdets_missing == SDL::Endcap) or (side_missing == SDL::Center)) ? 1.f : fmaxf(1.f, drdz_missing)/sqrtf(1+drdz_missing*drdz_missing); residual_missing = (moduleType_missing == SDL::PS) ? residual_missing / (0.15f*projection_missing) : residual_missing/5.f; - +/* float residual4 = (layer4 <= 6) ? ((z4 - z1) - slope * (rt4 - rt1)) : ((rt4 - rt1) - (z4 - z1)/slope); float residual5 = (layer4 <= 6) ? ((z5 - z1) - slope * (rt5 - rt1)) : ((rt5 - rt1) - (z5 - z1)/slope); +*/ + + residual4 = (layer4 <= 6) ? ((z4 - z1) - slope * (rt4 - rt1)) : ((rt4 - rt1) - (z4 - z1)/slope); + residual5 = (layer4 <= 6) ? ((z5 - z1) - slope * (rt5 - rt1)) : ((rt5 - rt1) - (z5 - z1)/slope); residual4 = (moduleLayer4 == 0) ? residual4/0.15f : residual4/5.0f; residual5 = (moduleLayer5 == 0) ? residual5/0.15f : residual5/5.0f; @@ -1320,9 +1345,9 @@ __global__ void SDL::createQuintupletsInGPUv2(struct SDL::modules& modulesInGPU, uint16_t lowerModule4 = tripletsInGPU.lowerModuleIndices[3 * outerTripletIndex + 1]; uint16_t lowerModule5 = tripletsInGPU.lowerModuleIndices[3 * outerTripletIndex + 2]; - float innerRadius, outerRadius, bridgeRadius, regressionG, regressionF, regressionRadius, rzChiSquared, chiSquared, nonAnchorChiSquared; //required for making distributions + float innerRadius, outerRadius, bridgeRadius, regressionG, regressionF, regressionRadius, rzChiSquared, chiSquared, nonAnchorChiSquared, residual_missing, residual4, residual; //required for making distributions - bool success = runQuintupletDefaultAlgo(modulesInGPU, mdsInGPU, segmentsInGPU, tripletsInGPU, lowerModule1, lowerModule2, lowerModule3, lowerModule4, lowerModule5, innerTripletIndex, outerTripletIndex, innerRadius, outerRadius, bridgeRadius, regressionG, regressionF, regressionRadius, rzChiSquared, chiSquared, nonAnchorChiSquared); + bool success = runQuintupletDefaultAlgo(modulesInGPU, mdsInGPU, segmentsInGPU, tripletsInGPU, lowerModule1, lowerModule2, lowerModule3, lowerModule4, lowerModule5, innerTripletIndex, outerTripletIndex, innerRadius, outerRadius, bridgeRadius, regressionG, regressionF, regressionRadius, rzChiSquared, chiSquared, nonAnchorChiSquared, residual_missing, residual4, residual5); if(success) { @@ -1362,7 +1387,7 @@ __global__ void SDL::createQuintupletsInGPUv2(struct SDL::modules& modulesInGPU, float eta = mdsInGPU.anchorEta[segmentsInGPU.mdIndices[2*tripletsInGPU.segmentIndices[2*innerTripletIndex+layer2_adjustment]]]; float pt = (innerRadius+outerRadius)*3.8f*1.602f/(2*100*5.39f); float scores = chiSquared + nonAnchorChiSquared; - addQuintupletToMemory(tripletsInGPU, quintupletsInGPU, innerTripletIndex, outerTripletIndex, lowerModule1, lowerModule2, lowerModule3, lowerModule4, lowerModule5, innerRadius, bridgeRadius, outerRadius, regressionG, regressionF, regressionRadius, rzChiSquared, chiSquared, nonAnchorChiSquared, pt,eta,phi,scores,layer,quintupletIndex); + addQuintupletToMemory(tripletsInGPU, quintupletsInGPU, innerTripletIndex, outerTripletIndex, lowerModule1, lowerModule2, lowerModule3, lowerModule4, lowerModule5, innerRadius, bridgeRadius, outerRadius, regressionG, regressionF, regressionRadius, rzChiSquared, chiSquared, nonAnchorChiSquared, rzChiSquared, residual_missing, residual4, residual5, pt,eta,phi,scores,layer,quintupletIndex); tripletsInGPU.partOfT5[quintupletsInGPU.tripletIndices[2 * quintupletIndex]] = true; tripletsInGPU.partOfT5[quintupletsInGPU.tripletIndices[2 * quintupletIndex + 1]] = true; diff --git a/SDL/Quintuplet.cuh b/SDL/Quintuplet.cuh index e397fae7..1a4e1adb 100644 --- a/SDL/Quintuplet.cuh +++ b/SDL/Quintuplet.cuh @@ -51,6 +51,10 @@ namespace SDL float* rzChiSquared; float* chiSquared; float* nonAnchorChiSquared; + float* residual_missing; + float* residual4; + float* residual5; +#endif quintuplets(); ~quintuplets(); @@ -65,12 +69,12 @@ namespace SDL // CUDA_DEV void rmQuintupletToMemory(struct SDL::quintuplets& quintupletsInGPU, unsigned int quintupletIndex); - CUDA_DEV void addQuintupletToMemory(struct SDL::triplets& tripletsInGPU, struct SDL::quintuplets& quintupletsInGPU, unsigned int innerTripletIndex, unsigned int outerTripletIndex, uint16_t& lowerModule1, uint16_t& lowerModule2, uint16_t& lowerModule3, uint16_t& lowerModule4, uint16_t& lowerModule5, float& innerRadius, float& bridgeRadius, float& outerRadius, float& regressionG, float& regressionF, float& regressionRadius, float& rzChiSquared, float& rPhiChiSquared, float& nonAnchorChiSquared, float pt, float eta, float phi, float scores, uint8_t layer, unsigned int quintupletIndex); + CUDA_DEV void addQuintupletToMemory(struct SDL::triplets& tripletsInGPU, struct SDL::quintuplets& quintupletsInGPU, unsigned int innerTripletIndex, unsigned int outerTripletIndex, uint16_t& lowerModule1, uint16_t& lowerModule2, uint16_t& lowerModule3, uint16_t& lowerModule4, uint16_t& lowerModule5, float& innerRadius, float& bridgeRadius, float& outerRadius, float& regressionG, float& regressionF, float& regressionRadius, float& rzChiSquared, float& rPhiChiSquared, float& nonAnchorChiSquared, float residual_missing, float residual4, float residual5, float pt, float eta, float phi, float scores, uint8_t layer, unsigned int quintupletIndex); - CUDA_DEV bool runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, struct SDL::segments& segmentsInGPU, struct SDL::triplets& tripletsInGPU, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, unsigned int& innerTripletIndex, unsigned int& outerTripletIndex, float& innerRadius, float& outerRadius, float& bridgeRadius, float& regressionG, float& regressionF, float& regressionRadius, float& rzChiSquared, float& chiSquared, float& nonAnchorChiSquared); + CUDA_DEV bool runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, struct SDL::segments& segmentsInGPU, struct SDL::triplets& tripletsInGPU, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, unsigned int& innerTripletIndex, unsigned int& outerTripletIndex, float& innerRadius, float& outerRadius, float& bridgeRadius, float& regressionG, float& regressionF, float& regressionRadius, float& rzChiSquared, float& chiSquared, float& nonAnchorChiSquared, float& residual_missing, float& residual4, float& residual5); - CUDA_DEV bool passT5RZConstraint(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex, unsigned int fifthMDIndex, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, float& rzChiSquared); + CUDA_DEV bool passT5RZConstraint(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex, unsigned int fifthMDIndex, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, float& rzChiSquared, float& residual_missing, float& residual4, float& residual5); CUDA_DEV bool T5HasCommonMiniDoublet(struct triplets& tripletsInGPU, struct segments& segmentsInGPU, unsigned int innerTripletIndex, unsigned int outerTripletIndex); From 30705b050b09cc2808168e5387d530eb7a0bef41 Mon Sep 17 00:00:00 2001 From: YonsiG Date: Thu, 1 Dec 2022 08:38:50 -0500 Subject: [PATCH 05/25] add helix approximation for T5 rz --- SDL/PixelTriplet.cu | 2 - SDL/Quintuplet.cu | 197 ++++++++++++++++++++++++++++---------------- SDL/Quintuplet.cuh | 2 +- 3 files changed, 128 insertions(+), 73 deletions(-) diff --git a/SDL/PixelTriplet.cu b/SDL/PixelTriplet.cu index 2e44ebf9..4ef569dc 100644 --- a/SDL/PixelTriplet.cu +++ b/SDL/PixelTriplet.cu @@ -474,10 +474,8 @@ __device__ bool SDL::passPT3RZChiSquaredCuts(struct modules& modulesInGPU, uint1 __device__ float SDL::computePT3RZChiSquared(struct modules& modulesInGPU, uint16_t* lowerModuleIndices, float* rtPix, float* xPix, float* yPix, float* zPix, float* rts, float* xs, float* ys, float* zs, float pixelSegmentPt, float pixelSegmentPx, float pixelSegmentPy, float pixelSegmentPz, int pixelSegmentCharge) { - float slope = (zPix[1] - zPix[0])/(rtPix[1] - rtPix[0]); float residual = 0; float error = 0; - //hardcoded array indices!!! float RMSE = 0; float Px=pixelSegmentPx, Py=pixelSegmentPy, Pz=pixelSegmentPz; diff --git a/SDL/Quintuplet.cu b/SDL/Quintuplet.cu index 026e020b..9fcc86e8 100644 --- a/SDL/Quintuplet.cu +++ b/SDL/Quintuplet.cu @@ -310,11 +310,6 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, pass = pass and runQuintupletDefaultAlgo(modulesInGPU, mdsInGPU, segmentsInGPU, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, secondSegmentIndex, fourthSegmentIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, zOut, rtOut, deltaPhiPos, deltaPhi, betaIn, betaOut, pt_beta, zLo, zHi, rtLo, rtHi, zLoPointed, zHiPointed, sdlCut, betaInCut, betaOutCut, deltaBetaCut, kZ); if(not pass) return pass; -// pass = pass and passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared); - passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5); - - if(not pass) return pass; - float x1 = mdsInGPU.anchorX[firstMDIndex]; float x2 = mdsInGPU.anchorX[secondMDIndex]; float x3 = mdsInGPU.anchorX[thirdMDIndex]; @@ -402,6 +397,11 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, outerRadius = computeRadiusFromThreeAnchorHits(x3, y3, x4, y4, x5, y5, g, f); bridgeRadius = computeRadiusFromThreeAnchorHits(x2, y2, x3, y3, x4, y4, g, f); +// pass = pass and passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, g, f); + float inner_pt = k2Rinv1GeVf * (innerRadius+outerRadius); + passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, g, f); + + if(not pass) return pass; pass = pass & (innerRadius >= 0.95f * ptCut/(2.f * k2Rinv1GeVf)); @@ -609,19 +609,20 @@ __device__ bool SDL::passChiSquaredConstraint(struct SDL::modules& modulesInGPU, } //bounds can be found at http://uaf-10.t2.ucsd.edu/~bsathian/SDL/T5_RZFix/t5_rz_thresholds.txt -__device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex, unsigned int fifthMDIndex, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, float& rzChiSquared, float& residual_missing, float& residual4, float& residual5) +__device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex, unsigned int fifthMDIndex, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, float& rzChiSquared, float& residual_missing, float& residual4, float& residual5, float inner_pt, float g, float f) { - const float& rt1 = mdsInGPU.anchorRt[firstMDIndex]; - const float& rt2 = mdsInGPU.anchorRt[secondMDIndex]; - const float& rt3 = mdsInGPU.anchorRt[thirdMDIndex]; - const float& rt4 = mdsInGPU.anchorRt[fourthMDIndex]; - const float& rt5 = mdsInGPU.anchorRt[fifthMDIndex]; - - const float& z1 = mdsInGPU.anchorZ[firstMDIndex]; - const float& z2 = mdsInGPU.anchorZ[secondMDIndex]; - const float& z3 = mdsInGPU.anchorZ[thirdMDIndex]; - const float& z4 = mdsInGPU.anchorZ[fourthMDIndex]; - const float& z5 = mdsInGPU.anchorZ[fifthMDIndex]; + //(g,f) is the center of the circle fitted by the innermost 3 points on x,y coordinates + const float& rt1 = mdsInGPU.anchorRt[firstMDIndex]/100; //in the unit of m instead of cm + const float& rt2 = mdsInGPU.anchorRt[secondMDIndex]/100; + const float& rt3 = mdsInGPU.anchorRt[thirdMDIndex]/100; + const float& rt4 = mdsInGPU.anchorRt[fourthMDIndex]/100; + const float& rt5 = mdsInGPU.anchorRt[fifthMDIndex]/100; + + const float& z1 = mdsInGPU.anchorZ[firstMDIndex]/100; + const float& z2 = mdsInGPU.anchorZ[secondMDIndex]/100; + const float& z3 = mdsInGPU.anchorZ[thirdMDIndex]/100; + const float& z4 = mdsInGPU.anchorZ[fourthMDIndex]/100; + const float& z5 = mdsInGPU.anchorZ[fifthMDIndex]/100; //following Philip's layer number prescription const int layer1 = modulesInGPU.layers[lowerModuleIndex1] + 6 * (modulesInGPU.subdets[lowerModuleIndex1] == SDL::Endcap) + 5 * (modulesInGPU.subdets[lowerModuleIndex1] == SDL::Endcap and modulesInGPU.moduleType[lowerModuleIndex1] == SDL::TwoS); @@ -631,63 +632,119 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc const int layer5 = modulesInGPU.layers[lowerModuleIndex5] + 6 * (modulesInGPU.subdets[lowerModuleIndex5] == SDL::Endcap) + 5 * (modulesInGPU.subdets[lowerModuleIndex5] == SDL::Endcap and modulesInGPU.moduleType[lowerModuleIndex5] == SDL::TwoS); //slope computed using the internal T3s - const int moduleLayer1 = modulesInGPU.moduleType[lowerModuleIndex1]; - const int moduleLayer2 = modulesInGPU.moduleType[lowerModuleIndex2]; - const int moduleLayer3 = modulesInGPU.moduleType[lowerModuleIndex3]; - const int moduleLayer4 = modulesInGPU.moduleType[lowerModuleIndex4]; - const int moduleLayer5 = modulesInGPU.moduleType[lowerModuleIndex5]; - - int missing_points = 0; - float slope; - short side_missing; - float drdz_missing; - short subdets_missing; - int moduleType_missing; - if(moduleLayer1 == 0 and moduleLayer2 == 0 and moduleLayer3 == 1) //PSPS2S - { - slope = (z2 -z1)/(rt2 - rt1); - missing_points=3; - side_missing=modulesInGPU.sides[lowerModuleIndex3]; - drdz_missing=modulesInGPU.drdzs[lowerModuleIndex3]; - subdets_missing=modulesInGPU.subdets[lowerModuleIndex3]; - moduleType_missing=modulesInGPU.moduleType[lowerModuleIndex3]; - } - else - { - slope = (z3 - z1)/(rt3 - rt1); - missing_points=2; - side_missing=modulesInGPU.sides[lowerModuleIndex2]; - drdz_missing=modulesInGPU.drdzs[lowerModuleIndex2]; - subdets_missing=modulesInGPU.subdets[lowerModuleIndex2]; - moduleType_missing=modulesInGPU.moduleType[lowerModuleIndex2]; - } + const int moduleType1 = modulesInGPU.moduleType[lowerModuleIndex1]; //0 is ps, 1 is 2s + const int moduleType2 = modulesInGPU.moduleType[lowerModuleIndex2]; + const int moduleType3 = modulesInGPU.moduleType[lowerModuleIndex3]; + const int moduleType4 = modulesInGPU.moduleType[lowerModuleIndex4]; + const int moduleType5 = modulesInGPU.moduleType[lowerModuleIndex5]; + + const float& x1 = mdsInGPU.anchorX[firstMDIndex]/100; + const float& x2 = mdsInGPU.anchorX[secondMDIndex]/100; + const float& x3 = mdsInGPU.anchorX[thirdMDIndex]/100; + const float& x4 = mdsInGPU.anchorX[fourthMDIndex]/100; + const float& x5 = mdsInGPU.anchorX[fifthMDIndex]/100; + const float& y1 = mdsInGPU.anchorY[firstMDIndex]/100; + const float& y2 = mdsInGPU.anchorY[secondMDIndex]/100; + const float& y3 = mdsInGPU.anchorY[thirdMDIndex]/100; + const float& y4 = mdsInGPU.anchorY[fourthMDIndex]/100; + const float& y5 = mdsInGPU.anchorY[fifthMDIndex]/100; + + float residual = 0; + float error = 0; + + float x_center=g, y_center=f; + float x_init=mdsInGPU.anchorX[secondMDIndex]/100; + float y_init=mdsInGPU.anchorY[secondMDIndex]/100; + float z_init=mdsInGPU.anchorZ[secondMDIndex]/100; + float rt_init=mdsInGPU.anchorRt[secondMDIndex]/100; //use the second MD as initial point + float pseudo_phi = atan((y_center-y_init)/(x_center-x_init)); //actually represent pi/2-phi, wrt helix axis z + + float Pt=inner_pt, Px = Pt*sin(pseudo_phi), Py=Pt*cos(pseudo_phi); + float Pz=(z_init-z1)/(y_init-y1)*Py; + float p = sqrt(Px*Px+Py*Py+Pz*Pz); + + int charge; + if (((y3-y2)/(x3-x2))<((y2-y1)/(x2-x1))) charge=1; + else if (((y3-y2)/(x3-x2))>((y2-y1)/(x2-x1))) charge=-1; + else return 0; + + float B = 3.8112; + float a = -0.299792*B*charge; + + float zsi, rtsi; + int layeri, moduleTypei; + float expectrt4=0,expectrt5=0,expectz4=0, expectz5=0; + for(size_t i = 4; i < 6; i++) + { + if (i==4){ + zsi = z4; + rtsi = rt4; + layeri=layer4; + moduleTypei=moduleType4; + } + else if (i==5){ + zsi = z5; + rtsi = rt5; + layeri=layer5; + moduleTypei=moduleType5; + } - //correction of tilted modules -// float residual_missing; - if (missing_points==2) - { - residual_missing = (layer2 <= 6 && ((side_missing == SDL::Center) or (drdz_missing < 1))) ? ((z2 - z1) - slope * (rt2 - rt1)) : ((rt2 - rt1) - (z2 - z1)/slope); - } - if (missing_points==3) - { - residual_missing = (layer3 <= 6 && ((side_missing == SDL::Center) or (drdz_missing < 1))) ? ((z3 - z1) - slope * (rt3 - rt1)) : ((rt3 - rt1) - (z3 - z1)/slope); - } + // calculation is copied from PixelTriplet.cu SDL::computePT3RZChiSquared + float diffr=0, diffz=0; + + float rou = a/p; + if (layeri>6){ + float s = (zsi-z_init)*p/Pz; + float x = x_init + Px/a*sin(rou*s)-Py/a*(1-cos(rou*s)); + float y = y_init + Py/a*sin(rou*s)+Px/a*(1-cos(rou*s)); + diffr = (rtsi-sqrt(x*x+y*y))*100; + if (i==4) expectrt4=sqrt(x*x+y*y); + if (i==5) expectrt5=sqrt(x*x+y*y); + } - float projection_missing = ((subdets_missing == SDL::Endcap) or (side_missing == SDL::Center)) ? 1.f : fmaxf(1.f, drdz_missing)/sqrtf(1+drdz_missing*drdz_missing); - residual_missing = (moduleType_missing == SDL::PS) ? residual_missing / (0.15f*projection_missing) : residual_missing/5.f; -/* - float residual4 = (layer4 <= 6) ? ((z4 - z1) - slope * (rt4 - rt1)) : ((rt4 - rt1) - (z4 - z1)/slope); - float residual5 = (layer4 <= 6) ? ((z5 - z1) - slope * (rt5 - rt1)) : ((rt5 - rt1) - (z5 - z1)/slope); -*/ + if (layeri<=6){ + float paraA = rt_init*rt_init + 2*(Px*Px+Py*Py)/(a*a) + 2*(y_init*Px-x_init*Py)/a - rtsi*rtsi; + float paraB = 2*(x_init*Px+y_init*Py)/a; + float paraC = 2*(y_init*Px-x_init*Py)/a+2*(Px*Px+Py*Py)/(a*a); + float A=paraB*paraB+paraC*paraC; + float B=2*paraA*paraB; + float C=paraA*paraA-paraC*paraC; + float sol1 = (-B+sqrt(B*B-4*A*C))/(2*A); + float sol2 = (-B-sqrt(B*B-4*A*C))/(2*A); + float solz1 = asin(sol1)/rou*Pz/p+z_init; + float solz2 = asin(sol2)/rou*Pz/p+z_init; + float diffz1 = (solz1-zsi)*100; + float diffz2 = (solz2-zsi)*100; + diffz = (fabs(diffz1)fabs(diffz2)) expectz4 = solz2; + if (i==5 && fabs(diffz1)fabs(diffz2)) expectz5 = solz2; + } - residual4 = (layer4 <= 6) ? ((z4 - z1) - slope * (rt4 - rt1)) : ((rt4 - rt1) - (z4 - z1)/slope); - residual5 = (layer4 <= 6) ? ((z5 - z1) - slope * (rt5 - rt1)) : ((rt5 - rt1) - (z5 - z1)/slope); + residual = (layeri>6) ? diffr : diffz ; + if (i==4) residual4=residual; + if (i==5) residual5=residual; - residual4 = (moduleLayer4 == 0) ? residual4/0.15f : residual4/5.0f; - residual5 = (moduleLayer5 == 0) ? residual5/0.15f : residual5/5.0f; + //PS Modules + if(moduleTypei == 0) + { + error = 0.15f; + } + else //2S modules + { + error = 5.0f; + } + rzChiSquared += 12*(residual * residual)/(error * error); + } + printf("rt1:%f, rt2:%f, rt3:%f, rt4:%f, rt5:%f\n", rt1, rt2, rt3, rt4, rt5); + printf("x1:%f, x2:%f, x3:%f, x4:%f, x5:%f\n", x1, x2, x3, x4, x5); + printf("y1:%f, y2:%f, y3:%f, y4:%f, y5:%f\n", y1, y2, y3, y4, y5); + printf("z1:%f, z2:%f, z3:%f, z4:%f, z5:%f\n", z1, z2, z3, z4, z5); + printf("rt4_ex:%f, rt5_ex:%f\n", expectrt4, expectrt5); + printf("z4_ex:%f, z5_ex:%f\n", expectz4, expectz5); + printf("Pt:%f, Px:%f, Py:%f, Pz:%f, charge: %i, residual4: %f, residual5:%f\n", Pt, Px, Py, Pz, charge, residual4, residual5); - rzChiSquared = 0.5*12*(residual4 * residual4 + residual5 * residual5 + residual_missing*residual_missing); // 12 is the factor for uniform random variable - //categories! if(layer1 == 1 and layer2 == 2 and layer3 == 3) { @@ -1059,7 +1116,7 @@ __device__ float SDL::computeRadiusFromThreeAnchorHits(float x1, float y1, float /* if((y1 - y3) * (x2 - x3) - (x1 - x3) * (y2 - y3) == 0) { - return -1; //WTF man three collinear points! + return -1; //WTF man, three collinear points! } */ diff --git a/SDL/Quintuplet.cuh b/SDL/Quintuplet.cuh index 1a4e1adb..e4f7d1db 100644 --- a/SDL/Quintuplet.cuh +++ b/SDL/Quintuplet.cuh @@ -74,7 +74,7 @@ namespace SDL CUDA_DEV bool runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, struct SDL::segments& segmentsInGPU, struct SDL::triplets& tripletsInGPU, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, unsigned int& innerTripletIndex, unsigned int& outerTripletIndex, float& innerRadius, float& outerRadius, float& bridgeRadius, float& regressionG, float& regressionF, float& regressionRadius, float& rzChiSquared, float& chiSquared, float& nonAnchorChiSquared, float& residual_missing, float& residual4, float& residual5); - CUDA_DEV bool passT5RZConstraint(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex, unsigned int fifthMDIndex, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, float& rzChiSquared, float& residual_missing, float& residual4, float& residual5); + CUDA_DEV bool passT5RZConstraint(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex, unsigned int fifthMDIndex, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, float& rzChiSquared, float& residual_missing, float& residual4, float& residual5, float inner_pt, float g, float f); CUDA_DEV bool T5HasCommonMiniDoublet(struct triplets& tripletsInGPU, struct segments& segmentsInGPU, unsigned int innerTripletIndex, unsigned int outerTripletIndex); From ff42be4ca77c788d48547a1ab2125fdb6c21f88f Mon Sep 17 00:00:00 2001 From: YonsiG Date: Sun, 4 Dec 2022 12:44:07 -0500 Subject: [PATCH 06/25] fix the sign of initial momentum --- SDL/Quintuplet.cu | 96 +++++++++++++++++++++++++++++++++++----------- SDL/Quintuplet.cuh | 2 +- 2 files changed, 74 insertions(+), 24 deletions(-) diff --git a/SDL/Quintuplet.cu b/SDL/Quintuplet.cu index 9fcc86e8..ff228252 100644 --- a/SDL/Quintuplet.cu +++ b/SDL/Quintuplet.cu @@ -393,13 +393,13 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, computeErrorInRadius(x3Vec, y3Vec, x1Vec, y1Vec, x2Vec, y2Vec, outerRadiusMin2S, outerRadiusMax2S); float g, f; - innerRadius = computeRadiusFromThreeAnchorHits(x1, y1, x2, y2, x3, y3, g, f); outerRadius = computeRadiusFromThreeAnchorHits(x3, y3, x4, y4, x5, y5, g, f); bridgeRadius = computeRadiusFromThreeAnchorHits(x2, y2, x3, y3, x4, y4, g, f); + innerRadius = computeRadiusFromThreeAnchorHits(x1, y1, x2, y2, x3, y3, g, f); -// pass = pass and passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, g, f); - float inner_pt = k2Rinv1GeVf * (innerRadius+outerRadius); - passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, g, f); +// pass = pass and passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f); + float inner_pt = 2 * k2Rinv1GeVf * innerRadius; + passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f); if(not pass) return pass; @@ -609,7 +609,7 @@ __device__ bool SDL::passChiSquaredConstraint(struct SDL::modules& modulesInGPU, } //bounds can be found at http://uaf-10.t2.ucsd.edu/~bsathian/SDL/T5_RZFix/t5_rz_thresholds.txt -__device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex, unsigned int fifthMDIndex, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, float& rzChiSquared, float& residual_missing, float& residual4, float& residual5, float inner_pt, float g, float f) +__device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex, unsigned int fifthMDIndex, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, float& rzChiSquared, float& residual_missing, float& residual4, float& residual5, float inner_pt, float innerRadius, float g, float f) { //(g,f) is the center of the circle fitted by the innermost 3 points on x,y coordinates const float& rt1 = mdsInGPU.anchorRt[firstMDIndex]/100; //in the unit of m instead of cm @@ -651,23 +651,57 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc float residual = 0; float error = 0; + float x_center=g/100, y_center=f/100; + float x_init=mdsInGPU.anchorX[thirdMDIndex]/100; + float y_init=mdsInGPU.anchorY[thirdMDIndex]/100; + float z_init=mdsInGPU.anchorZ[thirdMDIndex]/100; + float rt_init=mdsInGPU.anchorRt[thirdMDIndex]/100; //use the second MD as initial point - float x_center=g, y_center=f; - float x_init=mdsInGPU.anchorX[secondMDIndex]/100; - float y_init=mdsInGPU.anchorY[secondMDIndex]/100; - float z_init=mdsInGPU.anchorZ[secondMDIndex]/100; - float rt_init=mdsInGPU.anchorRt[secondMDIndex]/100; //use the second MD as initial point - float pseudo_phi = atan((y_center-y_init)/(x_center-x_init)); //actually represent pi/2-phi, wrt helix axis z - - float Pt=inner_pt, Px = Pt*sin(pseudo_phi), Py=Pt*cos(pseudo_phi); - float Pz=(z_init-z1)/(y_init-y1)*Py; - float p = sqrt(Px*Px+Py*Py+Pz*Pz); + //start from a circle of inner T3. + // to determine the charge + int charge=0; + float slope1c=(y3-y_center)/(x3-x_center); + float slope3c=(y1-y_center)/(x1-x_center); + if((y3-y_center)>0 && (y1-y_center)>0) + { + if (slope3c>slope1c) charge=-1; + else if (slope3cslope1c) charge=-1; + else if (slope3c5) charge=-1; + else if(slope3c>5 && slope1c<-5) charge=1; - int charge; if (((y3-y2)/(x3-x2))<((y2-y1)/(x2-x1))) charge=1; else if (((y3-y2)/(x3-x2))>((y2-y1)/(x2-x1))) charge=-1; else return 0; + float pseudo_phi = atan((y_init-y_center)/(x_init-x_center)); //actually represent pi/2-phi, wrt helix axis z + float Pt=inner_pt, Px, Py; + if (charge==1 && x_init<0 && y_init<0) {Px = -Pt*abs(sin(pseudo_phi)); Py=-Pt*abs(cos(pseudo_phi));} + else if (charge==1 && x_init<0 && y_init>0) {Px = -Pt*abs(sin(pseudo_phi)); Py=Pt*abs(cos(pseudo_phi));} + else if (charge==1 && x_init>0 && y_init>0) {Px = Pt*abs(sin(pseudo_phi)); Py=Pt*abs(cos(pseudo_phi));} + else if (charge==1 && x_init>0 && y_init<0) {Px = Pt*abs(sin(pseudo_phi)); Py=-Pt*abs(cos(pseudo_phi));} + else if (charge==-1 && x_init>0 && y_init>0) {Px = Pt*abs(sin(pseudo_phi)); Py=Pt*abs(cos(pseudo_phi));} + else if (charge==-1 && x_init>0 && y_init<0) {Px = Pt*abs(sin(pseudo_phi)); Py=-Pt*abs(cos(pseudo_phi));} + else if (charge==-1 && x_init<0 && y_init<0) {Px = -Pt*abs(sin(pseudo_phi)); Py=-Pt*abs(cos(pseudo_phi));} + else if (charge==-1 && x_init<0 && y_init>0) {Px = -Pt*abs(sin(pseudo_phi)); Py=Pt*abs(cos(pseudo_phi));} + else return 0; + + //to get Pz, we use pt/pz=ds/dz, ds is the arclength between MD1 and MD2. + float AO=sqrt((x1-x_center)*(x1-x_center)+(y1-y_center)*(y1-y_center)); + float BO=sqrt((x_init-x_center)*(x_init-x_center)+(y_init-y_center)*(y_init-y_center)); + float AB=sqrt((x1-x_init)*(x1-x_init)+(y1-y_init)*(y1-y_init)); + float dPhi = acos((AO*AO+BO*BO-AB*AB)/(2*AO*BO)); +// float ds=innerRadius/100*dPhi; + + float ds = sqrt((y_init-y1)*(y_init-y1)+(x_init-x1)*(x_init-x1)); //large ds->smallerPz->smaller residual + float Pz=(z_init-z1)/ds*Pt; + float p = sqrt(Px*Px+Py*Py+Pz*Pz); + float B = 3.8112; float a = -0.299792*B*charge; @@ -688,6 +722,18 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc layeri=layer5; moduleTypei=moduleType5; } + else if (i==2){ + zsi = z2; + rtsi = rt2; + layeri=layer2; + moduleTypei=moduleType2; + } + else if (i==3){ + zsi = z3; + rtsi = rt3; + layeri=layer3; + moduleTypei=moduleType3; + } // calculation is copied from PixelTriplet.cu SDL::computePT3RZChiSquared float diffr=0, diffz=0; @@ -737,13 +783,17 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc } rzChiSquared += 12*(residual * residual)/(error * error); } - printf("rt1:%f, rt2:%f, rt3:%f, rt4:%f, rt5:%f\n", rt1, rt2, rt3, rt4, rt5); - printf("x1:%f, x2:%f, x3:%f, x4:%f, x5:%f\n", x1, x2, x3, x4, x5); - printf("y1:%f, y2:%f, y3:%f, y4:%f, y5:%f\n", y1, y2, y3, y4, y5); - printf("z1:%f, z2:%f, z3:%f, z4:%f, z5:%f\n", z1, z2, z3, z4, z5); - printf("rt4_ex:%f, rt5_ex:%f\n", expectrt4, expectrt5); - printf("z4_ex:%f, z5_ex:%f\n", expectz4, expectz5); - printf("Pt:%f, Px:%f, Py:%f, Pz:%f, charge: %i, residual4: %f, residual5:%f\n", Pt, Px, Py, Pz, charge, residual4, residual5); + + if(layer1 == 1 and layer2 == 2 and layer3 == 3 and layer4 == 4 and layer5 == 5 and residual4>15){ +// printf("rt1:%f, rt2:%f, rt3:%f, rt4:%f, rt5:%f\n", rt1, rt2, rt3, rt4, rt5); +// printf("x1:%f, x2:%f, x3:%f, x4:%f, x5:%f\n", x1, x2, x3, x4, x5); +// printf("y1:%f, y2:%f, y3:%f, y4:%f, y5:%f\n", y1, y2, y3, y4, y5); +// printf("z1:%f, z2:%f, z3:%f, z4:%f, z5:%f\n", z1, z2, z3, z4, z5); +// printf("rt4_ex:%f, rt5_ex:%f\n", expectrt4, expectrt5); +// printf("z4_ex:%f, z5_ex:%f\n", expectz4, expectz5); +// printf("Pt:%f, Px:%f, Py:%f, Pz:%f, charge: %i, residual4: %f, residual5:%f\n", Pt, Px, Py, Pz, charge, residual4, residual5); + printf("x1:%f, x2:%f, x3:%f, y1:%f, y2:%f, y3:%f, CenterX:%f, CenterY:%f, innerRadius:%f, Pt:%f, Px:%f, Py:%f, charge:%i\n", x1, x2, x3, y1, y2, y3, x_center, y_center, innerRadius, Pt, Px, Py, charge); + } //categories! if(layer1 == 1 and layer2 == 2 and layer3 == 3) diff --git a/SDL/Quintuplet.cuh b/SDL/Quintuplet.cuh index e4f7d1db..acd5f1d3 100644 --- a/SDL/Quintuplet.cuh +++ b/SDL/Quintuplet.cuh @@ -74,7 +74,7 @@ namespace SDL CUDA_DEV bool runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, struct SDL::segments& segmentsInGPU, struct SDL::triplets& tripletsInGPU, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, unsigned int& innerTripletIndex, unsigned int& outerTripletIndex, float& innerRadius, float& outerRadius, float& bridgeRadius, float& regressionG, float& regressionF, float& regressionRadius, float& rzChiSquared, float& chiSquared, float& nonAnchorChiSquared, float& residual_missing, float& residual4, float& residual5); - CUDA_DEV bool passT5RZConstraint(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex, unsigned int fifthMDIndex, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, float& rzChiSquared, float& residual_missing, float& residual4, float& residual5, float inner_pt, float g, float f); + CUDA_DEV bool passT5RZConstraint(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex, unsigned int fifthMDIndex, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, float& rzChiSquared, float& residual_missing, float& residual4, float& residual5, float inner_pt, float innerRadius, float g, float f); CUDA_DEV bool T5HasCommonMiniDoublet(struct triplets& tripletsInGPU, struct segments& segmentsInGPU, unsigned int innerTripletIndex, unsigned int outerTripletIndex); From 315134c491b9b76edd83ce9392a28bb9fe46d061 Mon Sep 17 00:00:00 2001 From: YonsiG Date: Mon, 5 Dec 2022 17:27:40 -0500 Subject: [PATCH 07/25] add second MD for chi2 --- SDL/Quintuplet.cu | 55 ++++++++++++++++++++++++++--------------------- 1 file changed, 31 insertions(+), 24 deletions(-) diff --git a/SDL/Quintuplet.cu b/SDL/Quintuplet.cu index ff228252..fbc2c931 100644 --- a/SDL/Quintuplet.cu +++ b/SDL/Quintuplet.cu @@ -681,14 +681,10 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc float pseudo_phi = atan((y_init-y_center)/(x_init-x_center)); //actually represent pi/2-phi, wrt helix axis z float Pt=inner_pt, Px, Py; - if (charge==1 && x_init<0 && y_init<0) {Px = -Pt*abs(sin(pseudo_phi)); Py=-Pt*abs(cos(pseudo_phi));} - else if (charge==1 && x_init<0 && y_init>0) {Px = -Pt*abs(sin(pseudo_phi)); Py=Pt*abs(cos(pseudo_phi));} - else if (charge==1 && x_init>0 && y_init>0) {Px = Pt*abs(sin(pseudo_phi)); Py=Pt*abs(cos(pseudo_phi));} - else if (charge==1 && x_init>0 && y_init<0) {Px = Pt*abs(sin(pseudo_phi)); Py=-Pt*abs(cos(pseudo_phi));} - else if (charge==-1 && x_init>0 && y_init>0) {Px = Pt*abs(sin(pseudo_phi)); Py=Pt*abs(cos(pseudo_phi));} - else if (charge==-1 && x_init>0 && y_init<0) {Px = Pt*abs(sin(pseudo_phi)); Py=-Pt*abs(cos(pseudo_phi));} - else if (charge==-1 && x_init<0 && y_init<0) {Px = -Pt*abs(sin(pseudo_phi)); Py=-Pt*abs(cos(pseudo_phi));} - else if (charge==-1 && x_init<0 && y_init>0) {Px = -Pt*abs(sin(pseudo_phi)); Py=Pt*abs(cos(pseudo_phi));} + if (x_init<0 && y_init<0) {Px = -Pt*abs(sin(pseudo_phi)); Py=-Pt*abs(cos(pseudo_phi));} + else if (x_init<0 && y_init>0) {Px = -Pt*abs(sin(pseudo_phi)); Py=Pt*abs(cos(pseudo_phi));} + else if (x_init>0 && y_init>0) {Px = Pt*abs(sin(pseudo_phi)); Py=Pt*abs(cos(pseudo_phi));} + else if (x_init>0 && y_init<0) {Px = Pt*abs(sin(pseudo_phi)); Py=-Pt*abs(cos(pseudo_phi));} else return 0; //to get Pz, we use pt/pz=ds/dz, ds is the arclength between MD1 and MD2. @@ -696,9 +692,9 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc float BO=sqrt((x_init-x_center)*(x_init-x_center)+(y_init-y_center)*(y_init-y_center)); float AB=sqrt((x1-x_init)*(x1-x_init)+(y1-y_init)*(y1-y_init)); float dPhi = acos((AO*AO+BO*BO-AB*AB)/(2*AO*BO)); -// float ds=innerRadius/100*dPhi; + float ds=innerRadius/100*dPhi; - float ds = sqrt((y_init-y1)*(y_init-y1)+(x_init-x1)*(x_init-x1)); //large ds->smallerPz->smaller residual +// float ds = sqrt((y_init-y1)*(y_init-y1)+(x_init-x1)*(x_init-x1)); //large ds->smallerPz->smaller residual float Pz=(z_init-z1)/ds*Pt; float p = sqrt(Px*Px+Py*Py+Pz*Pz); @@ -708,9 +704,16 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc float zsi, rtsi; int layeri, moduleTypei; float expectrt4=0,expectrt5=0,expectz4=0, expectz5=0; - for(size_t i = 4; i < 6; i++) + for(size_t i = 2; i < 6; i++) { - if (i==4){ + if (i==2){ + zsi = z2; + rtsi = rt2; + layeri=layer2; + moduleTypei=moduleType2; + } + else if (i==3) continue; + else if (i==4){ zsi = z4; rtsi = rt4; layeri=layer4; @@ -722,18 +725,6 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc layeri=layer5; moduleTypei=moduleType5; } - else if (i==2){ - zsi = z2; - rtsi = rt2; - layeri=layer2; - moduleTypei=moduleType2; - } - else if (i==3){ - zsi = z3; - rtsi = rt3; - layeri=layer3; - moduleTypei=moduleType3; - } // calculation is copied from PixelTriplet.cu SDL::computePT3RZChiSquared float diffr=0, diffz=0; @@ -781,6 +772,22 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc { error = 5.0f; } + + //check the tilted module, side: PosZ, NegZ, Center(for not tilted) + if (i==2){ + float drdz2=modulesInGPU.drdzs[lowerModuleIndex2]; + short side2=modulesInGPU.sides[lowerModuleIndex2]; + short subdets=modulesInGPU.subdets[lowerModuleIndex2]; + residual = (layeri <= 6 && ((side2 == SDL::Center) or (drdz2 < 1))) ? diffz : diffr; + residual_missing=residual; + float projection_missing=1; + if (drdz2<1) + projection_missing = ((subdets == SDL::Endcap) or (side2 == SDL::Center)) ? 1.f : cos(atan(drdz2)); //if dr/dz<1 + if (drdz2>1) + projection_missing = ((subdets == SDL::Endcap) or (side2 == SDL::Center)) ? 1.f : sin(atan(drdz2)); //if dr/dz>1 + error=error*projection_missing; + } + rzChiSquared += 12*(residual * residual)/(error * error); } From e41642a25ef14cc3309e00945f900c9107eecf2d Mon Sep 17 00:00:00 2001 From: YonsiG Date: Tue, 6 Dec 2022 17:53:46 -0500 Subject: [PATCH 08/25] add T5 rz cut numbers --- SDL/Quintuplet.cu | 123 ++++++++++++++++++++++++---------------------- 1 file changed, 65 insertions(+), 58 deletions(-) diff --git a/SDL/Quintuplet.cu b/SDL/Quintuplet.cu index fbc2c931..e4c178e4 100644 --- a/SDL/Quintuplet.cu +++ b/SDL/Quintuplet.cu @@ -399,7 +399,7 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, // pass = pass and passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f); float inner_pt = 2 * k2Rinv1GeVf * innerRadius; - passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f); + pass = pass and passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f); if(not pass) return pass; @@ -447,7 +447,7 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, } //compute regression radius right here - this computation is expensive!!! -// pass = pass and tempPass; + pass = pass and tempPass; if(not pass) return pass; float xVec[] = {x1, x2, x3, x4, x5}; @@ -463,7 +463,7 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, //extra chi squared cuts! if(regressionRadius < 5.0f/(2.f * k2Rinv1GeVf)) { -// pass = pass and passChiSquaredConstraint(modulesInGPU, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, chiSquared); + pass = pass and passChiSquaredConstraint(modulesInGPU, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, chiSquared); if(not pass) return pass; } @@ -660,8 +660,8 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc //start from a circle of inner T3. // to determine the charge int charge=0; - float slope1c=(y3-y_center)/(x3-x_center); - float slope3c=(y1-y_center)/(x1-x_center); + float slope3c=(y3-y_center)/(x3-x_center); + float slope1c=(y1-y_center)/(x1-x_center); if((y3-y_center)>0 && (y1-y_center)>0) { if (slope3c>slope1c) charge=-1; @@ -675,9 +675,9 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc else if(slope3c<-5 && slope1c>5) charge=-1; else if(slope3c>5 && slope1c<-5) charge=1; - if (((y3-y2)/(x3-x2))<((y2-y1)/(x2-x1))) charge=1; - else if (((y3-y2)/(x3-x2))>((y2-y1)/(x2-x1))) charge=-1; - else return 0; +// if (((y3-y2)/(x3-x2))<((y2-y1)/(x2-x1))) charge=1; +// else if (((y3-y2)/(x3-x2))>((y2-y1)/(x2-x1))) charge=-1; + // else return 0; float pseudo_phi = atan((y_init-y_center)/(x_init-x_center)); //actually represent pi/2-phi, wrt helix axis z float Pt=inner_pt, Px, Py; @@ -704,6 +704,7 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc float zsi, rtsi; int layeri, moduleTypei; float expectrt4=0,expectrt5=0,expectz4=0, expectz5=0; + rzChiSquared=0; for(size_t i = 2; i < 6; i++) { if (i==2){ @@ -760,8 +761,6 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc } residual = (layeri>6) ? diffr : diffz ; - if (i==4) residual4=residual; - if (i==5) residual5=residual; //PS Modules if(moduleTypei == 0) @@ -772,26 +771,37 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc { error = 5.0f; } + if (i==4) residual4=residual/error; + if (i==5) residual5=residual/error; //check the tilted module, side: PosZ, NegZ, Center(for not tilted) + float drdz; + short side, subdets; if (i==2){ - float drdz2=modulesInGPU.drdzs[lowerModuleIndex2]; - short side2=modulesInGPU.sides[lowerModuleIndex2]; - short subdets=modulesInGPU.subdets[lowerModuleIndex2]; - residual = (layeri <= 6 && ((side2 == SDL::Center) or (drdz2 < 1))) ? diffz : diffr; - residual_missing=residual; + drdz=modulesInGPU.drdzs[lowerModuleIndex2]; + side=modulesInGPU.sides[lowerModuleIndex2]; + subdets=modulesInGPU.subdets[lowerModuleIndex2]; + } + if (i==3){ + drdz=modulesInGPU.drdzs[lowerModuleIndex3]; + side=modulesInGPU.sides[lowerModuleIndex3]; + subdets=modulesInGPU.subdets[lowerModuleIndex3]; + } + if (i==2 || i==3){ + residual = (layeri <= 6 && ((side == SDL::Center) or (drdz < 1))) ? diffz : diffr; +// residual_missing=residual; float projection_missing=1; - if (drdz2<1) - projection_missing = ((subdets == SDL::Endcap) or (side2 == SDL::Center)) ? 1.f : cos(atan(drdz2)); //if dr/dz<1 - if (drdz2>1) - projection_missing = ((subdets == SDL::Endcap) or (side2 == SDL::Center)) ? 1.f : sin(atan(drdz2)); //if dr/dz>1 + if (drdz<1) + projection_missing = ((subdets == SDL::Endcap) or (side == SDL::Center)) ? 1.f : 1/sqrt(1+drdz*drdz); // cos(atan(drdz)), if dr/dz<1 + if (drdz>1) + projection_missing = ((subdets == SDL::Endcap) or (side == SDL::Center)) ? 1.f : drdz/sqrt(1+drdz*drdz);//sin(atan(drdz)), if dr/dz>1 error=error*projection_missing; + residual_missing=residual/error; } - rzChiSquared += 12*(residual * residual)/(error * error); } - if(layer1 == 1 and layer2 == 2 and layer3 == 3 and layer4 == 4 and layer5 == 5 and residual4>15){ +// if(layer1 == 1 and layer2 == 2 and layer3 == 3 and layer4 == 4 and layer5 == 5 and rzChiSquared>1000){ // printf("rt1:%f, rt2:%f, rt3:%f, rt4:%f, rt5:%f\n", rt1, rt2, rt3, rt4, rt5); // printf("x1:%f, x2:%f, x3:%f, x4:%f, x5:%f\n", x1, x2, x3, x4, x5); // printf("y1:%f, y2:%f, y3:%f, y4:%f, y5:%f\n", y1, y2, y3, y4, y5); @@ -799,120 +809,117 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc // printf("rt4_ex:%f, rt5_ex:%f\n", expectrt4, expectrt5); // printf("z4_ex:%f, z5_ex:%f\n", expectz4, expectz5); // printf("Pt:%f, Px:%f, Py:%f, Pz:%f, charge: %i, residual4: %f, residual5:%f\n", Pt, Px, Py, Pz, charge, residual4, residual5); - printf("x1:%f, x2:%f, x3:%f, y1:%f, y2:%f, y3:%f, CenterX:%f, CenterY:%f, innerRadius:%f, Pt:%f, Px:%f, Py:%f, charge:%i\n", x1, x2, x3, y1, y2, y3, x_center, y_center, innerRadius, Pt, Px, Py, charge); - } +// printf("x1:%f, x2:%f, x3:%f, y1:%f, y2:%f, y3:%f, CenterX:%f, CenterY:%f, innerRadius:%f, Pt:%f, Px:%f, Py:%f, charge:%i\n", x1, x2, x3, y1, y2, y3, x_center, y_center, innerRadius, Pt, Px, Py, charge); +// printf("rzChi2: %f, residual2:%f, residual4: %f, residual5: %f\n", rzChiSquared, residual_missing, residual4, residual5); +// } //categories! if(layer1 == 1 and layer2 == 2 and layer3 == 3) { if(layer4 == 4 and layer5 == 5) { - return rzChiSquared < 0.545f; + return rzChiSquared < 136.5975f; } else if(layer4 == 4 and layer5 == 12) { - return rzChiSquared < 1.105f; + return rzChiSquared < 13.8784f; } else if(layer4 == 7 and layer5 == 13) { - return rzChiSquared < 0.775f; + return rzChiSquared < 197.4f; } else if(layer4 == 12 and layer5 == 13) { - return rzChiSquared < 0.625f; + return rzChiSquared < 9.9563f; } } else if(layer1 == 1 and layer2 == 2 and layer3 == 7) { if(layer4 == 8 and layer5 == 14) { - return rzChiSquared < 0.835f; + return rzChiSquared < 168.5f; } else if(layer4 == 13 and layer5 == 14) { - return rzChiSquared < 0.575f; + return rzChiSquared < 9.9181f; } } else if(layer1 == 1 and layer2 == 7 and layer3 == 8 and layer4 == 9 and layer5 == 15) { - return rzChiSquared < 0.825f; + return rzChiSquared < 99.5923f; } else if(layer1 == 2 and layer2 == 3 and layer3 == 4) { if(layer4 == 5 and layer5 == 6) { - return rzChiSquared < 0.845f; + return rzChiSquared < 933.1022f; } else if(layer4 == 5 and layer5 == 12) { - return rzChiSquared < 1.365f; + return rzChiSquared < 952.8f; } else if(layer4 == 12 and layer5 == 13) { - return rzChiSquared < 0.675f; + return rzChiSquared < 222.0667f; } } else if(layer1 == 2 and layer2 == 3 and layer3 == 7 and layer4 == 13 and layer5 == 14) { - return rzChiSquared < 0.495f; + return rzChiSquared < 9.9041f; } else if(layer1 == 2 and layer2 == 3 and layer3 == 12 and layer4 == 13 and layer5 == 14) { - return rzChiSquared < 0.695f; + return 1; +// return rzChiSquared < 0.695f; } else if(layer1 == 2 and layer2 == 7 and layer3 == 8) { if(layer4 == 9 and layer5 == 15) { - return rzChiSquared < 0.735f; + return rzChiSquared < 259.7f; } else if(layer4 == 14 and layer5 == 15) { - return rzChiSquared < 0.525f; + return rzChiSquared < 128.8375f; } else if(layer4 == 9 and layer5 == 10) { printf("here"); + return 1; } } - else if(layer1 == 2 and layer2 == 7 and layer3 == 13 and layer4 == 14 and layer5 == 15) - { - return rzChiSquared < 0.665f; - } + else if(layer1 == 3 and layer2 == 4 and layer3 == 5 and layer4 == 12 and layer5 == 13) { - return rzChiSquared < 0.995f; + return 1; +// return rzChiSquared < 0.995f; } else if(layer1 == 3 and layer2 == 4 and layer3 == 12 and layer4 == 13 and layer5 == 14) { - return rzChiSquared < 0.525f; + return 1; +// return rzChiSquared < 0.525f; } else if(layer1 == 3 and layer2 == 7 and layer3 == 8 and layer4 == 14 and layer5 == 15) { - return rzChiSquared < 0.525f; - } - else if(layer1 == 3 and layer2 == 7 and layer3 == 13 and layer4 == 14 and layer5 == 15) - { - return rzChiSquared < 0.745f; - } - else if(layer1 == 3 and layer2 == 12 and layer3 == 13 and layer4 == 14 and layer5 == 15) - { - return rzChiSquared < 0.555f; + return 1; +// return rzChiSquared < 0.525f; } + else if(layer1 == 7 and layer2 == 8 and layer3 == 9 and layer4 == 15 and layer5 == 16) { - return rzChiSquared < 0.525f; + return rzChiSquared < 157.05f; } - else if(layer1 == 7 and layer2 == 8 and layer3 == 14 and layer4 == 15 and layer5 == 16) + + else if(layer1 == 7 and layer2 == 8 and layer3 == 9 and layer4 == 10 and layer5 == 11) { - return rzChiSquared < 0.885f; + return rzChiSquared < 360.1125f; } - else if(layer1 == 7 and layer2 == 13 and layer3 == 14 and layer4 == 15 and layer5 == 16) + + else if(layer1 == 7 and layer2 == 8 and layer3 == 9 and layer4 == 10 and layer5 == 16) { - return rzChiSquared < 0.845f; + return rzChiSquared < 157.2273f; } - return true; } From 7a342ed0155384d921eff8a80321828fe8d8681b Mon Sep 17 00:00:00 2001 From: YonsiG Date: Wed, 7 Dec 2022 08:18:30 -0500 Subject: [PATCH 09/25] remove a few categories --- SDL/Quintuplet.cu | 12 ++---------- 1 file changed, 2 insertions(+), 10 deletions(-) diff --git a/SDL/Quintuplet.cu b/SDL/Quintuplet.cu index e4c178e4..e891eb45 100644 --- a/SDL/Quintuplet.cu +++ b/SDL/Quintuplet.cu @@ -868,11 +868,7 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc { return rzChiSquared < 9.9041f; } - else if(layer1 == 2 and layer2 == 3 and layer3 == 12 and layer4 == 13 and layer5 == 14) - { - return 1; -// return rzChiSquared < 0.695f; - } + else if(layer1 == 2 and layer2 == 7 and layer3 == 8) { if(layer4 == 9 and layer5 == 15) @@ -895,11 +891,7 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc return 1; // return rzChiSquared < 0.995f; } - else if(layer1 == 3 and layer2 == 4 and layer3 == 12 and layer4 == 13 and layer5 == 14) - { - return 1; -// return rzChiSquared < 0.525f; - } + else if(layer1 == 3 and layer2 == 7 and layer3 == 8 and layer4 == 14 and layer5 == 15) { return 1; From 17c125fffafb508415ac4258dd5896beeab43562 Mon Sep 17 00:00:00 2001 From: YonsiG Date: Mon, 19 Dec 2022 11:58:14 -0500 Subject: [PATCH 10/25] add charge correction --- SDL/Quintuplet.cu | 34 ++++++++++++++++------------------ 1 file changed, 16 insertions(+), 18 deletions(-) diff --git a/SDL/Quintuplet.cu b/SDL/Quintuplet.cu index e891eb45..fec19e9d 100644 --- a/SDL/Quintuplet.cu +++ b/SDL/Quintuplet.cu @@ -399,9 +399,9 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, // pass = pass and passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f); float inner_pt = 2 * k2Rinv1GeVf * innerRadius; - pass = pass and passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f); + passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f); - if(not pass) return pass; +// if(not pass) return pass; pass = pass & (innerRadius >= 0.95f * ptCut/(2.f * k2Rinv1GeVf)); @@ -447,8 +447,8 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, } //compute regression radius right here - this computation is expensive!!! - pass = pass and tempPass; - if(not pass) return pass; +// pass = pass and tempPass; +// if(not pass) return pass; float xVec[] = {x1, x2, x3, x4, x5}; float yVec[] = {y1, y2, y3, y4, y5}; @@ -463,8 +463,8 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, //extra chi squared cuts! if(regressionRadius < 5.0f/(2.f * k2Rinv1GeVf)) { - pass = pass and passChiSquaredConstraint(modulesInGPU, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, chiSquared); - if(not pass) return pass; +// pass = pass and passChiSquaredConstraint(modulesInGPU, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, chiSquared); +// if(not pass) return pass; } //compute the other chisquared @@ -666,18 +666,14 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc { if (slope3c>slope1c) charge=-1; else if (slope3c0 && slope3c<0) charge=-1; } else if((y3-y_center)<0 && (y1-y_center)<0) { if (slope3c>slope1c) charge=-1; else if (slope3c0 && slope3c<0) charge=-1; } - else if(slope3c<-5 && slope1c>5) charge=-1; - else if(slope3c>5 && slope1c<-5) charge=1; - -// if (((y3-y2)/(x3-x2))<((y2-y1)/(x2-x1))) charge=1; -// else if (((y3-y2)/(x3-x2))>((y2-y1)/(x2-x1))) charge=-1; - // else return 0; float pseudo_phi = atan((y_init-y_center)/(x_init-x_center)); //actually represent pi/2-phi, wrt helix axis z float Pt=inner_pt, Px, Py; @@ -778,12 +774,13 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc float drdz; short side, subdets; if (i==2){ - drdz=modulesInGPU.drdzs[lowerModuleIndex2]; + drdz=abs(modulesInGPU.drdzs[lowerModuleIndex2]); side=modulesInGPU.sides[lowerModuleIndex2]; subdets=modulesInGPU.subdets[lowerModuleIndex2]; } if (i==3){ - drdz=modulesInGPU.drdzs[lowerModuleIndex3]; + + drdz=abs(modulesInGPU.drdzs[lowerModuleIndex3]); side=modulesInGPU.sides[lowerModuleIndex3]; subdets=modulesInGPU.subdets[lowerModuleIndex3]; } @@ -800,8 +797,9 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc } rzChiSquared += 12*(residual * residual)/(error * error); } +// rzChiSquared = 12*(residual4 * residual4 + residual5 * residual5 + residual_missing * residual_missing); -// if(layer1 == 1 and layer2 == 2 and layer3 == 3 and layer4 == 4 and layer5 == 5 and rzChiSquared>1000){ + if(layer1 == 1 and layer2 == 2 and layer3 == 3 and layer4 == 4 and layer5 == 5 and rzChiSquared>1000){ // printf("rt1:%f, rt2:%f, rt3:%f, rt4:%f, rt5:%f\n", rt1, rt2, rt3, rt4, rt5); // printf("x1:%f, x2:%f, x3:%f, x4:%f, x5:%f\n", x1, x2, x3, x4, x5); // printf("y1:%f, y2:%f, y3:%f, y4:%f, y5:%f\n", y1, y2, y3, y4, y5); @@ -810,15 +808,15 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc // printf("z4_ex:%f, z5_ex:%f\n", expectz4, expectz5); // printf("Pt:%f, Px:%f, Py:%f, Pz:%f, charge: %i, residual4: %f, residual5:%f\n", Pt, Px, Py, Pz, charge, residual4, residual5); // printf("x1:%f, x2:%f, x3:%f, y1:%f, y2:%f, y3:%f, CenterX:%f, CenterY:%f, innerRadius:%f, Pt:%f, Px:%f, Py:%f, charge:%i\n", x1, x2, x3, y1, y2, y3, x_center, y_center, innerRadius, Pt, Px, Py, charge); -// printf("rzChi2: %f, residual2:%f, residual4: %f, residual5: %f\n", rzChiSquared, residual_missing, residual4, residual5); -// } + printf("rzChi2: %f, residual2: %f, inner_pt:%f, pseudo_phi: %f, charge: %i, Px:%f, Py:%f, x1:%f, x2:%f, x3:%f, x4:%f, x5:%f, y1:%f, y2:%f, y3:%f, y4:%f, y5:%f, z1:%f, z2:%f, z3:%f, z4:%f, z5:%f, x_center:%f, y_center:%f, slope1c:%f, slope3c:%f\n", rzChiSquared, residual_missing, inner_pt, pseudo_phi, charge, Px, Py, x1, x2, x3, x4, x5, y1, y2, y3, y4, y5, z1, z2, z3, z4, z5, x_center, y_center, slope1c, slope3c); + } //categories! if(layer1 == 1 and layer2 == 2 and layer3 == 3) { if(layer4 == 4 and layer5 == 5) { - return rzChiSquared < 136.5975f; + return rzChiSquared < 20.5139f; } else if(layer4 == 4 and layer5 == 12) { From 9c5c4e590c18c19aa67f53881d7c3d1a5eef74d6 Mon Sep 17 00:00:00 2001 From: YonsiG Date: Tue, 27 Dec 2022 00:55:32 -0500 Subject: [PATCH 11/25] add charge correction for nan --- SDL/Quintuplet.cu | 39 ++++++++++++++++++++++++++++++++++++--- 1 file changed, 36 insertions(+), 3 deletions(-) diff --git a/SDL/Quintuplet.cu b/SDL/Quintuplet.cu index fec19e9d..3cb56461 100644 --- a/SDL/Quintuplet.cu +++ b/SDL/Quintuplet.cu @@ -657,6 +657,14 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc float z_init=mdsInGPU.anchorZ[thirdMDIndex]/100; float rt_init=mdsInGPU.anchorRt[thirdMDIndex]/100; //use the second MD as initial point + if (moduleType3==1) // 1: if MD3 is in 2s layer + { + x_init=mdsInGPU.anchorX[secondMDIndex]/100; + y_init=mdsInGPU.anchorY[secondMDIndex]/100; + z_init=mdsInGPU.anchorZ[secondMDIndex]/100; + rt_init=mdsInGPU.anchorRt[secondMDIndex]/100; + } + //start from a circle of inner T3. // to determine the charge int charge=0; @@ -667,13 +675,25 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc if (slope3c>slope1c) charge=-1; else if (slope3c0 && slope3c<0) charge=-1; + if (slope1c<0 && slope3c>0) charge=1; } else if((y3-y_center)<0 && (y1-y_center)<0) { if (slope3c>slope1c) charge=-1; else if (slope3c0) charge=1; if (slope1c>0 && slope3c<0) charge=-1; } + else if ((y3-y_center)<0 && (y1-y_center)>0) + { + if ((x3-x_center)>0 && (x1-x_center)>0) charge = 1; + else if ((x3-x_center)<0 && (x1-x_center)<0) charge = -1; + } + else if ((y3-y_center)>0 && (y1-y_center)<0) + { + if ((x3-x_center)>0 && (x1-x_center)>0) charge = -1; + else if ((x3-x_center)<0 && (x1-x_center)<0) charge = 1; + } float pseudo_phi = atan((y_init-y_center)/(x_init-x_center)); //actually represent pi/2-phi, wrt helix axis z float Pt=inner_pt, Px, Py; @@ -683,7 +703,7 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc else if (x_init>0 && y_init<0) {Px = Pt*abs(sin(pseudo_phi)); Py=-Pt*abs(cos(pseudo_phi));} else return 0; - //to get Pz, we use pt/pz=ds/dz, ds is the arclength between MD1 and MD2. + //to get Pz, we use pt/pz=ds/dz, ds is the arclength between MD1 and MD3. float AO=sqrt((x1-x_center)*(x1-x_center)+(y1-y_center)*(y1-y_center)); float BO=sqrt((x_init-x_center)*(x_init-x_center)+(y_init-y_center)*(y_init-y_center)); float AB=sqrt((x1-x_init)*(x1-x_init)+(y1-y_init)*(y1-y_init)); @@ -709,7 +729,12 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc layeri=layer2; moduleTypei=moduleType2; } - else if (i==3) continue; + else if (i==3) { + zsi = z3; + rtsi = rt3; + layeri=layer3; + moduleTypei=moduleType3; + } else if (i==4){ zsi = z4; rtsi = rt4; @@ -723,6 +748,13 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc moduleTypei=moduleType5; } + if (moduleType3==0) { //0: ps + if (i==3) continue; + } + else{ + if (i==2) continue; + } + // calculation is copied from PixelTriplet.cu SDL::computePT3RZChiSquared float diffr=0, diffz=0; @@ -799,6 +831,8 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc } // rzChiSquared = 12*(residual4 * residual4 + residual5 * residual5 + residual_missing * residual_missing); +// if (isnan(rzChiSquared)) printf("rzChi2: %f, residual2: %f, inner_pt:%f, pseudo_phi: %f, charge: %i, Px:%f, Py:%f, x1:%f, x2:%f, x3:%f, x4:%f, x5:%f, y1:%f, y2:%f, y3:%f, y4:%f, y5:%f, z1:%f, z2:%f, z3:%f, z4:%f, z5:%f, x_center:%f, y_center:%f, slope1c:%f, slope3c:%f\n", rzChiSquared, residual_missing, inner_pt, pseudo_phi, charge, Px, Py, x1, x2, x3, x4, x5, y1, y2, y3, y4, y5, z1, z2, z3, z4, z5, x_center, y_center, slope1c, slope3c); + if(layer1 == 1 and layer2 == 2 and layer3 == 3 and layer4 == 4 and layer5 == 5 and rzChiSquared>1000){ // printf("rt1:%f, rt2:%f, rt3:%f, rt4:%f, rt5:%f\n", rt1, rt2, rt3, rt4, rt5); // printf("x1:%f, x2:%f, x3:%f, x4:%f, x5:%f\n", x1, x2, x3, x4, x5); @@ -807,7 +841,6 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc // printf("rt4_ex:%f, rt5_ex:%f\n", expectrt4, expectrt5); // printf("z4_ex:%f, z5_ex:%f\n", expectz4, expectz5); // printf("Pt:%f, Px:%f, Py:%f, Pz:%f, charge: %i, residual4: %f, residual5:%f\n", Pt, Px, Py, Pz, charge, residual4, residual5); -// printf("x1:%f, x2:%f, x3:%f, y1:%f, y2:%f, y3:%f, CenterX:%f, CenterY:%f, innerRadius:%f, Pt:%f, Px:%f, Py:%f, charge:%i\n", x1, x2, x3, y1, y2, y3, x_center, y_center, innerRadius, Pt, Px, Py, charge); printf("rzChi2: %f, residual2: %f, inner_pt:%f, pseudo_phi: %f, charge: %i, Px:%f, Py:%f, x1:%f, x2:%f, x3:%f, x4:%f, x5:%f, y1:%f, y2:%f, y3:%f, y4:%f, y5:%f, z1:%f, z2:%f, z3:%f, z4:%f, z5:%f, x_center:%f, y_center:%f, slope1c:%f, slope3c:%f\n", rzChiSquared, residual_missing, inner_pt, pseudo_phi, charge, Px, Py, x1, x2, x3, x4, x5, y1, y2, y3, y4, y5, z1, z2, z3, z4, z5, x_center, y_center, slope1c, slope3c); } From e2205a25b77224ec17b6870498fec2e004aeb26f Mon Sep 17 00:00:00 2001 From: YonsiG Date: Tue, 27 Dec 2022 20:39:11 -0500 Subject: [PATCH 12/25] update cut values --- SDL/Quintuplet.cu | 122 ++++++++++++++++++++++++---------------------- 1 file changed, 65 insertions(+), 57 deletions(-) diff --git a/SDL/Quintuplet.cu b/SDL/Quintuplet.cu index 3cb56461..a6f7b8cc 100644 --- a/SDL/Quintuplet.cu +++ b/SDL/Quintuplet.cu @@ -696,12 +696,15 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc } float pseudo_phi = atan((y_init-y_center)/(x_init-x_center)); //actually represent pi/2-phi, wrt helix axis z - float Pt=inner_pt, Px, Py; - if (x_init<0 && y_init<0) {Px = -Pt*abs(sin(pseudo_phi)); Py=-Pt*abs(cos(pseudo_phi));} - else if (x_init<0 && y_init>0) {Px = -Pt*abs(sin(pseudo_phi)); Py=Pt*abs(cos(pseudo_phi));} - else if (x_init>0 && y_init>0) {Px = Pt*abs(sin(pseudo_phi)); Py=Pt*abs(cos(pseudo_phi));} - else if (x_init>0 && y_init<0) {Px = Pt*abs(sin(pseudo_phi)); Py=-Pt*abs(cos(pseudo_phi));} - else return 0; + float Pt=inner_pt, Px=Pt*abs(sin(pseudo_phi)), Py=Pt*abs(cos(pseudo_phi)); +// if (x_init<0 && y_init<0) {Px = -Pt*abs(sin(pseudo_phi)); Py=-Pt*abs(cos(pseudo_phi));} +// else if (x_init<0 && y_init>0) {Px = -Pt*abs(sin(pseudo_phi)); Py=Pt*abs(cos(pseudo_phi));} +// else if (x_init>0 && y_init>0) {Px = Pt*abs(sin(pseudo_phi)); Py=Pt*abs(cos(pseudo_phi));} +// else if (x_init>0 && y_init<0) {Px = Pt*abs(sin(pseudo_phi)); Py=-Pt*abs(cos(pseudo_phi));} + if (x4-x2<0) Px = -Px; + if (y4-y2<0) Py = -Py; +// else return 0; + //to get Pz, we use pt/pz=ds/dz, ds is the arclength between MD1 and MD3. float AO=sqrt((x1-x_center)*(x1-x_center)+(y1-y_center)*(y1-y_center)); @@ -831,9 +834,9 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc } // rzChiSquared = 12*(residual4 * residual4 + residual5 * residual5 + residual_missing * residual_missing); -// if (isnan(rzChiSquared)) printf("rzChi2: %f, residual2: %f, inner_pt:%f, pseudo_phi: %f, charge: %i, Px:%f, Py:%f, x1:%f, x2:%f, x3:%f, x4:%f, x5:%f, y1:%f, y2:%f, y3:%f, y4:%f, y5:%f, z1:%f, z2:%f, z3:%f, z4:%f, z5:%f, x_center:%f, y_center:%f, slope1c:%f, slope3c:%f\n", rzChiSquared, residual_missing, inner_pt, pseudo_phi, charge, Px, Py, x1, x2, x3, x4, x5, y1, y2, y3, y4, y5, z1, z2, z3, z4, z5, x_center, y_center, slope1c, slope3c); + if (isnan(rzChiSquared)) printf("rzChi2: %f, residual2: %f, inner_pt:%f, pseudo_phi: %f, charge: %i, Px:%f, Py:%f, x1:%f, x2:%f, x3:%f, x4:%f, x5:%f, y1:%f, y2:%f, y3:%f, y4:%f, y5:%f, z1:%f, z2:%f, z3:%f, z4:%f, z5:%f, x_center:%f, y_center:%f, slope1c:%f, slope3c:%f\n", rzChiSquared, residual_missing, inner_pt, pseudo_phi, charge, Px, Py, x1, x2, x3, x4, x5, y1, y2, y3, y4, y5, z1, z2, z3, z4, z5, x_center, y_center, slope1c, slope3c); - if(layer1 == 1 and layer2 == 2 and layer3 == 3 and layer4 == 4 and layer5 == 5 and rzChiSquared>1000){ + if(layer1 == 1 and layer2 == 2 and layer3 == 3 and layer4 == 4 and layer5 == 5 and rzChiSquared>100){ // printf("rt1:%f, rt2:%f, rt3:%f, rt4:%f, rt5:%f\n", rt1, rt2, rt3, rt4, rt5); // printf("x1:%f, x2:%f, x3:%f, x4:%f, x5:%f\n", x1, x2, x3, x4, x5); // printf("y1:%f, y2:%f, y3:%f, y4:%f, y5:%f\n", y1, y2, y3, y4, y5); @@ -841,7 +844,7 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc // printf("rt4_ex:%f, rt5_ex:%f\n", expectrt4, expectrt5); // printf("z4_ex:%f, z5_ex:%f\n", expectz4, expectz5); // printf("Pt:%f, Px:%f, Py:%f, Pz:%f, charge: %i, residual4: %f, residual5:%f\n", Pt, Px, Py, Pz, charge, residual4, residual5); - printf("rzChi2: %f, residual2: %f, inner_pt:%f, pseudo_phi: %f, charge: %i, Px:%f, Py:%f, x1:%f, x2:%f, x3:%f, x4:%f, x5:%f, y1:%f, y2:%f, y3:%f, y4:%f, y5:%f, z1:%f, z2:%f, z3:%f, z4:%f, z5:%f, x_center:%f, y_center:%f, slope1c:%f, slope3c:%f\n", rzChiSquared, residual_missing, inner_pt, pseudo_phi, charge, Px, Py, x1, x2, x3, x4, x5, y1, y2, y3, y4, y5, z1, z2, z3, z4, z5, x_center, y_center, slope1c, slope3c); +// printf("rzChi2: %f, residual2: %f, inner_pt:%f, pseudo_phi: %f, charge: %i, Px:%f, Py:%f, x1:%f, x2:%f, x3:%f, x4:%f, x5:%f, y1:%f, y2:%f, y3:%f, y4:%f, y5:%f, z1:%f, z2:%f, z3:%f, z4:%f, z5:%f, x_center:%f, y_center:%f, slope1c:%f, slope3c:%f\n", rzChiSquared, residual_missing, inner_pt, pseudo_phi, charge, Px, Py, x1, x2, x3, x4, x5, y1, y2, y3, y4, y5, z1, z2, z3, z4, z5, x_center, y_center, slope1c, slope3c); } //categories! @@ -849,99 +852,104 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc { if(layer4 == 4 and layer5 == 5) { - return rzChiSquared < 20.5139f; + return rzChiSquared < 28.3f; //11 } - else if(layer4 == 4 and layer5 == 12) + else if(layer4 == 4 and layer5 == 12) //12 { - return rzChiSquared < 13.8784f; + return rzChiSquared < 4.475f; } - else if(layer4 == 7 and layer5 == 13) + else if(layer4 == 7 and layer5 == 8) //8 + { + return rzChiSquared < 44.243f; + } + else if(layer4 == 7 and layer5 == 13) //9 { - return rzChiSquared < 197.4f; + return rzChiSquared < 30.405f; } - else if(layer4 == 12 and layer5 == 13) + else if(layer4 == 12 and layer5 == 13) //10 { - return rzChiSquared < 9.9563f; + return rzChiSquared < 4.331f; } } else if(layer1 == 1 and layer2 == 2 and layer3 == 7) { - if(layer4 == 8 and layer5 == 14) + if(layer4 == 8 and layer5 == 9) //5 { - return rzChiSquared < 168.5f; + return rzChiSquared < 113.184f; } - else if(layer4 == 13 and layer5 == 14) + if(layer4 == 8 and layer5 == 14) //6 + { + return rzChiSquared < 38.839f; + } + else if(layer4 == 13 and layer5 == 14) //7 { - return rzChiSquared < 9.9181f; + return rzChiSquared < 4.868f; } } - else if(layer1 == 1 and layer2 == 7 and layer3 == 8 and layer4 == 9 and layer5 == 15) + else if(layer1 == 1 and layer2 == 7 and layer3 == 8 and layer4 == 9) { - return rzChiSquared < 99.5923f; + if (layer5 == 10) //3 + { + return rzChiSquared < 109.577f; + } + if (layer5 == 15) //4 + { + return rzChiSquared < 37.933f; + } } else if(layer1 == 2 and layer2 == 3 and layer3 == 4) { - if(layer4 == 5 and layer5 == 6) + if(layer4 == 5 and layer5 == 6) //18 { - return rzChiSquared < 933.1022f; + return rzChiSquared < 7.951f; } - else if(layer4 == 5 and layer5 == 12) + else if(layer4 == 5 and layer5 == 12) //19 { - return rzChiSquared < 952.8f; + return rzChiSquared < 14.413f; } - else if(layer4 == 12 and layer5 == 13) + else if(layer4 == 12 and layer5 == 13) //20 { - return rzChiSquared < 222.0667f; + return rzChiSquared < 12.310f; } } - else if(layer1 == 2 and layer2 == 3 and layer3 == 7 and layer4 == 13 and layer5 == 14) + else if(layer1 == 2 and layer2 == 3 and layer3 == 7) { - return rzChiSquared < 9.9041f; + if(layer4 == 8 and layer5 == 14) //16 + { + return rzChiSquared < 18.649f; + } + if(layer4 == 13 and layer5 == 14) //17 + { + return rzChiSquared < 4.661f; + } } else if(layer1 == 2 and layer2 == 7 and layer3 == 8) { - if(layer4 == 9 and layer5 == 15) - { - return rzChiSquared < 259.7f; - } - else if(layer4 == 14 and layer5 == 15) + if(layer4 == 9 and layer5 == 15) //14 { - return rzChiSquared < 128.8375f; + return rzChiSquared < 41.036f; } - else if(layer4 == 9 and layer5 == 10) + else if(layer4 == 14 and layer5 == 15) //15 { - printf("here"); - return 1; + return rzChiSquared < 13.821f; } } - else if(layer1 == 3 and layer2 == 4 and layer3 == 5 and layer4 == 12 and layer5 == 13) - { - return 1; -// return rzChiSquared < 0.995f; - } - - else if(layer1 == 3 and layer2 == 7 and layer3 == 8 and layer4 == 14 and layer5 == 15) - { - return 1; -// return rzChiSquared < 0.525f; - } - - else if(layer1 == 7 and layer2 == 8 and layer3 == 9 and layer4 == 15 and layer5 == 16) + else if(layer1 == 7 and layer2 == 8 and layer3 == 9 and layer4 == 15 and layer5 == 16) //2 { - return rzChiSquared < 157.05f; + return rzChiSquared < 12.223f; } - else if(layer1 == 7 and layer2 == 8 and layer3 == 9 and layer4 == 10 and layer5 == 11) + else if(layer1 == 7 and layer2 == 8 and layer3 == 9 and layer4 == 10 and layer5 == 11) //0 { - return rzChiSquared < 360.1125f; + return rzChiSquared < 93.893f; } - else if(layer1 == 7 and layer2 == 8 and layer3 == 9 and layer4 == 10 and layer5 == 16) + else if(layer1 == 7 and layer2 == 8 and layer3 == 9 and layer4 == 10 and layer5 == 16) //1 { - return rzChiSquared < 157.2273f; + return rzChiSquared < 37.087f; } return true; } From f44a7fec6c440f9c470f1783a6fc061db14dd6e3 Mon Sep 17 00:00:00 2001 From: YonsiG Date: Wed, 28 Dec 2022 01:20:55 -0500 Subject: [PATCH 13/25] add back pass cuts (deleted for testing) --- SDL/Quintuplet.cu | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) diff --git a/SDL/Quintuplet.cu b/SDL/Quintuplet.cu index a6f7b8cc..d9b097ea 100644 --- a/SDL/Quintuplet.cu +++ b/SDL/Quintuplet.cu @@ -399,9 +399,9 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, // pass = pass and passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f); float inner_pt = 2 * k2Rinv1GeVf * innerRadius; - passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f); + pass = pass and passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f); -// if(not pass) return pass; + if(not pass) return pass; pass = pass & (innerRadius >= 0.95f * ptCut/(2.f * k2Rinv1GeVf)); @@ -447,8 +447,8 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, } //compute regression radius right here - this computation is expensive!!! -// pass = pass and tempPass; -// if(not pass) return pass; + pass = pass and tempPass; + if(not pass) return pass; float xVec[] = {x1, x2, x3, x4, x5}; float yVec[] = {y1, y2, y3, y4, y5}; @@ -463,8 +463,8 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, //extra chi squared cuts! if(regressionRadius < 5.0f/(2.f * k2Rinv1GeVf)) { -// pass = pass and passChiSquaredConstraint(modulesInGPU, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, chiSquared); -// if(not pass) return pass; + pass = pass and passChiSquaredConstraint(modulesInGPU, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, chiSquared); + if(not pass) return pass; } //compute the other chisquared @@ -834,7 +834,7 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc } // rzChiSquared = 12*(residual4 * residual4 + residual5 * residual5 + residual_missing * residual_missing); - if (isnan(rzChiSquared)) printf("rzChi2: %f, residual2: %f, inner_pt:%f, pseudo_phi: %f, charge: %i, Px:%f, Py:%f, x1:%f, x2:%f, x3:%f, x4:%f, x5:%f, y1:%f, y2:%f, y3:%f, y4:%f, y5:%f, z1:%f, z2:%f, z3:%f, z4:%f, z5:%f, x_center:%f, y_center:%f, slope1c:%f, slope3c:%f\n", rzChiSquared, residual_missing, inner_pt, pseudo_phi, charge, Px, Py, x1, x2, x3, x4, x5, y1, y2, y3, y4, y5, z1, z2, z3, z4, z5, x_center, y_center, slope1c, slope3c); +// if (isnan(rzChiSquared)) printf("rzChi2: %f, residual2: %f, inner_pt:%f, pseudo_phi: %f, charge: %i, Px:%f, Py:%f, x1:%f, x2:%f, x3:%f, x4:%f, x5:%f, y1:%f, y2:%f, y3:%f, y4:%f, y5:%f, z1:%f, z2:%f, z3:%f, z4:%f, z5:%f, x_center:%f, y_center:%f, slope1c:%f, slope3c:%f\n", rzChiSquared, residual_missing, inner_pt, pseudo_phi, charge, Px, Py, x1, x2, x3, x4, x5, y1, y2, y3, y4, y5, z1, z2, z3, z4, z5, x_center, y_center, slope1c, slope3c); if(layer1 == 1 and layer2 == 2 and layer3 == 3 and layer4 == 4 and layer5 == 5 and rzChiSquared>100){ // printf("rt1:%f, rt2:%f, rt3:%f, rt4:%f, rt5:%f\n", rt1, rt2, rt3, rt4, rt5); From e5b38a325e04d2beaa6d428d6785f0bc4afcae3b Mon Sep 17 00:00:00 2001 From: YonsiG Date: Fri, 6 Jan 2023 03:55:37 -0500 Subject: [PATCH 14/25] fix 0 for tilted modules with drdz>1 --- SDL/Quintuplet.cu | 84 +++++++++++++++++++++++------------------------ 1 file changed, 41 insertions(+), 43 deletions(-) diff --git a/SDL/Quintuplet.cu b/SDL/Quintuplet.cu index d9b097ea..1a9d0a1e 100644 --- a/SDL/Quintuplet.cu +++ b/SDL/Quintuplet.cu @@ -397,9 +397,9 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, bridgeRadius = computeRadiusFromThreeAnchorHits(x2, y2, x3, y3, x4, y4, g, f); innerRadius = computeRadiusFromThreeAnchorHits(x1, y1, x2, y2, x3, y3, g, f); -// pass = pass and passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f); float inner_pt = 2 * k2Rinv1GeVf * innerRadius; pass = pass and passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f); +// passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f); if(not pass) return pass; @@ -762,34 +762,32 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc float diffr=0, diffz=0; float rou = a/p; - if (layeri>6){ - float s = (zsi-z_init)*p/Pz; - float x = x_init + Px/a*sin(rou*s)-Py/a*(1-cos(rou*s)); - float y = y_init + Py/a*sin(rou*s)+Px/a*(1-cos(rou*s)); - diffr = (rtsi-sqrt(x*x+y*y))*100; - if (i==4) expectrt4=sqrt(x*x+y*y); - if (i==5) expectrt5=sqrt(x*x+y*y); - } - - if (layeri<=6){ - float paraA = rt_init*rt_init + 2*(Px*Px+Py*Py)/(a*a) + 2*(y_init*Px-x_init*Py)/a - rtsi*rtsi; - float paraB = 2*(x_init*Px+y_init*Py)/a; - float paraC = 2*(y_init*Px-x_init*Py)/a+2*(Px*Px+Py*Py)/(a*a); - float A=paraB*paraB+paraC*paraC; - float B=2*paraA*paraB; - float C=paraA*paraA-paraC*paraC; - float sol1 = (-B+sqrt(B*B-4*A*C))/(2*A); - float sol2 = (-B-sqrt(B*B-4*A*C))/(2*A); - float solz1 = asin(sol1)/rou*Pz/p+z_init; - float solz2 = asin(sol2)/rou*Pz/p+z_init; - float diffz1 = (solz1-zsi)*100; - float diffz2 = (solz2-zsi)*100; - diffz = (fabs(diffz1)fabs(diffz2)) expectz4 = solz2; - if (i==5 && fabs(diffz1)fabs(diffz2)) expectz5 = solz2; - } + // for barrel + float s = (zsi-z_init)*p/Pz; + float x = x_init + Px/a*sin(rou*s)-Py/a*(1-cos(rou*s)); + float y = y_init + Py/a*sin(rou*s)+Px/a*(1-cos(rou*s)); + diffr = (rtsi-sqrt(x*x+y*y))*100; + if (i==4) expectrt4=sqrt(x*x+y*y); + if (i==5) expectrt5=sqrt(x*x+y*y); + + // for endcap + float paraA = rt_init*rt_init + 2*(Px*Px+Py*Py)/(a*a) + 2*(y_init*Px-x_init*Py)/a - rtsi*rtsi; + float paraB = 2*(x_init*Px+y_init*Py)/a; + float paraC = 2*(y_init*Px-x_init*Py)/a+2*(Px*Px+Py*Py)/(a*a); + float A=paraB*paraB+paraC*paraC; + float B=2*paraA*paraB; + float C=paraA*paraA-paraC*paraC; + float sol1 = (-B+sqrt(B*B-4*A*C))/(2*A); + float sol2 = (-B-sqrt(B*B-4*A*C))/(2*A); + float solz1 = asin(sol1)/rou*Pz/p+z_init; + float solz2 = asin(sol2)/rou*Pz/p+z_init; + float diffz1 = (solz1-zsi)*100; + float diffz2 = (solz2-zsi)*100; + diffz = (fabs(diffz1)fabs(diffz2)) expectz4 = solz2; + if (i==5 && fabs(diffz1)fabs(diffz2)) expectz5 = solz2; residual = (layeri>6) ? diffr : diffz ; @@ -814,7 +812,6 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc subdets=modulesInGPU.subdets[lowerModuleIndex2]; } if (i==3){ - drdz=abs(modulesInGPU.drdzs[lowerModuleIndex3]); side=modulesInGPU.sides[lowerModuleIndex3]; subdets=modulesInGPU.subdets[lowerModuleIndex3]; @@ -836,54 +833,55 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc // if (isnan(rzChiSquared)) printf("rzChi2: %f, residual2: %f, inner_pt:%f, pseudo_phi: %f, charge: %i, Px:%f, Py:%f, x1:%f, x2:%f, x3:%f, x4:%f, x5:%f, y1:%f, y2:%f, y3:%f, y4:%f, y5:%f, z1:%f, z2:%f, z3:%f, z4:%f, z5:%f, x_center:%f, y_center:%f, slope1c:%f, slope3c:%f\n", rzChiSquared, residual_missing, inner_pt, pseudo_phi, charge, Px, Py, x1, x2, x3, x4, x5, y1, y2, y3, y4, y5, z1, z2, z3, z4, z5, x_center, y_center, slope1c, slope3c); - if(layer1 == 1 and layer2 == 2 and layer3 == 3 and layer4 == 4 and layer5 == 5 and rzChiSquared>100){ +// if(layer1 == 1 and layer2 == 2 and layer3 == 3 and layer4 == 12 and layer5 == 13){ // printf("rt1:%f, rt2:%f, rt3:%f, rt4:%f, rt5:%f\n", rt1, rt2, rt3, rt4, rt5); // printf("x1:%f, x2:%f, x3:%f, x4:%f, x5:%f\n", x1, x2, x3, x4, x5); // printf("y1:%f, y2:%f, y3:%f, y4:%f, y5:%f\n", y1, y2, y3, y4, y5); // printf("z1:%f, z2:%f, z3:%f, z4:%f, z5:%f\n", z1, z2, z3, z4, z5); // printf("rt4_ex:%f, rt5_ex:%f\n", expectrt4, expectrt5); // printf("z4_ex:%f, z5_ex:%f\n", expectz4, expectz5); -// printf("Pt:%f, Px:%f, Py:%f, Pz:%f, charge: %i, residual4: %f, residual5:%f\n", Pt, Px, Py, Pz, charge, residual4, residual5); +// printf("residual_missing:%f\n", residual_missing); +// printf("Pt:%f, Px:%f, Py:%f, Pz:%f, charge: %i, residual_missing: %f, residual4: %f, residual5:%f, moduleType3:%i\n", Pt, Px, Py, Pz, charge, residual_missing, residual4, residual5, moduleType3); // printf("rzChi2: %f, residual2: %f, inner_pt:%f, pseudo_phi: %f, charge: %i, Px:%f, Py:%f, x1:%f, x2:%f, x3:%f, x4:%f, x5:%f, y1:%f, y2:%f, y3:%f, y4:%f, y5:%f, z1:%f, z2:%f, z3:%f, z4:%f, z5:%f, x_center:%f, y_center:%f, slope1c:%f, slope3c:%f\n", rzChiSquared, residual_missing, inner_pt, pseudo_phi, charge, Px, Py, x1, x2, x3, x4, x5, y1, y2, y3, y4, y5, z1, z2, z3, z4, z5, x_center, y_center, slope1c, slope3c); - } +// } //categories! if(layer1 == 1 and layer2 == 2 and layer3 == 3) { if(layer4 == 4 and layer5 == 5) { - return rzChiSquared < 28.3f; //11 + return rzChiSquared < 30.545f; //11 } else if(layer4 == 4 and layer5 == 12) //12 { - return rzChiSquared < 4.475f; + return rzChiSquared < 14.765f; } else if(layer4 == 7 and layer5 == 8) //8 { - return rzChiSquared < 44.243f; + return rzChiSquared < 44.247f; } else if(layer4 == 7 and layer5 == 13) //9 { - return rzChiSquared < 30.405f; + return rzChiSquared < 43.578f; } else if(layer4 == 12 and layer5 == 13) //10 { - return rzChiSquared < 4.331f; + return rzChiSquared < 14.377f; } } else if(layer1 == 1 and layer2 == 2 and layer3 == 7) { if(layer4 == 8 and layer5 == 9) //5 { - return rzChiSquared < 113.184f; + return rzChiSquared < 115.989f; } if(layer4 == 8 and layer5 == 14) //6 { - return rzChiSquared < 38.839f; + return rzChiSquared < 56.447f; } else if(layer4 == 13 and layer5 == 14) //7 { - return rzChiSquared < 4.868f; + return rzChiSquared < 10.938f; } } else if(layer1 == 1 and layer2 == 7 and layer3 == 8 and layer4 == 9) @@ -917,11 +915,11 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc { if(layer4 == 8 and layer5 == 14) //16 { - return rzChiSquared < 18.649f; + return rzChiSquared < 23.748f; } if(layer4 == 13 and layer5 == 14) //17 { - return rzChiSquared < 4.661f; + return rzChiSquared < 17.975f; } } From c0c81d2397d4827e51438cd9d6a9891245c1c5b0 Mon Sep 17 00:00:00 2001 From: YonsiG Date: Sat, 7 Jan 2023 02:42:47 -0500 Subject: [PATCH 15/25] add 95% cut on TC collection for T5 --- SDL/Quintuplet.cu | 65 +++++++++++++++++++++++++++++++------------ SDL/Quintuplet.cuh | 7 +++-- SDL/TrackCandidate.cu | 4 ++- 3 files changed, 54 insertions(+), 22 deletions(-) diff --git a/SDL/Quintuplet.cu b/SDL/Quintuplet.cu index 1a9d0a1e..74f0ba3f 100644 --- a/SDL/Quintuplet.cu +++ b/SDL/Quintuplet.cu @@ -15,6 +15,7 @@ SDL::quintuplets::quintuplets() outerRadius = nullptr; regressionRadius = nullptr; isDup = nullptr; + TightCutFlag = nullptr; partOfPT5 = nullptr; pt = nullptr; layer = nullptr; @@ -48,6 +49,7 @@ void SDL::quintuplets::freeMemoryCache() cms::cuda::free_device(dev, outerRadius); cms::cuda::free_device(dev, partOfPT5); cms::cuda::free_device(dev, isDup); + cms::cuda::free_device(dev, TightCutFlag); cms::cuda::free_device(dev, pt); cms::cuda::free_device(dev, layer); cms::cuda::free_device(dev, regressionG); @@ -77,6 +79,7 @@ void SDL::quintuplets::freeMemory(cudaStream_t stream) cudaFree(regressionRadius); cudaFree(partOfPT5); cudaFree(isDup); + cudaFree(TightCutFlag); cudaFree(pt); cudaFree(layer); cudaFree(regressionG); @@ -174,6 +177,7 @@ void SDL::createQuintupletsInExplicitMemory(struct SDL::quintuplets& quintuplets quintupletsInGPU.pt = (FPX*)cms::cuda::allocate_device(dev, nTotalQuintuplets *4* sizeof(FPX), stream); quintupletsInGPU.layer = (uint8_t*)cms::cuda::allocate_device(dev, nTotalQuintuplets * sizeof(uint8_t), stream); quintupletsInGPU.isDup = (bool*)cms::cuda::allocate_device(dev, nTotalQuintuplets * sizeof(bool), stream); + quintupletsInGPU.TightCutFlag = (bool*)cms::cuda::allocate_device(dev, nTotalQuintuplets * sizeof(bool), stream); quintupletsInGPU.partOfPT5 = (bool*)cms::cuda::allocate_device(dev, nTotalQuintuplets * sizeof(bool), stream); quintupletsInGPU.regressionRadius = (float*)cms::cuda::allocate_device(dev, nTotalQuintuplets * sizeof(float), stream); quintupletsInGPU.regressionG = (float*)cms::cuda::allocate_device(dev, nTotalQuintuplets * sizeof(float), stream); @@ -198,6 +202,7 @@ void SDL::createQuintupletsInExplicitMemory(struct SDL::quintuplets& quintuplets cudaMalloc(&quintupletsInGPU.outerRadius, nTotalQuintuplets * sizeof(FPX)); cudaMalloc(&quintupletsInGPU.pt, nTotalQuintuplets *4* sizeof(FPX)); cudaMalloc(&quintupletsInGPU.isDup, nTotalQuintuplets * sizeof(bool)); + cudaMalloc(&quintupletsInGPU.TightCutFlag, nTotalQuintuplets * sizeof(bool)); cudaMalloc(&quintupletsInGPU.partOfPT5, nTotalQuintuplets * sizeof(bool)); cudaMalloc(&quintupletsInGPU.layer, nTotalQuintuplets * sizeof(uint8_t)); cudaMalloc(&quintupletsInGPU.regressionRadius, nTotalQuintuplets * sizeof(float)); @@ -218,6 +223,7 @@ void SDL::createQuintupletsInExplicitMemory(struct SDL::quintuplets& quintuplets cudaMemsetAsync(quintupletsInGPU.nQuintuplets,0,nLowerModules * sizeof(unsigned int),stream); cudaMemsetAsync(quintupletsInGPU.totOccupancyQuintuplets,0,nLowerModules * sizeof(unsigned int),stream); cudaMemsetAsync(quintupletsInGPU.isDup,0,nTotalQuintuplets * sizeof(bool),stream); + cudaMemsetAsync(quintupletsInGPU.TightCutFlag,0,nTotalQuintuplets * sizeof(bool),stream); cudaMemsetAsync(quintupletsInGPU.partOfPT5,0,nTotalQuintuplets * sizeof(bool),stream); cudaStreamSynchronize(stream); quintupletsInGPU.eta = quintupletsInGPU.pt + nTotalQuintuplets; @@ -227,7 +233,7 @@ void SDL::createQuintupletsInExplicitMemory(struct SDL::quintuplets& quintuplets __device__ void SDL::addQuintupletToMemory(struct SDL::triplets& tripletsInGPU, struct SDL::quintuplets& quintupletsInGPU, unsigned int innerTripletIndex, unsigned int outerTripletIndex, uint16_t& lowerModule1, uint16_t& lowerModule2, uint16_t& lowerModule3, uint16_t& lowerModule4, uint16_t& lowerModule5, float& innerRadius, float& bridgeRadius, float& outerRadius, float& regressionG, float& regressionF, float& regressionRadius, float& rzChiSquared, float& rPhiChiSquared, float& - nonAnchorChiSquared, float residual_missing, float residual4, float residual5, float pt, float eta, float phi, float scores, uint8_t layer, unsigned int quintupletIndex) + nonAnchorChiSquared, float residual_missing, float residual4, float residual5, float pt, float eta, float phi, float scores, uint8_t layer, unsigned int quintupletIndex, bool TightCutFlag) { quintupletsInGPU.tripletIndices[2 * quintupletIndex] = innerTripletIndex; @@ -246,6 +252,7 @@ __device__ void SDL::addQuintupletToMemory(struct SDL::triplets& tripletsInGPU, quintupletsInGPU.score_rphisum[quintupletIndex] = __F2H(scores); quintupletsInGPU.layer[quintupletIndex] = layer; quintupletsInGPU.isDup[quintupletIndex] = false; + quintupletsInGPU.TightCutFlag[quintupletIndex] = TightCutFlag; quintupletsInGPU.regressionRadius[quintupletIndex] = regressionRadius; quintupletsInGPU.regressionG[quintupletIndex] = regressionG; quintupletsInGPU.regressionF[quintupletIndex] = regressionF; @@ -274,7 +281,7 @@ __device__ void SDL::addQuintupletToMemory(struct SDL::triplets& tripletsInGPU, quintupletsInGPU.residual5[quintupletIndex] = residual5; } -__device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, struct SDL::segments& segmentsInGPU, struct SDL::triplets& tripletsInGPU, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, unsigned int& innerTripletIndex, unsigned int& outerTripletIndex, float& innerRadius, float& outerRadius, float& bridgeRadius, float& regressionG, float& regressionF, float& regressionRadius, float& rzChiSquared, float& chiSquared, float& nonAnchorChiSquared, float& residual_missing, float& residual4, float& residual5) +__device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, struct SDL::segments& segmentsInGPU, struct SDL::triplets& tripletsInGPU, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, unsigned int& innerTripletIndex, unsigned int& outerTripletIndex, float& innerRadius, float& outerRadius, float& bridgeRadius, float& regressionG, float& regressionF, float& regressionRadius, float& rzChiSquared, float& chiSquared, float& nonAnchorChiSquared, float& residual_missing, float& residual4, float& residual5, bool& TightCutFlag) { bool pass = true; unsigned int firstSegmentIndex = tripletsInGPU.segmentIndices[2 * innerTripletIndex]; @@ -398,7 +405,7 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, innerRadius = computeRadiusFromThreeAnchorHits(x1, y1, x2, y2, x3, y3, g, f); float inner_pt = 2 * k2Rinv1GeVf * innerRadius; - pass = pass and passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f); + pass = pass and passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f, TightCutFlag); // passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f); if(not pass) return pass; @@ -609,7 +616,7 @@ __device__ bool SDL::passChiSquaredConstraint(struct SDL::modules& modulesInGPU, } //bounds can be found at http://uaf-10.t2.ucsd.edu/~bsathian/SDL/T5_RZFix/t5_rz_thresholds.txt -__device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex, unsigned int fifthMDIndex, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, float& rzChiSquared, float& residual_missing, float& residual4, float& residual5, float inner_pt, float innerRadius, float g, float f) +__device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex, unsigned int fifthMDIndex, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, float& rzChiSquared, float& residual_missing, float& residual4, float& residual5, float inner_pt, float innerRadius, float g, float f, bool& TightCutFlag) { //(g,f) is the center of the circle fitted by the innermost 3 points on x,y coordinates const float& rt1 = mdsInGPU.anchorRt[firstMDIndex]/100; //in the unit of m instead of cm @@ -845,53 +852,65 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc // printf("rzChi2: %f, residual2: %f, inner_pt:%f, pseudo_phi: %f, charge: %i, Px:%f, Py:%f, x1:%f, x2:%f, x3:%f, x4:%f, x5:%f, y1:%f, y2:%f, y3:%f, y4:%f, y5:%f, z1:%f, z2:%f, z3:%f, z4:%f, z5:%f, x_center:%f, y_center:%f, slope1c:%f, slope3c:%f\n", rzChiSquared, residual_missing, inner_pt, pseudo_phi, charge, Px, Py, x1, x2, x3, x4, x5, y1, y2, y3, y4, y5, z1, z2, z3, z4, z5, x_center, y_center, slope1c, slope3c); // } + // when building T5, apply 99% chi2 cuts as default, and add to pT5 collection. But when adding T5 to TC collections, appy 95% cut to reduce the fake rate + TightCutFlag = false; //categories! if(layer1 == 1 and layer2 == 2 and layer3 == 3) { - if(layer4 == 4 and layer5 == 5) + if(layer4 == 4 and layer5 == 5) //11 { - return rzChiSquared < 30.545f; //11 + if (rzChiSquared < 15.754f) TightCutFlag = 1; + return rzChiSquared < 30.545f; } else if(layer4 == 4 and layer5 == 12) //12 { - return rzChiSquared < 14.765f; + if (rzChiSquared < 14.765f) TightCutFlag = 1; + return rzChiSquared < 23.704f; } else if(layer4 == 7 and layer5 == 8) //8 { + if (rzChiSquared < 27.824f) TightCutFlag = 1; return rzChiSquared < 44.247f; } else if(layer4 == 7 and layer5 == 13) //9 { + if (rzChiSquared < 19.658f) TightCutFlag = 1; return rzChiSquared < 43.578f; } else if(layer4 == 12 and layer5 == 13) //10 { - return rzChiSquared < 14.377f; + if (rzChiSquared < 14.377f) TightCutFlag = 1; + return rzChiSquared < 43.099f; } } else if(layer1 == 1 and layer2 == 2 and layer3 == 7) { if(layer4 == 8 and layer5 == 9) //5 { + if (rzChiSquared < 60.409f) TightCutFlag = 1; return rzChiSquared < 115.989f; } if(layer4 == 8 and layer5 == 14) //6 { + if (rzChiSquared < 19.603f) TightCutFlag = 1; return rzChiSquared < 56.447f; } else if(layer4 == 13 and layer5 == 14) //7 { - return rzChiSquared < 10.938f; + if (rzChiSquared < 10.938f) TightCutFlag = 1; + return rzChiSquared < 50.731f; } } else if(layer1 == 1 and layer2 == 7 and layer3 == 8 and layer4 == 9) { if (layer5 == 10) //3 { + if (rzChiSquared < 63.724f) TightCutFlag = 1; return rzChiSquared < 109.577f; } if (layer5 == 15) //4 { + if (rzChiSquared < 18.351f) TightCutFlag = 1; return rzChiSquared < 37.933f; } } @@ -899,26 +918,31 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc { if(layer4 == 5 and layer5 == 6) //18 { - return rzChiSquared < 7.951f; + if (rzChiSquared < 7.951f) TightCutFlag = 1; + return rzChiSquared < 40.849f; } else if(layer4 == 5 and layer5 == 12) //19 { - return rzChiSquared < 14.413f; + if (rzChiSquared < 14.413f) TightCutFlag = 1; + return rzChiSquared < 74.155f; } else if(layer4 == 12 and layer5 == 13) //20 { - return rzChiSquared < 12.310f; + if (rzChiSquared < 12.310f) TightCutFlag = 1; + return rzChiSquared < 57.153f; } } else if(layer1 == 2 and layer2 == 3 and layer3 == 7) { if(layer4 == 8 and layer5 == 14) //16 { + if (rzChiSquared < 23.730f) TightCutFlag = 1; return rzChiSquared < 23.748f; } if(layer4 == 13 and layer5 == 14) //17 { + if (rzChiSquared < 10.585f) TightCutFlag = 1; return rzChiSquared < 17.975f; } } @@ -927,27 +951,32 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc { if(layer4 == 9 and layer5 == 15) //14 { + if (rzChiSquared < 24.576f) TightCutFlag = 1; return rzChiSquared < 41.036f; } else if(layer4 == 14 and layer5 == 15) //15 { + if (rzChiSquared < 8.785f) TightCutFlag = 1; return rzChiSquared < 13.821f; } } else if(layer1 == 7 and layer2 == 8 and layer3 == 9 and layer4 == 15 and layer5 == 16) //2 { - return rzChiSquared < 12.223f; + if (rzChiSquared < 8.002f) TightCutFlag = 1; + return rzChiSquared < 12.223f; } else if(layer1 == 7 and layer2 == 8 and layer3 == 9 and layer4 == 10 and layer5 == 11) //0 { - return rzChiSquared < 93.893f; + if (rzChiSquared < 56.138f) TightCutFlag = 1; + return rzChiSquared < 93.893f; } else if(layer1 == 7 and layer2 == 8 and layer3 == 9 and layer4 == 10 and layer5 == 16) //1 { - return rzChiSquared < 37.087f; + if (rzChiSquared < 22.041f) TightCutFlag = 1; + return rzChiSquared < 37.087f; } return true; } @@ -1496,8 +1525,8 @@ __global__ void SDL::createQuintupletsInGPUv2(struct SDL::modules& modulesInGPU, uint16_t lowerModule5 = tripletsInGPU.lowerModuleIndices[3 * outerTripletIndex + 2]; float innerRadius, outerRadius, bridgeRadius, regressionG, regressionF, regressionRadius, rzChiSquared, chiSquared, nonAnchorChiSquared, residual_missing, residual4, residual; //required for making distributions - - bool success = runQuintupletDefaultAlgo(modulesInGPU, mdsInGPU, segmentsInGPU, tripletsInGPU, lowerModule1, lowerModule2, lowerModule3, lowerModule4, lowerModule5, innerTripletIndex, outerTripletIndex, innerRadius, outerRadius, bridgeRadius, regressionG, regressionF, regressionRadius, rzChiSquared, chiSquared, nonAnchorChiSquared, residual_missing, residual4, residual5); + bool TightCutFlag; + bool success = runQuintupletDefaultAlgo(modulesInGPU, mdsInGPU, segmentsInGPU, tripletsInGPU, lowerModule1, lowerModule2, lowerModule3, lowerModule4, lowerModule5, innerTripletIndex, outerTripletIndex, innerRadius, outerRadius, bridgeRadius, regressionG, regressionF, regressionRadius, rzChiSquared, chiSquared, nonAnchorChiSquared, residual_missing, residual4, residual5, TightCutFlag); if(success) { @@ -1537,7 +1566,7 @@ __global__ void SDL::createQuintupletsInGPUv2(struct SDL::modules& modulesInGPU, float eta = mdsInGPU.anchorEta[segmentsInGPU.mdIndices[2*tripletsInGPU.segmentIndices[2*innerTripletIndex+layer2_adjustment]]]; float pt = (innerRadius+outerRadius)*3.8f*1.602f/(2*100*5.39f); float scores = chiSquared + nonAnchorChiSquared; - addQuintupletToMemory(tripletsInGPU, quintupletsInGPU, innerTripletIndex, outerTripletIndex, lowerModule1, lowerModule2, lowerModule3, lowerModule4, lowerModule5, innerRadius, bridgeRadius, outerRadius, regressionG, regressionF, regressionRadius, rzChiSquared, chiSquared, nonAnchorChiSquared, rzChiSquared, residual_missing, residual4, residual5, pt,eta,phi,scores,layer,quintupletIndex); + addQuintupletToMemory(tripletsInGPU, quintupletsInGPU, innerTripletIndex, outerTripletIndex, lowerModule1, lowerModule2, lowerModule3, lowerModule4, lowerModule5, innerRadius, bridgeRadius, outerRadius, regressionG, regressionF, regressionRadius, rzChiSquared, chiSquared, nonAnchorChiSquared, rzChiSquared, residual_missing, residual4, residual5, pt,eta,phi,scores,layer,quintupletIndex, TightCutFlag); tripletsInGPU.partOfT5[quintupletsInGPU.tripletIndices[2 * quintupletIndex]] = true; tripletsInGPU.partOfT5[quintupletsInGPU.tripletIndices[2 * quintupletIndex + 1]] = true; diff --git a/SDL/Quintuplet.cuh b/SDL/Quintuplet.cuh index acd5f1d3..84e26145 100644 --- a/SDL/Quintuplet.cuh +++ b/SDL/Quintuplet.cuh @@ -40,6 +40,7 @@ namespace SDL FPX* score_rphisum; uint8_t* layer; bool* isDup; + bool* TightCutFlag; bool* partOfPT5; float* regressionRadius; @@ -69,12 +70,12 @@ namespace SDL // CUDA_DEV void rmQuintupletToMemory(struct SDL::quintuplets& quintupletsInGPU, unsigned int quintupletIndex); - CUDA_DEV void addQuintupletToMemory(struct SDL::triplets& tripletsInGPU, struct SDL::quintuplets& quintupletsInGPU, unsigned int innerTripletIndex, unsigned int outerTripletIndex, uint16_t& lowerModule1, uint16_t& lowerModule2, uint16_t& lowerModule3, uint16_t& lowerModule4, uint16_t& lowerModule5, float& innerRadius, float& bridgeRadius, float& outerRadius, float& regressionG, float& regressionF, float& regressionRadius, float& rzChiSquared, float& rPhiChiSquared, float& nonAnchorChiSquared, float residual_missing, float residual4, float residual5, float pt, float eta, float phi, float scores, uint8_t layer, unsigned int quintupletIndex); + CUDA_DEV void addQuintupletToMemory(struct SDL::triplets& tripletsInGPU, struct SDL::quintuplets& quintupletsInGPU, unsigned int innerTripletIndex, unsigned int outerTripletIndex, uint16_t& lowerModule1, uint16_t& lowerModule2, uint16_t& lowerModule3, uint16_t& lowerModule4, uint16_t& lowerModule5, float& innerRadius, float& bridgeRadius, float& outerRadius, float& regressionG, float& regressionF, float& regressionRadius, float& rzChiSquared, float& rPhiChiSquared, float& nonAnchorChiSquared, float residual_missing, float residual4, float residual5, float pt, float eta, float phi, float scores, uint8_t layer, unsigned int quintupletIndex, bool TightCutFlag); - CUDA_DEV bool runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, struct SDL::segments& segmentsInGPU, struct SDL::triplets& tripletsInGPU, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, unsigned int& innerTripletIndex, unsigned int& outerTripletIndex, float& innerRadius, float& outerRadius, float& bridgeRadius, float& regressionG, float& regressionF, float& regressionRadius, float& rzChiSquared, float& chiSquared, float& nonAnchorChiSquared, float& residual_missing, float& residual4, float& residual5); + CUDA_DEV bool runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, struct SDL::segments& segmentsInGPU, struct SDL::triplets& tripletsInGPU, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, unsigned int& innerTripletIndex, unsigned int& outerTripletIndex, float& innerRadius, float& outerRadius, float& bridgeRadius, float& regressionG, float& regressionF, float& regressionRadius, float& rzChiSquared, float& chiSquared, float& nonAnchorChiSquared, float& residual_missing, float& residual4, float& residual5, bool& TightCutFlag); - CUDA_DEV bool passT5RZConstraint(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex, unsigned int fifthMDIndex, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, float& rzChiSquared, float& residual_missing, float& residual4, float& residual5, float inner_pt, float innerRadius, float g, float f); + CUDA_DEV bool passT5RZConstraint(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex, unsigned int fifthMDIndex, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, float& rzChiSquared, float& residual_missing, float& residual4, float& residual5, float inner_pt, float innerRadius, float g, float f, bool& TightCutFlag); CUDA_DEV bool T5HasCommonMiniDoublet(struct triplets& tripletsInGPU, struct segments& segmentsInGPU, unsigned int innerTripletIndex, unsigned int outerTripletIndex); diff --git a/SDL/TrackCandidate.cu b/SDL/TrackCandidate.cu index 16752b8e..f4b2202e 100644 --- a/SDL/TrackCandidate.cu +++ b/SDL/TrackCandidate.cu @@ -355,7 +355,9 @@ __global__ void SDL::addT5asTrackCandidateInGPU(struct SDL::modules& modulesInGP { int quintupletIndex = rangesInGPU.quintupletModuleIndices[idx] + jdx; - if(quintupletsInGPU.isDup[quintupletIndex] or quintupletsInGPU.partOfPT5[quintupletIndex]) continue; + if (quintupletsInGPU.isDup[quintupletIndex] or quintupletsInGPU.partOfPT5[quintupletIndex]) continue; + + if (!(quintupletsInGPU.TightCutFlag[quintupletIndex])) continue; unsigned int trackCandidateIdx = atomicAdd(trackCandidatesInGPU.nTrackCandidates,1); atomicAdd(trackCandidatesInGPU.nTrackCandidatesT5,1); From 7d9fdcfe699a5260b26270465c04493c458a0f2f Mon Sep 17 00:00:00 2001 From: YonsiG Date: Sat, 7 Jan 2023 08:08:06 -0500 Subject: [PATCH 16/25] new method for Px Py initial for t5 --- SDL/Quintuplet.cu | 102 ++++++++++++++++++++++++++-------------------- 1 file changed, 57 insertions(+), 45 deletions(-) diff --git a/SDL/Quintuplet.cu b/SDL/Quintuplet.cu index 74f0ba3f..6dcda094 100644 --- a/SDL/Quintuplet.cu +++ b/SDL/Quintuplet.cu @@ -405,8 +405,8 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, innerRadius = computeRadiusFromThreeAnchorHits(x1, y1, x2, y2, x3, y3, g, f); float inner_pt = 2 * k2Rinv1GeVf * innerRadius; - pass = pass and passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f, TightCutFlag); -// passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f); +// pass = pass and passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f, TightCutFlag); + passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f, TightCutFlag); if(not pass) return pass; @@ -454,7 +454,7 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, } //compute regression radius right here - this computation is expensive!!! - pass = pass and tempPass; +// pass = pass and tempPass; if(not pass) return pass; float xVec[] = {x1, x2, x3, x4, x5}; @@ -470,7 +470,7 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, //extra chi squared cuts! if(regressionRadius < 5.0f/(2.f * k2Rinv1GeVf)) { - pass = pass and passChiSquaredConstraint(modulesInGPU, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, chiSquared); +// pass = pass and passChiSquaredConstraint(modulesInGPU, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, chiSquared); if(not pass) return pass; } @@ -704,14 +704,25 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc float pseudo_phi = atan((y_init-y_center)/(x_init-x_center)); //actually represent pi/2-phi, wrt helix axis z float Pt=inner_pt, Px=Pt*abs(sin(pseudo_phi)), Py=Pt*abs(cos(pseudo_phi)); -// if (x_init<0 && y_init<0) {Px = -Pt*abs(sin(pseudo_phi)); Py=-Pt*abs(cos(pseudo_phi));} -// else if (x_init<0 && y_init>0) {Px = -Pt*abs(sin(pseudo_phi)); Py=Pt*abs(cos(pseudo_phi));} -// else if (x_init>0 && y_init>0) {Px = Pt*abs(sin(pseudo_phi)); Py=Pt*abs(cos(pseudo_phi));} -// else if (x_init>0 && y_init<0) {Px = Pt*abs(sin(pseudo_phi)); Py=-Pt*abs(cos(pseudo_phi));} - if (x4-x2<0) Px = -Px; - if (y4-y2<0) Py = -Py; -// else return 0; - + + if (x_init>x_center && y_init>y_center) //1st quad + { + if (charge==1) Py=-Py; + if (charge==-1) Px=-Px; + } + if (x_inity_center) //2nd quad + { + if (charge==-1) {Px=-Px; Py=-Py;} + } + if (x_initx_center && y_init Date: Sun, 8 Jan 2023 21:56:20 -0500 Subject: [PATCH 17/25] small fix --- SDL/Quintuplet.cu | 54 +++++++++++++++++++++++++++++------------------ 1 file changed, 34 insertions(+), 20 deletions(-) diff --git a/SDL/Quintuplet.cu b/SDL/Quintuplet.cu index 6dcda094..658cb789 100644 --- a/SDL/Quintuplet.cu +++ b/SDL/Quintuplet.cu @@ -724,6 +724,20 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc if (charge==1) {Px=-Px; Py=-Py;} } + if (moduleType3==0){ + if (x4x3 && x3>x2) Px=abs(Px); + if (y4y3 && y3>y2) Py=abs(Py); + } + else (moduleType3==1) + { + if (x3x2 && x2>x1) Px=abs(Px); + if (y3y2 && y2>y1) Py=abs(Py); + } + //to get Pz, we use pt/pz=ds/dz, ds is the arclength between MD1 and MD3. float AO=sqrt((x1-x_center)*(x1-x_center)+(y1-y_center)*(y1-y_center)); float BO=sqrt((x_init-x_center)*(x_init-x_center)+(y_init-y_center)*(y_init-y_center)); @@ -780,7 +794,7 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc float diffr=0, diffz=0; float rou = a/p; - // for barrel + // for endcap float s = (zsi-z_init)*p/Pz; float x = x_init + Px/a*sin(rou*s)-Py/a*(1-cos(rou*s)); float y = y_init + Py/a*sin(rou*s)+Px/a*(1-cos(rou*s)); @@ -788,25 +802,25 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc if (i==4) expectrt4=sqrt(x*x+y*y); if (i==5) expectrt5=sqrt(x*x+y*y); - // for endcap - float paraA = rt_init*rt_init + 2*(Px*Px+Py*Py)/(a*a) + 2*(y_init*Px-x_init*Py)/a - rtsi*rtsi; - float paraB = 2*(x_init*Px+y_init*Py)/a; - float paraC = 2*(y_init*Px-x_init*Py)/a+2*(Px*Px+Py*Py)/(a*a); - float A=paraB*paraB+paraC*paraC; - float B=2*paraA*paraB; - float C=paraA*paraA-paraC*paraC; - float sol1 = (-B+sqrt(B*B-4*A*C))/(2*A); - float sol2 = (-B-sqrt(B*B-4*A*C))/(2*A); - float solz1 = asin(sol1)/rou*Pz/p+z_init; - float solz2 = asin(sol2)/rou*Pz/p+z_init; - float diffz1 = (solz1-zsi)*100; - float diffz2 = (solz2-zsi)*100; - diffz = (fabs(diffz1)fabs(diffz2)) expectz4 = solz2; - if (i==5 && fabs(diffz1)fabs(diffz2)) expectz5 = solz2; - + // for barrel + if (layeri<=6) + { + float paraA = rt_init*rt_init + 2*(Px*Px+Py*Py)/(a*a) + 2*(y_init*Px-x_init*Py)/a - rtsi*rtsi; + float paraB = 2*(x_init*Px+y_init*Py)/a; + float paraC = 2*(y_init*Px-x_init*Py)/a+2*(Px*Px+Py*Py)/(a*a); + float A=paraB*paraB+paraC*paraC; + float B=2*paraA*paraB; + float C=paraA*paraA-paraC*paraC; + float sol1 = (-B+sqrt(B*B-4*A*C))/(2*A); + float sol2 = (-B-sqrt(B*B-4*A*C))/(2*A); + float solz1 = asin(sol1)/rou*Pz/p+z_init; + float solz2 = asin(sol2)/rou*Pz/p+z_init; + float diffz1 = (solz1-zsi)*100; + float diffz2 = (solz2-zsi)*100; + if (isnan(diffz1)) diffz = diffz2; + else if (isnan(diffz2)) diffz = diffz1; + else {diffz = (fabs(diffz1)6) ? diffr : diffz ; //PS Modules From fd725fa9fcb9e529c99d2404f0aa910533e2b5a8 Mon Sep 17 00:00:00 2001 From: YonsiG Date: Tue, 10 Jan 2023 04:15:19 -0500 Subject: [PATCH 18/25] add tight/loose cut set for T5 --- SDL/Quintuplet.cu | 86 +++++++++++++++++++++++------------------------ 1 file changed, 43 insertions(+), 43 deletions(-) diff --git a/SDL/Quintuplet.cu b/SDL/Quintuplet.cu index 658cb789..6b60f328 100644 --- a/SDL/Quintuplet.cu +++ b/SDL/Quintuplet.cu @@ -311,10 +311,10 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, pass = pass and runQuintupletDefaultAlgo(modulesInGPU, mdsInGPU, segmentsInGPU, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex4, lowerModuleIndex5, firstSegmentIndex, fourthSegmentIndex, firstMDIndex, secondMDIndex, fourthMDIndex, fifthMDIndex, zOut, rtOut, deltaPhiPos, deltaPhi, betaIn, betaOut, pt_beta, zLo, zHi, rtLo, rtHi, zLoPointed, zHiPointed, sdlCut, betaInCut, betaOutCut, deltaBetaCut, kZ); if(not pass) return pass; - pass = pass and runQuintupletDefaultAlgo(modulesInGPU, mdsInGPU, segmentsInGPU, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex3, lowerModuleIndex4, secondSegmentIndex, thirdSegmentIndex, secondMDIndex, thirdMDIndex, thirdMDIndex, fourthMDIndex, zOut, rtOut, deltaPhiPos, deltaPhi, betaIn, betaOut, pt_beta, zLo, zHi, rtLo, rtHi, zLoPointed, zHiPointed, sdlCut, betaInCut, betaOutCut, deltaBetaCut, kZ); +// pass = pass and runQuintupletDefaultAlgo(modulesInGPU, mdsInGPU, segmentsInGPU, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex3, lowerModuleIndex4, secondSegmentIndex, thirdSegmentIndex, secondMDIndex, thirdMDIndex, thirdMDIndex, fourthMDIndex, zOut, rtOut, deltaPhiPos, deltaPhi, betaIn, betaOut, pt_beta, zLo, zHi, rtLo, rtHi, zLoPointed, zHiPointed, sdlCut, betaInCut, betaOutCut, deltaBetaCut, kZ); if(not pass) return pass; - pass = pass and runQuintupletDefaultAlgo(modulesInGPU, mdsInGPU, segmentsInGPU, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, secondSegmentIndex, fourthSegmentIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, zOut, rtOut, deltaPhiPos, deltaPhi, betaIn, betaOut, pt_beta, zLo, zHi, rtLo, rtHi, zLoPointed, zHiPointed, sdlCut, betaInCut, betaOutCut, deltaBetaCut, kZ); +// pass = pass and runQuintupletDefaultAlgo(modulesInGPU, mdsInGPU, segmentsInGPU, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, secondSegmentIndex, fourthSegmentIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, zOut, rtOut, deltaPhiPos, deltaPhi, betaIn, betaOut, pt_beta, zLo, zHi, rtLo, rtHi, zLoPointed, zHiPointed, sdlCut, betaInCut, betaOutCut, deltaBetaCut, kZ); if(not pass) return pass; float x1 = mdsInGPU.anchorX[firstMDIndex]; @@ -405,8 +405,8 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, innerRadius = computeRadiusFromThreeAnchorHits(x1, y1, x2, y2, x3, y3, g, f); float inner_pt = 2 * k2Rinv1GeVf * innerRadius; -// pass = pass and passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f, TightCutFlag); - passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f, TightCutFlag); + pass = pass and passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f, TightCutFlag); +// passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f, TightCutFlag); if(not pass) return pass; @@ -454,7 +454,7 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, } //compute regression radius right here - this computation is expensive!!! -// pass = pass and tempPass; + pass = pass and tempPass; if(not pass) return pass; float xVec[] = {x1, x2, x3, x4, x5}; @@ -470,7 +470,7 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, //extra chi squared cuts! if(regressionRadius < 5.0f/(2.f * k2Rinv1GeVf)) { -// pass = pass and passChiSquaredConstraint(modulesInGPU, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, chiSquared); + pass = pass and passChiSquaredConstraint(modulesInGPU, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, chiSquared); if(not pass) return pass; } @@ -730,14 +730,14 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc if (y4y3 && y3>y2) Py=abs(Py); } - else (moduleType3==1) + else if(moduleType3==1) { if (x3x2 && x2>x1) Px=abs(Px); if (y3y2 && y2>y1) Py=abs(Py); } - + //to get Pz, we use pt/pz=ds/dz, ds is the arclength between MD1 and MD3. float AO=sqrt((x1-x_center)*(x1-x_center)+(y1-y_center)*(y1-y_center)); float BO=sqrt((x_init-x_center)*(x_init-x_center)+(y_init-y_center)*(y_init-y_center)); @@ -865,7 +865,7 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc // if (isnan(rzChiSquared)) printf("rzChi2: %f, residual2: %f, inner_pt:%f, pseudo_phi: %f, charge: %i, Px:%f, Py:%f, x1:%f, x2:%f, x3:%f, x4:%f, x5:%f, y1:%f, y2:%f, y3:%f, y4:%f, y5:%f, z1:%f, z2:%f, z3:%f, z4:%f, z5:%f, x_center:%f, y_center:%f, slope1c:%f, slope3c:%f\n", rzChiSquared, residual_missing, inner_pt, pseudo_phi, charge, Px, Py, x1, x2, x3, x4, x5, y1, y2, y3, y4, y5, z1, z2, z3, z4, z5, x_center, y_center, slope1c, slope3c); -// if(layer1 == 1 and layer2 == 2 and layer3 == 3 and layer4 == 12 and layer5 == 13){ +// if(layer1 == 7 and layer2 == 8 and layer3 == 9 and layer4 == 10 and layer5 == 11 and rzChiSquared>100){ // printf("rt1:%f, rt2:%f, rt3:%f, rt4:%f, rt5:%f\n", rt1, rt2, rt3, rt4, rt5); // printf("x1:%f, x2:%f, x3:%f, x4:%f, x5:%f\n", x1, x2, x3, x4, x5); // printf("y1:%f, y2:%f, y3:%f, y4:%f, y5:%f\n", y1, y2, y3, y4, y5); @@ -874,7 +874,7 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc // printf("z4_ex:%f, z5_ex:%f\n", expectz4, expectz5); // printf("residual_missing:%f\n", residual_missing); // printf("Pt:%f, Px:%f, Py:%f, Pz:%f, charge: %i, residual_missing: %f, residual4: %f, residual5:%f, moduleType3:%i\n", Pt, Px, Py, Pz, charge, residual_missing, residual4, residual5, moduleType3); -// printf("rzChi2: %f, residual2: %f, inner_pt:%f, pseudo_phi: %f, charge: %i, Px:%f, Py:%f, x1:%f, x2:%f, x3:%f, x4:%f, x5:%f, y1:%f, y2:%f, y3:%f, y4:%f, y5:%f, z1:%f, z2:%f, z3:%f, z4:%f, z5:%f, x_center:%f, y_center:%f, slope1c:%f, slope3c:%f\n", rzChiSquared, residual_missing, inner_pt, pseudo_phi, charge, Px, Py, x1, x2, x3, x4, x5, y1, y2, y3, y4, y5, z1, z2, z3, z4, z5, x_center, y_center, slope1c, slope3c); +// if (fabs(rzChiSquared-434.901)<0.01) printf("rzChi2: %f, residual2: %f, residual4: %f, residual5:%f, inner_pt:%f, pseudo_phi: %f, charge: %i, Px:%f, Py:%f, x1:%f, x2:%f, x3:%f, x4:%f, x5:%f, y1:%f, y2:%f, y3:%f, y4:%f, y5:%f, z1:%f, z2:%f, z3:%f, z4:%f, z5:%f, x_center:%f, y_center:%f, slope1c:%f, slope3c:%f\n", rzChiSquared, residual_missing, residual4, residual5, inner_pt, pseudo_phi, charge, Px, Py, x1, x2, x3, x4, x5, y1, y2, y3, y4, y5, z1, z2, z3, z4, z5, x_center, y_center, slope1c, slope3c); // printf("residual_missing:%f\n", residual_missing); // } @@ -885,12 +885,12 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc { if(layer4 == 4 and layer5 == 5) //11 { - if (rzChiSquared < 15.595f) TightCutFlag = 1; - return rzChiSquared < 28.902f; + if (rzChiSquared < 15.627f) TightCutFlag = 1; + return rzChiSquared < 29.035f; } else if(layer4 == 4 and layer5 == 12) //12 { - if (rzChiSquared < 14.614f) TightCutFlag = 1; + if (rzChiSquared < 14.64f) TightCutFlag = 1; return rzChiSquared < 23.037f; } else if(layer4 == 7 and layer5 == 8) //8 @@ -900,63 +900,63 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc } else if(layer4 == 7 and layer5 == 13) //9 { - if (rzChiSquared < 18.085f) TightCutFlag = 1; - return rzChiSquared < 33.023f; + if (rzChiSquared < 18.145f) TightCutFlag = 1; + return rzChiSquared < 33.752f; } else if(layer4 == 12 and layer5 == 13) //10 { - if (rzChiSquared < 13.267f) TightCutFlag = 1; - return rzChiSquared < 21.186f; + if (rzChiSquared < 13.308f) TightCutFlag = 1; + return rzChiSquared < 21.213f; } } else if(layer1 == 1 and layer2 == 2 and layer3 == 7) { if(layer4 == 8 and layer5 == 9) //5 { - if (rzChiSquared < 60.195f) TightCutFlag = 1; - return rzChiSquared < 117.118f; + if (rzChiSquared < 116.148f) TightCutFlag = 1; + return true; } if(layer4 == 8 and layer5 == 14) //6 { - if (rzChiSquared < 19.490f) TightCutFlag = 1; - return rzChiSquared < 55.322f; + if (rzChiSquared < 19.352f) TightCutFlag = 1; + return rzChiSquared < 52.561f; } else if(layer4 == 13 and layer5 == 14) //7 { - if (rzChiSquared < 10.157f) TightCutFlag = 1; - return rzChiSquared < 14.217f; + if (rzChiSquared < 10.392f) TightCutFlag = 1; + return rzChiSquared < 13.76f; } } else if(layer1 == 1 and layer2 == 7 and layer3 == 8 and layer4 == 9) { if (layer5 == 10) //3 { - if (rzChiSquared < 63.697f) TightCutFlag = 1; - return rzChiSquared < 109.584f; + if (rzChiSquared < 111.390f) TightCutFlag = 1; + return true; } if (layer5 == 15) //4 { - if (rzChiSquared < 18.346f) TightCutFlag = 1; - return rzChiSquared < 34.941f; + if (rzChiSquared < 18.351f) TightCutFlag = 1; + return rzChiSquared < 37.941f; } } else if(layer1 == 2 and layer2 == 3 and layer3 == 4) { if(layer4 == 5 and layer5 == 6) //18 { - if (rzChiSquared < 6.053f) TightCutFlag = 1; - return rzChiSquared < 8.629f; + if (rzChiSquared < 6.065f) TightCutFlag = 1; + return rzChiSquared < 8.803f; } else if(layer4 == 5 and layer5 == 12) //19 { if (rzChiSquared < 5.693f) TightCutFlag = 1; - return rzChiSquared < 7.929f; + return rzChiSquared < 7.930f; } else if(layer4 == 12 and layer5 == 13) //20 { - if (rzChiSquared < 5.44f) TightCutFlag = 1; - return rzChiSquared < 7.627f; + if (rzChiSquared < 5.473f) TightCutFlag = 1; + return rzChiSquared < 7.626f; } } else if(layer1 == 2 and layer2 == 3 and layer3 == 7) @@ -968,8 +968,8 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc } if(layer4 == 13 and layer5 == 14) //17 { - if (rzChiSquared < 10.55f) TightCutFlag = 1; - return rzChiSquared < 17.817f; + if (rzChiSquared < 10.772f) TightCutFlag = 1; + return rzChiSquared < 17.945f; } } @@ -977,32 +977,32 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc { if(layer4 == 9 and layer5 == 15) //14 { - if (rzChiSquared < 24.558f) TightCutFlag = 1; - return rzChiSquared < 40.918f; + if (rzChiSquared < 24.662f) TightCutFlag = 1; + return rzChiSquared < 41.036f; } else if(layer4 == 14 and layer5 == 15) //15 { - if (rzChiSquared < 8.752f) TightCutFlag = 1; - return rzChiSquared < 13.678f; + if (rzChiSquared < 8.866f) TightCutFlag = 1; + return rzChiSquared < 14.092f; } } else if(layer1 == 7 and layer2 == 8 and layer3 == 9 and layer4 == 15 and layer5 == 16) //2 { - if (rzChiSquared < 7.994f) TightCutFlag = 1; + if (rzChiSquared < 7.992f) TightCutFlag = 1; return rzChiSquared < 11.622f; } else if(layer1 == 7 and layer2 == 8 and layer3 == 9 and layer4 == 10 and layer5 == 11) //0 { - if (rzChiSquared < 56.313f) TightCutFlag = 1; - return rzChiSquared < 93.893f; + if (rzChiSquared < 94.470f) TightCutFlag = 1; + return true; } else if(layer1 == 7 and layer2 == 8 and layer3 == 9 and layer4 == 10 and layer5 == 16) //1 { - if (rzChiSquared < 22.041f) TightCutFlag = 1; - return rzChiSquared < 37.087f; + if (rzChiSquared < 22.099f) TightCutFlag = 1; + return rzChiSquared < 37.956f; } return true; } From 1ffa5d2e8447a4bc6ed5e880cae64a122b84d36a Mon Sep 17 00:00:00 2001 From: YonsiG Date: Fri, 27 Jan 2023 01:31:31 -0500 Subject: [PATCH 19/25] add linear T5 rz cuts for pT>100 --- SDL/Quintuplet.cu | 26 ++++++++++++++++++++++++++ 1 file changed, 26 insertions(+) diff --git a/SDL/Quintuplet.cu b/SDL/Quintuplet.cu index 6b60f328..784f5420 100644 --- a/SDL/Quintuplet.cu +++ b/SDL/Quintuplet.cu @@ -861,6 +861,32 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc } rzChiSquared += 12*(residual * residual)/(error * error); } + // for set rzchi2 cut +// if (!isnan(rzChiSquared)) rzChiSquared=100; + // if the 5 points are linear, helix calculation gives nan + if (inner_pt>100 || isnan(rzChiSquared)){ + float slope; + if(moduleType1 == 0 and moduleType2 == 0 and moduleType3 == 1) //PSPS2S + { + slope = (z2 -z1)/(rt2 - rt1); + } + else + { + slope = (z3 - z1)/(rt3 - rt1); + } + float residual4_linear = (layer4 <= 6)? ((z4 - z1) - slope * (rt4 - rt1)) : ((rt4 - rt1) - (z4 - z1)/slope); + float residual5_linear = (layer4 <= 6) ? ((z5 - z1) - slope * (rt5 - rt1)) : ((rt5 - rt1) - (z5 - z1)/slope); + + // creating a chi squared type quantity + // 0-> PS, 1->2S + residual4_linear = (moduleType4 == 0) ? residual4_linear/0.15f : residual4_linear/5.0f; + residual5_linear = (moduleType5 == 0) ? residual5_linear/0.15f : residual5_linear/5.0f; + residual4_linear = residual4_linear*100; + residual5_linear = residual5_linear*100; + + rzChiSquared = 12 * (residual4_linear * residual4_linear + residual5_linear * residual5_linear); + return rzChiSquared < 4.677f; + } // rzChiSquared = 12*(residual4 * residual4 + residual5 * residual5 + residual_missing * residual_missing); // if (isnan(rzChiSquared)) printf("rzChi2: %f, residual2: %f, inner_pt:%f, pseudo_phi: %f, charge: %i, Px:%f, Py:%f, x1:%f, x2:%f, x3:%f, x4:%f, x5:%f, y1:%f, y2:%f, y3:%f, y4:%f, y5:%f, z1:%f, z2:%f, z3:%f, z4:%f, z5:%f, x_center:%f, y_center:%f, slope1c:%f, slope3c:%f\n", rzChiSquared, residual_missing, inner_pt, pseudo_phi, charge, Px, Py, x1, x2, x3, x4, x5, y1, y2, y3, y4, y5, z1, z2, z3, z4, z5, x_center, y_center, slope1c, slope3c); From 4eeb31c2402efea49590b6f6dae228f20e7c0d01 Mon Sep 17 00:00:00 2001 From: YonsiG Date: Tue, 25 Apr 2023 10:30:37 -0700 Subject: [PATCH 20/25] bug fix --- SDL/Event.cu | 1 - SDL/Quintuplet.cu | 11 ++--------- SDL/Quintuplet.cuh | 1 - 3 files changed, 2 insertions(+), 11 deletions(-) diff --git a/SDL/Event.cu b/SDL/Event.cu index fb8250c8..c7e1254f 100644 --- a/SDL/Event.cu +++ b/SDL/Event.cu @@ -176,7 +176,6 @@ SDL::Event::~Event() delete[] quintupletsInCPU->residual_missing; delete[] quintupletsInCPU->residual4; delete[] quintupletsInCPU->residual5; -#endif delete quintupletsInCPU; } diff --git a/SDL/Quintuplet.cu b/SDL/Quintuplet.cu index 784f5420..eb43535d 100644 --- a/SDL/Quintuplet.cu +++ b/SDL/Quintuplet.cu @@ -30,7 +30,6 @@ SDL::quintuplets::quintuplets() residual_missing = nullptr; residual4 = nullptr; residual5 = nullptr; -#endif } SDL::quintuplets::~quintuplets() @@ -65,7 +64,6 @@ void SDL::quintuplets::freeMemoryCache() cms::cuda::free_device(dev, residual_missing); cms::cuda::free_device(dev, residual4); cms::cuda::free_device(dev, residual5); -#endif } void SDL::quintuplets::freeMemory(cudaStream_t stream) @@ -94,7 +92,6 @@ void SDL::quintuplets::freeMemory(cudaStream_t stream) cudaFree(residual_missing); cudaFree(residual4); cudaFree(residual5); -#endif cudaStreamSynchronize(stream); } //TODO:Reuse the track candidate one instead of this! @@ -192,7 +189,6 @@ void SDL::createQuintupletsInExplicitMemory(struct SDL::quintuplets& quintuplets quintupletsInGPU.residual_missing = (float*)cms::cuda::allocate_device(dev, nTotalQuintuplets * sizeof(float), stream); quintupletsInGPU.residual4 = (float*)cms::cuda::allocate_device(dev, nTotalQuintuplets * sizeof(float), stream); quintupletsInGPU.residual5 = (float*)cms::cuda::allocate_device(dev, nTotalQuintuplets * sizeof(float), stream); -#endif #else cudaMalloc(&quintupletsInGPU.tripletIndices, 2 * nTotalQuintuplets * sizeof(unsigned int)); cudaMalloc(&quintupletsInGPU.lowerModuleIndices, 5 * nTotalQuintuplets * sizeof(uint16_t)); @@ -754,7 +750,6 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc float zsi, rtsi; int layeri, moduleTypei; - float expectrt4=0,expectrt5=0,expectz4=0, expectz5=0; rzChiSquared=0; for(size_t i = 2; i < 6; i++) { @@ -799,8 +794,6 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc float x = x_init + Px/a*sin(rou*s)-Py/a*(1-cos(rou*s)); float y = y_init + Py/a*sin(rou*s)+Px/a*(1-cos(rou*s)); diffr = (rtsi-sqrt(x*x+y*y))*100; - if (i==4) expectrt4=sqrt(x*x+y*y); - if (i==5) expectrt5=sqrt(x*x+y*y); // for barrel if (layeri<=6) @@ -1576,7 +1569,7 @@ __global__ void SDL::createQuintupletsInGPUv2(struct SDL::modules& modulesInGPU, uint16_t lowerModule4 = tripletsInGPU.lowerModuleIndices[3 * outerTripletIndex + 1]; uint16_t lowerModule5 = tripletsInGPU.lowerModuleIndices[3 * outerTripletIndex + 2]; - float innerRadius, outerRadius, bridgeRadius, regressionG, regressionF, regressionRadius, rzChiSquared, chiSquared, nonAnchorChiSquared, residual_missing, residual4, residual; //required for making distributions + float innerRadius, outerRadius, bridgeRadius, regressionG, regressionF, regressionRadius, rzChiSquared, chiSquared, nonAnchorChiSquared, residual_missing, residual4, residual5; //required for making distributions bool TightCutFlag; bool success = runQuintupletDefaultAlgo(modulesInGPU, mdsInGPU, segmentsInGPU, tripletsInGPU, lowerModule1, lowerModule2, lowerModule3, lowerModule4, lowerModule5, innerTripletIndex, outerTripletIndex, innerRadius, outerRadius, bridgeRadius, regressionG, regressionF, regressionRadius, rzChiSquared, chiSquared, nonAnchorChiSquared, residual_missing, residual4, residual5, TightCutFlag); @@ -1618,7 +1611,7 @@ __global__ void SDL::createQuintupletsInGPUv2(struct SDL::modules& modulesInGPU, float eta = mdsInGPU.anchorEta[segmentsInGPU.mdIndices[2*tripletsInGPU.segmentIndices[2*innerTripletIndex+layer2_adjustment]]]; float pt = (innerRadius+outerRadius)*3.8f*1.602f/(2*100*5.39f); float scores = chiSquared + nonAnchorChiSquared; - addQuintupletToMemory(tripletsInGPU, quintupletsInGPU, innerTripletIndex, outerTripletIndex, lowerModule1, lowerModule2, lowerModule3, lowerModule4, lowerModule5, innerRadius, bridgeRadius, outerRadius, regressionG, regressionF, regressionRadius, rzChiSquared, chiSquared, nonAnchorChiSquared, rzChiSquared, residual_missing, residual4, residual5, pt,eta,phi,scores,layer,quintupletIndex, TightCutFlag); + addQuintupletToMemory(tripletsInGPU, quintupletsInGPU, innerTripletIndex, outerTripletIndex, lowerModule1, lowerModule2, lowerModule3, lowerModule4, lowerModule5, innerRadius, bridgeRadius, outerRadius, regressionG, regressionF, regressionRadius, rzChiSquared, chiSquared, nonAnchorChiSquared, residual_missing, residual4, residual5, pt,eta,phi,scores,layer,quintupletIndex, TightCutFlag); tripletsInGPU.partOfT5[quintupletsInGPU.tripletIndices[2 * quintupletIndex]] = true; tripletsInGPU.partOfT5[quintupletsInGPU.tripletIndices[2 * quintupletIndex + 1]] = true; diff --git a/SDL/Quintuplet.cuh b/SDL/Quintuplet.cuh index 84e26145..cb79ba55 100644 --- a/SDL/Quintuplet.cuh +++ b/SDL/Quintuplet.cuh @@ -55,7 +55,6 @@ namespace SDL float* residual_missing; float* residual4; float* residual5; -#endif quintuplets(); ~quintuplets(); From f37d2504c9a6db2bedbf4e1996830305261f878c Mon Sep 17 00:00:00 2001 From: YonsiG Date: Sun, 7 May 2023 07:10:29 -0700 Subject: [PATCH 21/25] remove comments and redundant variables --- SDL/Event.cu | 6 ---- SDL/Quintuplet.cu | 78 ++++++++++------------------------------------ SDL/Quintuplet.cuh | 9 ++---- 3 files changed, 19 insertions(+), 74 deletions(-) diff --git a/SDL/Event.cu b/SDL/Event.cu index c7e1254f..9a3075b1 100644 --- a/SDL/Event.cu +++ b/SDL/Event.cu @@ -173,9 +173,6 @@ SDL::Event::~Event() delete[] quintupletsInCPU->chiSquared; delete[] quintupletsInCPU->rzChiSquared; delete[] quintupletsInCPU->nonAnchorChiSquared; - delete[] quintupletsInCPU->residual_missing; - delete[] quintupletsInCPU->residual4; - delete[] quintupletsInCPU->residual5; delete quintupletsInCPU; } @@ -1965,9 +1962,6 @@ SDL::quintuplets* SDL::Event::getQuintuplets() quintupletsInCPU->chiSquared = new float[nMemoryLocations]; quintupletsInCPU->nonAnchorChiSquared = new float[nMemoryLocations]; quintupletsInCPU->rzChiSquared = new float[nMemoryLocations]; - quintupletsInCPU->residual_missing = new float[nMemoryLocations]; - quintupletsInCPU->residual4 = new float[nMemoryLocations]; - quintupletsInCPU->residual5 = new float[nMemoryLocations]; cudaMemcpyAsync(quintupletsInCPU->nQuintuplets, quintupletsInGPU->nQuintuplets, nLowerModules * sizeof(unsigned int), cudaMemcpyDeviceToHost,stream); cudaMemcpyAsync(quintupletsInCPU->totOccupancyQuintuplets, quintupletsInGPU->totOccupancyQuintuplets, nLowerModules * sizeof(unsigned int), cudaMemcpyDeviceToHost,stream); diff --git a/SDL/Quintuplet.cu b/SDL/Quintuplet.cu index eb43535d..54da08af 100644 --- a/SDL/Quintuplet.cu +++ b/SDL/Quintuplet.cu @@ -27,9 +27,6 @@ SDL::quintuplets::quintuplets() chiSquared = nullptr; rzChiSquared = nullptr; nonAnchorChiSquared = nullptr; - residual_missing = nullptr; - residual4 = nullptr; - residual5 = nullptr; } SDL::quintuplets::~quintuplets() @@ -61,9 +58,6 @@ void SDL::quintuplets::freeMemoryCache() cms::cuda::free_device(dev, rzChiSquared); cms::cuda::free_device(dev, chiSquared); cms::cuda::free_device(dev, nonAnchorChiSquared); - cms::cuda::free_device(dev, residual_missing); - cms::cuda::free_device(dev, residual4); - cms::cuda::free_device(dev, residual5); } void SDL::quintuplets::freeMemory(cudaStream_t stream) @@ -89,10 +83,7 @@ void SDL::quintuplets::freeMemory(cudaStream_t stream) cudaFree(rzChiSquared); cudaFree(chiSquared); cudaFree(nonAnchorChiSquared); - cudaFree(residual_missing); - cudaFree(residual4); - cudaFree(residual5); -cudaStreamSynchronize(stream); + cudaStreamSynchronize(stream); } //TODO:Reuse the track candidate one instead of this! __global__ void SDL::createEligibleModulesListForQuintupletsGPU(struct modules& modulesInGPU,struct triplets& tripletsInGPU, struct objectRanges& rangesInGPU) @@ -186,9 +177,6 @@ void SDL::createQuintupletsInExplicitMemory(struct SDL::quintuplets& quintuplets quintupletsInGPU.rzChiSquared = (float*)cms::cuda::allocate_device(dev, nTotalQuintuplets * sizeof(float), stream); quintupletsInGPU.chiSquared = (float*)cms::cuda::allocate_device(dev, nTotalQuintuplets * sizeof(float), stream); quintupletsInGPU.nonAnchorChiSquared = (float*)cms::cuda::allocate_device(dev, nTotalQuintuplets * sizeof(float), stream); - quintupletsInGPU.residual_missing = (float*)cms::cuda::allocate_device(dev, nTotalQuintuplets * sizeof(float), stream); - quintupletsInGPU.residual4 = (float*)cms::cuda::allocate_device(dev, nTotalQuintuplets * sizeof(float), stream); - quintupletsInGPU.residual5 = (float*)cms::cuda::allocate_device(dev, nTotalQuintuplets * sizeof(float), stream); #else cudaMalloc(&quintupletsInGPU.tripletIndices, 2 * nTotalQuintuplets * sizeof(unsigned int)); cudaMalloc(&quintupletsInGPU.lowerModuleIndices, 5 * nTotalQuintuplets * sizeof(uint16_t)); @@ -211,9 +199,6 @@ void SDL::createQuintupletsInExplicitMemory(struct SDL::quintuplets& quintuplets cudaMalloc(&quintupletsInGPU.rzChiSquared, nTotalQuintuplets * sizeof(float)); cudaMalloc(&quintupletsInGPU.chiSquared, nTotalQuintuplets * sizeof(float)); cudaMalloc(&quintupletsInGPU.nonAnchorChiSquared, nTotalQuintuplets * sizeof(float)); - cudaMalloc(&quintupletsInGPU.residual_missing, nTotalQuintuplets * sizeof(float)); - cudaMalloc(&quintupletsInGPU.residual4, nTotalQuintuplets * sizeof(float)); - cudaMalloc(&quintupletsInGPU.residual5, nTotalQuintuplets * sizeof(float)); cudaMalloc(&quintupletsInGPU.nMemoryLocations, sizeof(unsigned int)); #endif cudaMemsetAsync(quintupletsInGPU.nQuintuplets,0,nLowerModules * sizeof(unsigned int),stream); @@ -229,7 +214,7 @@ void SDL::createQuintupletsInExplicitMemory(struct SDL::quintuplets& quintuplets __device__ void SDL::addQuintupletToMemory(struct SDL::triplets& tripletsInGPU, struct SDL::quintuplets& quintupletsInGPU, unsigned int innerTripletIndex, unsigned int outerTripletIndex, uint16_t& lowerModule1, uint16_t& lowerModule2, uint16_t& lowerModule3, uint16_t& lowerModule4, uint16_t& lowerModule5, float& innerRadius, float& bridgeRadius, float& outerRadius, float& regressionG, float& regressionF, float& regressionRadius, float& rzChiSquared, float& rPhiChiSquared, float& - nonAnchorChiSquared, float residual_missing, float residual4, float residual5, float pt, float eta, float phi, float scores, uint8_t layer, unsigned int quintupletIndex, bool TightCutFlag) + nonAnchorChiSquared, float pt, float eta, float phi, float scores, uint8_t layer, unsigned int quintupletIndex, bool TightCutFlag) { quintupletsInGPU.tripletIndices[2 * quintupletIndex] = innerTripletIndex; @@ -272,12 +257,9 @@ __device__ void SDL::addQuintupletToMemory(struct SDL::triplets& tripletsInGPU, quintupletsInGPU.rzChiSquared[quintupletIndex] = rzChiSquared; quintupletsInGPU.chiSquared[quintupletIndex] = rPhiChiSquared; quintupletsInGPU.nonAnchorChiSquared[quintupletIndex] = nonAnchorChiSquared; - quintupletsInGPU.residual_missing[quintupletIndex] = residual_missing; - quintupletsInGPU.residual4[quintupletIndex] = residual4; - quintupletsInGPU.residual5[quintupletIndex] = residual5; } -__device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, struct SDL::segments& segmentsInGPU, struct SDL::triplets& tripletsInGPU, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, unsigned int& innerTripletIndex, unsigned int& outerTripletIndex, float& innerRadius, float& outerRadius, float& bridgeRadius, float& regressionG, float& regressionF, float& regressionRadius, float& rzChiSquared, float& chiSquared, float& nonAnchorChiSquared, float& residual_missing, float& residual4, float& residual5, bool& TightCutFlag) +__device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, struct SDL::segments& segmentsInGPU, struct SDL::triplets& tripletsInGPU, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, unsigned int& innerTripletIndex, unsigned int& outerTripletIndex, float& innerRadius, float& outerRadius, float& bridgeRadius, float& regressionG, float& regressionF, float& regressionRadius, float& rzChiSquared, float& chiSquared, float& nonAnchorChiSquared, bool& TightCutFlag) { bool pass = true; unsigned int firstSegmentIndex = tripletsInGPU.segmentIndices[2 * innerTripletIndex]; @@ -307,12 +289,6 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, pass = pass and runQuintupletDefaultAlgo(modulesInGPU, mdsInGPU, segmentsInGPU, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex4, lowerModuleIndex5, firstSegmentIndex, fourthSegmentIndex, firstMDIndex, secondMDIndex, fourthMDIndex, fifthMDIndex, zOut, rtOut, deltaPhiPos, deltaPhi, betaIn, betaOut, pt_beta, zLo, zHi, rtLo, rtHi, zLoPointed, zHiPointed, sdlCut, betaInCut, betaOutCut, deltaBetaCut, kZ); if(not pass) return pass; -// pass = pass and runQuintupletDefaultAlgo(modulesInGPU, mdsInGPU, segmentsInGPU, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex3, lowerModuleIndex4, secondSegmentIndex, thirdSegmentIndex, secondMDIndex, thirdMDIndex, thirdMDIndex, fourthMDIndex, zOut, rtOut, deltaPhiPos, deltaPhi, betaIn, betaOut, pt_beta, zLo, zHi, rtLo, rtHi, zLoPointed, zHiPointed, sdlCut, betaInCut, betaOutCut, deltaBetaCut, kZ); - if(not pass) return pass; - -// pass = pass and runQuintupletDefaultAlgo(modulesInGPU, mdsInGPU, segmentsInGPU, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, secondSegmentIndex, fourthSegmentIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, zOut, rtOut, deltaPhiPos, deltaPhi, betaIn, betaOut, pt_beta, zLo, zHi, rtLo, rtHi, zLoPointed, zHiPointed, sdlCut, betaInCut, betaOutCut, deltaBetaCut, kZ); - if(not pass) return pass; - float x1 = mdsInGPU.anchorX[firstMDIndex]; float x2 = mdsInGPU.anchorX[secondMDIndex]; float x3 = mdsInGPU.anchorX[thirdMDIndex]; @@ -401,8 +377,7 @@ __device__ bool SDL::runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, innerRadius = computeRadiusFromThreeAnchorHits(x1, y1, x2, y2, x3, y3, g, f); float inner_pt = 2 * k2Rinv1GeVf * innerRadius; - pass = pass and passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f, TightCutFlag); -// passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, residual_missing, residual4, residual5, inner_pt, innerRadius, g, f, TightCutFlag); + pass = pass and passT5RZConstraint(modulesInGPU, mdsInGPU, firstMDIndex, secondMDIndex, thirdMDIndex, fourthMDIndex, fifthMDIndex, lowerModuleIndex1, lowerModuleIndex2, lowerModuleIndex3, lowerModuleIndex4, lowerModuleIndex5, rzChiSquared, inner_pt, innerRadius, g, f, TightCutFlag); if(not pass) return pass; @@ -612,7 +587,7 @@ __device__ bool SDL::passChiSquaredConstraint(struct SDL::modules& modulesInGPU, } //bounds can be found at http://uaf-10.t2.ucsd.edu/~bsathian/SDL/T5_RZFix/t5_rz_thresholds.txt -__device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex, unsigned int fifthMDIndex, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, float& rzChiSquared, float& residual_missing, float& residual4, float& residual5, float inner_pt, float innerRadius, float g, float f, bool& TightCutFlag) +__device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex, unsigned int fifthMDIndex, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, float& rzChiSquared, float inner_pt, float innerRadius, float g, float f, bool& TightCutFlag) { //(g,f) is the center of the circle fitted by the innermost 3 points on x,y coordinates const float& rt1 = mdsInGPU.anchorRt[firstMDIndex]/100; //in the unit of m instead of cm @@ -708,7 +683,10 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc } if (x_inity_center) //2nd quad { - if (charge==-1) {Px=-Px; Py=-Py;} + if (charge==-1) { + Px=-Px; + Py=-Py; + } } if (x_initx_center && y_initsmallerPz->smaller residual float Pz=(z_init-z1)/ds*Pt; float p = sqrt(Px*Px+Py*Py+Pz*Pz); @@ -843,19 +823,16 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc } if (i==2 || i==3){ residual = (layeri <= 6 && ((side == SDL::Center) or (drdz < 1))) ? diffz : diffr; -// residual_missing=residual; float projection_missing=1; if (drdz<1) projection_missing = ((subdets == SDL::Endcap) or (side == SDL::Center)) ? 1.f : 1/sqrt(1+drdz*drdz); // cos(atan(drdz)), if dr/dz<1 if (drdz>1) projection_missing = ((subdets == SDL::Endcap) or (side == SDL::Center)) ? 1.f : drdz/sqrt(1+drdz*drdz);//sin(atan(drdz)), if dr/dz>1 error=error*projection_missing; - residual_missing=residual/error; } rzChiSquared += 12*(residual * residual)/(error * error); } // for set rzchi2 cut -// if (!isnan(rzChiSquared)) rzChiSquared=100; // if the 5 points are linear, helix calculation gives nan if (inner_pt>100 || isnan(rzChiSquared)){ float slope; @@ -880,22 +857,6 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc rzChiSquared = 12 * (residual4_linear * residual4_linear + residual5_linear * residual5_linear); return rzChiSquared < 4.677f; } -// rzChiSquared = 12*(residual4 * residual4 + residual5 * residual5 + residual_missing * residual_missing); - -// if (isnan(rzChiSquared)) printf("rzChi2: %f, residual2: %f, inner_pt:%f, pseudo_phi: %f, charge: %i, Px:%f, Py:%f, x1:%f, x2:%f, x3:%f, x4:%f, x5:%f, y1:%f, y2:%f, y3:%f, y4:%f, y5:%f, z1:%f, z2:%f, z3:%f, z4:%f, z5:%f, x_center:%f, y_center:%f, slope1c:%f, slope3c:%f\n", rzChiSquared, residual_missing, inner_pt, pseudo_phi, charge, Px, Py, x1, x2, x3, x4, x5, y1, y2, y3, y4, y5, z1, z2, z3, z4, z5, x_center, y_center, slope1c, slope3c); - -// if(layer1 == 7 and layer2 == 8 and layer3 == 9 and layer4 == 10 and layer5 == 11 and rzChiSquared>100){ -// printf("rt1:%f, rt2:%f, rt3:%f, rt4:%f, rt5:%f\n", rt1, rt2, rt3, rt4, rt5); -// printf("x1:%f, x2:%f, x3:%f, x4:%f, x5:%f\n", x1, x2, x3, x4, x5); -// printf("y1:%f, y2:%f, y3:%f, y4:%f, y5:%f\n", y1, y2, y3, y4, y5); -// printf("z1:%f, z2:%f, z3:%f, z4:%f, z5:%f\n", z1, z2, z3, z4, z5); -// printf("rt4_ex:%f, rt5_ex:%f\n", expectrt4, expectrt5); -// printf("z4_ex:%f, z5_ex:%f\n", expectz4, expectz5); -// printf("residual_missing:%f\n", residual_missing); -// printf("Pt:%f, Px:%f, Py:%f, Pz:%f, charge: %i, residual_missing: %f, residual4: %f, residual5:%f, moduleType3:%i\n", Pt, Px, Py, Pz, charge, residual_missing, residual4, residual5, moduleType3); -// if (fabs(rzChiSquared-434.901)<0.01) printf("rzChi2: %f, residual2: %f, residual4: %f, residual5:%f, inner_pt:%f, pseudo_phi: %f, charge: %i, Px:%f, Py:%f, x1:%f, x2:%f, x3:%f, x4:%f, x5:%f, y1:%f, y2:%f, y3:%f, y4:%f, y5:%f, z1:%f, z2:%f, z3:%f, z4:%f, z5:%f, x_center:%f, y_center:%f, slope1c:%f, slope3c:%f\n", rzChiSquared, residual_missing, residual4, residual5, inner_pt, pseudo_phi, charge, Px, Py, x1, x2, x3, x4, x5, y1, y2, y3, y4, y5, z1, z2, z3, z4, z5, x_center, y_center, slope1c, slope3c); -// printf("residual_missing:%f\n", residual_missing); -// } // when building T5, apply 99% chi2 cuts as default, and add to pT5 collection. But when adding T5 to TC collections, appy 95% cut to reduce the fake rate TightCutFlag = false; @@ -1280,13 +1241,6 @@ __device__ float SDL::computeRadiusFromThreeAnchorHits(float x1, float y1, float //(g,f) -> center //first anchor hit - (x1,y1), second anchor hit - (x2,y2), third anchor hit - (x3, y3) - /* - if((y1 - y3) * (x2 - x3) - (x1 - x3) * (y2 - y3) == 0) - { - return -1; //WTF man, three collinear points! - } - */ - float denomInv = 1.0f/((y1 - y3) * (x2 - x3) - (x1 - x3) * (y2 - y3)); float xy1sqr = x1 * x1 + y1 * y1; @@ -1569,9 +1523,9 @@ __global__ void SDL::createQuintupletsInGPUv2(struct SDL::modules& modulesInGPU, uint16_t lowerModule4 = tripletsInGPU.lowerModuleIndices[3 * outerTripletIndex + 1]; uint16_t lowerModule5 = tripletsInGPU.lowerModuleIndices[3 * outerTripletIndex + 2]; - float innerRadius, outerRadius, bridgeRadius, regressionG, regressionF, regressionRadius, rzChiSquared, chiSquared, nonAnchorChiSquared, residual_missing, residual4, residual5; //required for making distributions + float innerRadius, outerRadius, bridgeRadius, regressionG, regressionF, regressionRadius, rzChiSquared, chiSquared, nonAnchorChiSquared; //required for making distributions bool TightCutFlag; - bool success = runQuintupletDefaultAlgo(modulesInGPU, mdsInGPU, segmentsInGPU, tripletsInGPU, lowerModule1, lowerModule2, lowerModule3, lowerModule4, lowerModule5, innerTripletIndex, outerTripletIndex, innerRadius, outerRadius, bridgeRadius, regressionG, regressionF, regressionRadius, rzChiSquared, chiSquared, nonAnchorChiSquared, residual_missing, residual4, residual5, TightCutFlag); + bool success = runQuintupletDefaultAlgo(modulesInGPU, mdsInGPU, segmentsInGPU, tripletsInGPU, lowerModule1, lowerModule2, lowerModule3, lowerModule4, lowerModule5, innerTripletIndex, outerTripletIndex, innerRadius, outerRadius, bridgeRadius, regressionG, regressionF, regressionRadius, rzChiSquared, chiSquared, nonAnchorChiSquared, TightCutFlag); if(success) { @@ -1611,7 +1565,7 @@ __global__ void SDL::createQuintupletsInGPUv2(struct SDL::modules& modulesInGPU, float eta = mdsInGPU.anchorEta[segmentsInGPU.mdIndices[2*tripletsInGPU.segmentIndices[2*innerTripletIndex+layer2_adjustment]]]; float pt = (innerRadius+outerRadius)*3.8f*1.602f/(2*100*5.39f); float scores = chiSquared + nonAnchorChiSquared; - addQuintupletToMemory(tripletsInGPU, quintupletsInGPU, innerTripletIndex, outerTripletIndex, lowerModule1, lowerModule2, lowerModule3, lowerModule4, lowerModule5, innerRadius, bridgeRadius, outerRadius, regressionG, regressionF, regressionRadius, rzChiSquared, chiSquared, nonAnchorChiSquared, residual_missing, residual4, residual5, pt,eta,phi,scores,layer,quintupletIndex, TightCutFlag); + addQuintupletToMemory(tripletsInGPU, quintupletsInGPU, innerTripletIndex, outerTripletIndex, lowerModule1, lowerModule2, lowerModule3, lowerModule4, lowerModule5, innerRadius, bridgeRadius, outerRadius, regressionG, regressionF, regressionRadius, rzChiSquared, chiSquared, nonAnchorChiSquared, pt,eta,phi,scores,layer,quintupletIndex, TightCutFlag); tripletsInGPU.partOfT5[quintupletsInGPU.tripletIndices[2 * quintupletIndex]] = true; tripletsInGPU.partOfT5[quintupletsInGPU.tripletIndices[2 * quintupletIndex + 1]] = true; diff --git a/SDL/Quintuplet.cuh b/SDL/Quintuplet.cuh index cb79ba55..1da92a4c 100644 --- a/SDL/Quintuplet.cuh +++ b/SDL/Quintuplet.cuh @@ -52,9 +52,6 @@ namespace SDL float* rzChiSquared; float* chiSquared; float* nonAnchorChiSquared; - float* residual_missing; - float* residual4; - float* residual5; quintuplets(); ~quintuplets(); @@ -69,12 +66,12 @@ namespace SDL // CUDA_DEV void rmQuintupletToMemory(struct SDL::quintuplets& quintupletsInGPU, unsigned int quintupletIndex); - CUDA_DEV void addQuintupletToMemory(struct SDL::triplets& tripletsInGPU, struct SDL::quintuplets& quintupletsInGPU, unsigned int innerTripletIndex, unsigned int outerTripletIndex, uint16_t& lowerModule1, uint16_t& lowerModule2, uint16_t& lowerModule3, uint16_t& lowerModule4, uint16_t& lowerModule5, float& innerRadius, float& bridgeRadius, float& outerRadius, float& regressionG, float& regressionF, float& regressionRadius, float& rzChiSquared, float& rPhiChiSquared, float& nonAnchorChiSquared, float residual_missing, float residual4, float residual5, float pt, float eta, float phi, float scores, uint8_t layer, unsigned int quintupletIndex, bool TightCutFlag); + CUDA_DEV void addQuintupletToMemory(struct SDL::triplets& tripletsInGPU, struct SDL::quintuplets& quintupletsInGPU, unsigned int innerTripletIndex, unsigned int outerTripletIndex, uint16_t& lowerModule1, uint16_t& lowerModule2, uint16_t& lowerModule3, uint16_t& lowerModule4, uint16_t& lowerModule5, float& innerRadius, float& bridgeRadius, float& outerRadius, float& regressionG, float& regressionF, float& regressionRadius, float& rzChiSquared, float& rPhiChiSquared, float& nonAnchorChiSquared, float pt, float eta, float phi, float scores, uint8_t layer, unsigned int quintupletIndex, bool TightCutFlag); - CUDA_DEV bool runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, struct SDL::segments& segmentsInGPU, struct SDL::triplets& tripletsInGPU, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, unsigned int& innerTripletIndex, unsigned int& outerTripletIndex, float& innerRadius, float& outerRadius, float& bridgeRadius, float& regressionG, float& regressionF, float& regressionRadius, float& rzChiSquared, float& chiSquared, float& nonAnchorChiSquared, float& residual_missing, float& residual4, float& residual5, bool& TightCutFlag); + CUDA_DEV bool runQuintupletDefaultAlgo(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, struct SDL::segments& segmentsInGPU, struct SDL::triplets& tripletsInGPU, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, unsigned int& innerTripletIndex, unsigned int& outerTripletIndex, float& innerRadius, float& outerRadius, float& bridgeRadius, float& regressionG, float& regressionF, float& regressionRadius, float& rzChiSquared, float& chiSquared, float& nonAnchorChiSquared, bool& TightCutFlag); - CUDA_DEV bool passT5RZConstraint(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex, unsigned int fifthMDIndex, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, float& rzChiSquared, float& residual_missing, float& residual4, float& residual5, float inner_pt, float innerRadius, float g, float f, bool& TightCutFlag); + CUDA_DEV bool passT5RZConstraint(struct SDL::modules& modulesInGPU, struct SDL::miniDoublets& mdsInGPU, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex, unsigned int fifthMDIndex, uint16_t& lowerModuleIndex1, uint16_t& lowerModuleIndex2, uint16_t& lowerModuleIndex3, uint16_t& lowerModuleIndex4, uint16_t& lowerModuleIndex5, float& rzChiSquared, float inner_pt, float innerRadius, float g, float f, bool& TightCutFlag); CUDA_DEV bool T5HasCommonMiniDoublet(struct triplets& tripletsInGPU, struct segments& segmentsInGPU, unsigned int innerTripletIndex, unsigned int outerTripletIndex); From be78b4e499a99d7738747650000357f24a8a06b6 Mon Sep 17 00:00:00 2001 From: YonsiG Date: Mon, 8 May 2023 04:54:09 -0700 Subject: [PATCH 22/25] remove the residual4 and 5 variable --- SDL/Quintuplet.cu | 2 -- 1 file changed, 2 deletions(-) diff --git a/SDL/Quintuplet.cu b/SDL/Quintuplet.cu index 54da08af..4e69e0c9 100644 --- a/SDL/Quintuplet.cu +++ b/SDL/Quintuplet.cu @@ -805,8 +805,6 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc { error = 5.0f; } - if (i==4) residual4=residual/error; - if (i==5) residual5=residual/error; //check the tilted module, side: PosZ, NegZ, Center(for not tilted) float drdz; From 0edef6d3f8307bf311555f8789cfb3453cb305e6 Mon Sep 17 00:00:00 2001 From: YonsiG Date: Mon, 8 May 2023 05:26:17 -0700 Subject: [PATCH 23/25] reordering of the regions --- SDL/Quintuplet.cu | 136 ++++++++++++++++++++++++---------------------- 1 file changed, 70 insertions(+), 66 deletions(-) diff --git a/SDL/Quintuplet.cu b/SDL/Quintuplet.cu index 4e69e0c9..4d6e3761 100644 --- a/SDL/Quintuplet.cu +++ b/SDL/Quintuplet.cu @@ -859,34 +859,41 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc // when building T5, apply 99% chi2 cuts as default, and add to pT5 collection. But when adding T5 to TC collections, appy 95% cut to reduce the fake rate TightCutFlag = false; //categories! - if(layer1 == 1 and layer2 == 2 and layer3 == 3) + // The category numbers are related to module regions and layers, decoding of the region numbers can be found here in slide 2 table. https://github.com/SegmentLinking/TrackLooper/files/11420927/part.2.pdf + // The commented numbers after each case is the region code, and can look it up from the table to see which category it belongs to. For example, //0 means T5 built with Endcap 1,2,3,4,5 ps modules + + if(layer1 == 7 and layer2 == 8 and layer3 == 9 and layer4 == 10 and layer5 == 11) //0 { - if(layer4 == 4 and layer5 == 5) //11 - { - if (rzChiSquared < 15.627f) TightCutFlag = 1; - return rzChiSquared < 29.035f; - } - else if(layer4 == 4 and layer5 == 12) //12 - { - if (rzChiSquared < 14.64f) TightCutFlag = 1; - return rzChiSquared < 23.037f; - } - else if(layer4 == 7 and layer5 == 8) //8 - { - if (rzChiSquared < 27.824f) TightCutFlag = 1; - return rzChiSquared < 44.247f; - } - else if(layer4 == 7 and layer5 == 13) //9 + if (rzChiSquared < 94.470f) TightCutFlag = 1; + return true; + } + + else if(layer1 == 7 and layer2 == 8 and layer3 == 9 and layer4 == 10 and layer5 == 16) //1 + { + if (rzChiSquared < 22.099f) TightCutFlag = 1; + return rzChiSquared < 37.956f; + } + + else if(layer1 == 7 and layer2 == 8 and layer3 == 9 and layer4 == 15 and layer5 == 16) //2 + { + if (rzChiSquared < 7.992f) TightCutFlag = 1; + return rzChiSquared < 11.622f; + } + + else if(layer1 == 1 and layer2 == 7 and layer3 == 8 and layer4 == 9) + { + if (layer5 == 10) //3 { - if (rzChiSquared < 18.145f) TightCutFlag = 1; - return rzChiSquared < 33.752f; + if (rzChiSquared < 111.390f) TightCutFlag = 1; + return true; } - else if(layer4 == 12 and layer5 == 13) //10 + if (layer5 == 15) //4 { - if (rzChiSquared < 13.308f) TightCutFlag = 1; - return rzChiSquared < 21.213f; + if (rzChiSquared < 18.351f) TightCutFlag = 1; + return rzChiSquared < 37.941f; } } + else if(layer1 == 1 and layer2 == 2 and layer3 == 7) { if(layer4 == 8 and layer5 == 9) //5 @@ -905,36 +912,46 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc return rzChiSquared < 13.76f; } } - else if(layer1 == 1 and layer2 == 7 and layer3 == 8 and layer4 == 9) + else if(layer1 == 1 and layer2 == 2 and layer3 == 3) { - if (layer5 == 10) //3 + if(layer4 == 7 and layer5 == 8) //8 + { + if (rzChiSquared < 27.824f) TightCutFlag = 1; + return rzChiSquared < 44.247f; + } + else if(layer4 == 7 and layer5 == 13) //9 { - if (rzChiSquared < 111.390f) TightCutFlag = 1; - return true; + if (rzChiSquared < 18.145f) TightCutFlag = 1; + return rzChiSquared < 33.752f; } - if (layer5 == 15) //4 + else if(layer4 == 12 and layer5 == 13) //10 { - if (rzChiSquared < 18.351f) TightCutFlag = 1; - return rzChiSquared < 37.941f; + if (rzChiSquared < 13.308f) TightCutFlag = 1; + return rzChiSquared < 21.213f; } - } - else if(layer1 == 2 and layer2 == 3 and layer3 == 4) - { - if(layer4 == 5 and layer5 == 6) //18 + else if(layer4 == 4 and layer5 == 5) //11 { - if (rzChiSquared < 6.065f) TightCutFlag = 1; - return rzChiSquared < 8.803f; + if (rzChiSquared < 15.627f) TightCutFlag = 1; + return rzChiSquared < 29.035f; } - else if(layer4 == 5 and layer5 == 12) //19 + else if(layer4 == 4 and layer5 == 12) //12 { - if (rzChiSquared < 5.693f) TightCutFlag = 1; - return rzChiSquared < 7.930f; + if (rzChiSquared < 14.64f) TightCutFlag = 1; + return rzChiSquared < 23.037f; } + } - else if(layer4 == 12 and layer5 == 13) //20 + else if(layer1 == 2 and layer2 == 7 and layer3 == 8) + { + if(layer4 == 9 and layer5 == 15) //14 { - if (rzChiSquared < 5.473f) TightCutFlag = 1; - return rzChiSquared < 7.626f; + if (rzChiSquared < 24.662f) TightCutFlag = 1; + return rzChiSquared < 41.036f; + } + else if(layer4 == 14 and layer5 == 15) //15 + { + if (rzChiSquared < 8.866f) TightCutFlag = 1; + return rzChiSquared < 14.092f; } } else if(layer1 == 2 and layer2 == 3 and layer3 == 7) @@ -950,37 +967,24 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc return rzChiSquared < 17.945f; } } - - else if(layer1 == 2 and layer2 == 7 and layer3 == 8) + else if(layer1 == 2 and layer2 == 3 and layer3 == 4) { - if(layer4 == 9 and layer5 == 15) //14 + if(layer4 == 5 and layer5 == 6) //18 { - if (rzChiSquared < 24.662f) TightCutFlag = 1; - return rzChiSquared < 41.036f; + if (rzChiSquared < 6.065f) TightCutFlag = 1; + return rzChiSquared < 8.803f; } - else if(layer4 == 14 and layer5 == 15) //15 + else if(layer4 == 5 and layer5 == 12) //19 { - if (rzChiSquared < 8.866f) TightCutFlag = 1; - return rzChiSquared < 14.092f; + if (rzChiSquared < 5.693f) TightCutFlag = 1; + return rzChiSquared < 7.930f; } - } - else if(layer1 == 7 and layer2 == 8 and layer3 == 9 and layer4 == 15 and layer5 == 16) //2 - { - if (rzChiSquared < 7.992f) TightCutFlag = 1; - return rzChiSquared < 11.622f; - } - - else if(layer1 == 7 and layer2 == 8 and layer3 == 9 and layer4 == 10 and layer5 == 11) //0 - { - if (rzChiSquared < 94.470f) TightCutFlag = 1; - return true; - } - - else if(layer1 == 7 and layer2 == 8 and layer3 == 9 and layer4 == 10 and layer5 == 16) //1 - { - if (rzChiSquared < 22.099f) TightCutFlag = 1; - return rzChiSquared < 37.956f; + else if(layer4 == 12 and layer5 == 13) //20 + { + if (rzChiSquared < 5.473f) TightCutFlag = 1; + return rzChiSquared < 7.626f; + } } return true; } From 9c9d618d30d78e529168672c35476873b8dd6c6c Mon Sep 17 00:00:00 2001 From: YonsiG Date: Tue, 9 May 2023 11:55:09 -0700 Subject: [PATCH 24/25] add comments and explanations on charge and momentum sign determination --- SDL/Quintuplet.cu | 22 +++++++++++++--------- 1 file changed, 13 insertions(+), 9 deletions(-) diff --git a/SDL/Quintuplet.cu b/SDL/Quintuplet.cu index 4d6e3761..a9a5a5e0 100644 --- a/SDL/Quintuplet.cu +++ b/SDL/Quintuplet.cu @@ -643,24 +643,25 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc rt_init=mdsInGPU.anchorRt[secondMDIndex]/100; } - //start from a circle of inner T3. + // start from a circle of inner T3. // to determine the charge int charge=0; float slope3c=(y3-y_center)/(x3-x_center); float slope1c=(y1-y_center)/(x1-x_center); + // these 4 "if"s basically separate the x-y plane into 4 quarters. It determines geometrically how a circle and line slope goes and their positions, and we can get the charges correspondingly. if((y3-y_center)>0 && (y1-y_center)>0) { - if (slope3c>slope1c) charge=-1; + if (slope1c>0 && slope3c<0) charge=-1; // on x axis of a quarter, 3 hits go anti-clockwise + else if (slope1c<0 && slope3c>0) charge=1; // on x axis of a quarter, 3 hits go clockwise + else if (slope3c>slope1c) charge=-1; else if (slope3c0 && slope3c<0) charge=-1; - if (slope1c<0 && slope3c>0) charge=1; } else if((y3-y_center)<0 && (y1-y_center)<0) { - if (slope3c>slope1c) charge=-1; - else if (slope3c0) charge=1; - if (slope1c>0 && slope3c<0) charge=-1; + else if (slope1c>0 && slope3c<0) charge=-1; + else if (slope3c>slope1c) charge=-1; + else if (slope3c0) { @@ -676,6 +677,8 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc float pseudo_phi = atan((y_init-y_center)/(x_init-x_center)); //actually represent pi/2-phi, wrt helix axis z float Pt=inner_pt, Px=Pt*abs(sin(pseudo_phi)), Py=Pt*abs(cos(pseudo_phi)); + // Above line only gives you the correct value of Px and Py, but signs of Px and Py calculated below. + // We look at if the circle is clockwise or anti-clock wise, to make it simpler, we separate the x-y plane into 4 quarters. if (x_init>x_center && y_init>y_center) //1st quad { if (charge==1) Py=-Py; @@ -701,13 +704,14 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc } } - if (moduleType3==0){ + // But if the initial T5 curve goes across quarters(i.e. cross axis to separate the quarters), need special redeclaration of Px,Py signs on these to avoid errors + if (moduleType3==0){ // 0 is ps if (x4x3 && x3>x2) Px=abs(Px); if (y4y3 && y3>y2) Py=abs(Py); } - else if(moduleType3==1) + else if(moduleType3==1) // 1 is 2s { if (x3x2 && x2>x1) Px=abs(Px); From 9d5b64a74789572a2b477f2906482b954f0a4bdd Mon Sep 17 00:00:00 2001 From: YonsiG Date: Tue, 9 May 2023 12:04:59 -0700 Subject: [PATCH 25/25] put magnetic field into constant.cu file --- SDL/Constants.cu | 1 + SDL/Constants.cuh | 1 + SDL/Quintuplet.cu | 2 +- 3 files changed, 3 insertions(+), 1 deletion(-) diff --git a/SDL/Constants.cu b/SDL/Constants.cu index 70323eae..310e7d04 100644 --- a/SDL/Constants.cu +++ b/SDL/Constants.cu @@ -17,3 +17,4 @@ CUDA_CONST_VAR const float SDL::deltaZLum = 15.0; CUDA_CONST_VAR const float SDL::pixelPSZpitch = 0.15; CUDA_CONST_VAR const float SDL::strip2SZpitch = 5.0; CUDA_CONST_VAR const float SDL::pt_betaMax = 7.0f; +CUDA_CONST_VAR const float SDL::magnetic_field = 3.8112; diff --git a/SDL/Constants.cuh b/SDL/Constants.cuh index b926a4b7..4a2ff0b9 100644 --- a/SDL/Constants.cuh +++ b/SDL/Constants.cuh @@ -73,5 +73,6 @@ namespace SDL extern CUDA_CONST_VAR const float pixelPSZpitch; extern CUDA_CONST_VAR const float strip2SZpitch; extern CUDA_CONST_VAR const float pt_betaMax; + extern CUDA_CONST_VAR const float magnetic_field; } #endif diff --git a/SDL/Quintuplet.cu b/SDL/Quintuplet.cu index a9a5a5e0..81d2371f 100644 --- a/SDL/Quintuplet.cu +++ b/SDL/Quintuplet.cu @@ -729,7 +729,7 @@ __device__ bool SDL::passT5RZConstraint(struct SDL::modules& modulesInGPU, struc float Pz=(z_init-z1)/ds*Pt; float p = sqrt(Px*Px+Py*Py+Pz*Pz); - float B = 3.8112; + float B = SDL::magnetic_field; float a = -0.299792*B*charge; float zsi, rtsi;