38 square(
const double x) {
43 distance(
const std::vector<Amg::Vector2D> &vectorOfPositions,
44 const std::vector<Amg::Vector2D> &allLocalPositions,
45 const std::vector<Amg::MatrixX> &allErrorMatrix,
46 const int i,
const int j,
const int k) {
48 square(vectorOfPositions[
i][0] - allLocalPositions[0][0]) / allErrorMatrix[0](0, 0) +
49 square(vectorOfPositions[j][0] - allLocalPositions[1][0]) / allErrorMatrix[1](0, 0) +
50 square(vectorOfPositions[
k][0] - allLocalPositions[2][0]) / allErrorMatrix[2](0, 0) +
51 square(vectorOfPositions[
i][1] - allLocalPositions[0][1]) / allErrorMatrix[0](1, 1) +
52 square(vectorOfPositions[j][1] - allLocalPositions[1][1]) / allErrorMatrix[1](1, 1) +
53 square(vectorOfPositions[
k][1] - allLocalPositions[2][1]) / allErrorMatrix[2](1, 1);
58 (
const std::string &
t,
const std::string &
n,
const IInterface *
p) :
60 m_disableDistortions(false),
62 m_NNIBLcorrection(false),
64 m_NnClusterizationFactory(
"InDet::NnClusterizationFactory/NnClusterizationFactory", this),
65 m_IBLParameterSvc(
"IBLParameterSvc",
n),
66 m_doNotRecalibrateNN(false),
67 m_noNNandBroadErrors(false),
68 m_usingTIDE_Ambi(false),
69 m_splitClusterMapKey(
"") {
70 declareInterface<IRIO_OnTrackCreator>(
this);
128 constexpr
double phimin=-0.27, phimax=0.27;
130 constexpr
double etacen[
s_nbineta]={-0.,1.,1.55,1.9,2.15,2.35};
139 return StatusCode::SUCCESS;
150 return StatusCode::SUCCESS;
201 const double TOPHAT_SIGMA = 1. / std::sqrt(12.);
230 if (
pix->rdoList().empty()) {
231 ATH_MSG_WARNING(
"Pixel RDO-list size is 0, check integrity of pixel clusters! stop ROT creation.");
242 float trkphicomp = my_track.dot(my_phiax);
243 float trketacomp = my_track.dot(my_etaax);
244 float trknormcomp = my_track.dot(my_normal);
245 double bowphi = std::atan2(trkphicomp, trknormcomp);
246 double boweta = std::atan2(trketacomp, trknormcomp);
247 float etatrack = trackPar.
eta();
253 if (bowphi >
M_PI *0.5) {
256 if (bowphi < -
M_PI *0.5) {
268 double PixTrkPt = trackPar.
pT();
269 double PixTrkEta = trackPar.
eta();
270 ATH_MSG_VERBOSE(
"tanl = " << tanl <<
" readout side is " << readoutside <<
271 " module " << PixEtaModule <<
" " << PixPhiModule <<
272 " track pt, eta = " << PixTrkPt <<
" " << PixTrkEta <<
273 " track momentum phi, norm = " << trkphicomp <<
" " <<
274 trknormcomp <<
" bowphi = " << bowphi <<
" angle = " <<
angle);
276 float omegaphi =
pix->omegax();
277 float omegaeta =
pix->omegay();
278 double localphi = -9999.;
279 double localeta = -9999.;
281 const std::vector<Identifier> & rdos =
pix->rdoList();
287 for (
const auto & rId:rdos) {
296 meanpos = meanpos / rdos.size();
308 int nrows = rowmax - rowmin + 1;
309 int ncol = colmax - colmin + 1;
320 ang = 180 *
angle * M_1_PI;
335 ATH_MSG_ERROR(
"bin out of range in line " << __LINE__ <<
" of PixelClusterOnTrackTool.cxx.");
343 localphi += delta * (omegaphi - 0.5);
345 double thetaloc = -999.;
346 if (boweta > -0.5 *
M_PI && boweta <
M_PI / 2.) {
347 thetaloc = M_PI_2 - boweta;
348 }
else if (boweta > M_PI_2 && boweta <
M_PI) {
349 thetaloc = 1.5 *
M_PI - boweta;
351 thetaloc = -M_PI_2 - boweta;
353 double etaloc = -1 *
log(
tan(thetaloc * 0.5));
357 etaloc = std::abs(etaloc);
368 }
else if (ncol ==
bin + 1) {
370 }
else if (ncol ==
bin + 2) {
376 ATH_MSG_ERROR(
"bin out of range in line " << __LINE__ <<
" of PixelClusterOnTrackTool.cxx.");
383 localeta += delta * (omegaeta - 0.5);
389 localphi += deltax * (omegaphi - 0.5);
390 localeta += deltay * (omegaeta - 0.5);
398 localphi += deltax * (omegaphi - 0.5);
399 localeta += deltay * (omegaeta - 0.5);
405 localphi = meanpos.
xPhi() + shift;
406 localeta = meanpos.
xEta();
418 if (std::abs(
angle) > 1) {
422 ATH_MSG_VERBOSE(
"Shallow track with tanl = " << tanl <<
" bowphi = " <<
423 bowphi <<
" angle = " <<
angle <<
" width.z = " <<
width.z() <<
424 " errphi = " << errphi <<
" erreta = " << erreta);
426 errphi =
width.phiR() * TOPHAT_SIGMA;
427 erreta =
width.z() * TOPHAT_SIGMA;
429 errphi = (
width.phiR() / nrows) * TOPHAT_SIGMA;
430 erreta = (
width.z() / ncol) * TOPHAT_SIGMA;
438 errphi =
width.phiR() * TOPHAT_SIGMA;
447 }
else if (nrows == 2) {
453 ATH_MSG_ERROR(
"bin out of range in line " << __LINE__ <<
" of PixelClusterOnTrackTool.cxx.");
462 double etaloc = std::abs(etatrack);
464 erreta =
width.z() * TOPHAT_SIGMA;
471 ATH_MSG_ERROR(
"bin out of range in line " << __LINE__ <<
" of PixelClusterOnTrackTool.cxx.");
475 }
else if (ncol ==
bin + 1) {
477 }
else if (ncol ==
bin + 2) {
480 erreta =
width.z() * TOPHAT_SIGMA;
490 if (errphi > erreta) {
491 erreta =
width.z() * TOPHAT_SIGMA;
512 cov(0, 0) = errphi * errphi;
515 cov(1, 1) = erreta * erreta;
523 cov = Trk::ErrorScalingCast<PixelRIO_OnTrackErrorScaling>(*error_scaling)
531 glob,
pix->gangedPixel(), isbroad);
536 const double theta)
const {
547 int initial_errorStrategy;
555 newROT =
correct(rio, trackPar);
575 if (pixelPrepCluster ==
nullptr) {
582 ATH_MSG_WARNING(
"Cannot access detector element. Aborting cluster correction...");
595 float trkphicomp = my_track.dot(my_phiax);
596 float trketacomp = my_track.dot(my_etaax);
597 float trknormcomp = my_track.dot(my_normal);
598 double bowphi = std::atan2(trkphicomp, trknormcomp);
599 double boweta = std::atan2(trketacomp, trknormcomp);
619 if (!
getErrorsTIDE_Ambi(pixelPrepCluster, trackPar, finalposition, finalerrormatrix)) {
631 <<
" +/- " << std::sqrt(pixelPrepCluster->
localCovariance()(1, 1)) <<
"\n"
632 <<
" Final position x: " << finalposition[0]
633 <<
" +/- " << std::sqrt(finalerrormatrix(0, 0))
634 <<
" y: " << finalposition[1] <<
" +/- "
635 <<std::sqrt(finalerrormatrix(1, 1)) );
641 float trkphicomp = my_track.dot(my_phiax);
642 float trketacomp = my_track.dot(my_etaax);
643 float trknormcomp = my_track.dot(my_normal);
644 double bowphi = std::atan2(trkphicomp, trknormcomp);
645 double boweta = std::atan2(trketacomp, trknormcomp);
656 cov = Trk::ErrorScalingCast<PixelRIO_OnTrackErrorScaling>(*error_scaling)
681 std::vector<Amg::Vector2D> vectorOfPositions;
682 int numberOfSubclusters = 1;
683 vectorOfPositions.push_back(pixelPrepCluster->
localPosition());
687 InDet::PixelGangedClusterAmbiguities::const_iterator mapBegin = splitClusterMap->begin();
688 InDet::PixelGangedClusterAmbiguities::const_iterator mapEnd = splitClusterMap->end();
689 for (InDet::PixelGangedClusterAmbiguities::const_iterator mapIter = mapBegin; mapIter != mapEnd; ++mapIter) {
692 if (
first == pixelPrepCluster &&
second != pixelPrepCluster) {
693 ATH_MSG_DEBUG(
"Found additional split cluster in ambiguity map (+=1).");
694 numberOfSubclusters += 1;
700 if (pixelAddCluster ==
nullptr) {
701 ATH_MSG_WARNING(
"Pixel ambiguity map has empty pixel cluster. Please DEBUG!");
704 vectorOfPositions.push_back(pixelAddCluster->
localPosition());
717 "Parameters are not at a plane ! Aborting cluster correction... ");
721 std::vector<Amg::Vector2D> allLocalPositions;
722 std::vector<Amg::MatrixX> allErrorMatrix;
728 numberOfSubclusters);
730 if (allLocalPositions.empty()) {
731 ATH_MSG_DEBUG(
" Cluster cannot be treated by NN. Giving back to default clusterization " );
736 if (allLocalPositions.size() !=
size_t(numberOfSubclusters)) {
737 ATH_MSG_WARNING(
"Returned position vector size " << allLocalPositions.size() <<
738 " not according to expected number of subclusters: " << numberOfSubclusters <<
". Abort cluster correction..." );
746 if (numberOfSubclusters == 1) {
747 finalposition = allLocalPositions[0];
748 finalerrormatrix = allErrorMatrix[0];
751 else if (numberOfSubclusters == 2) {
753 square(vectorOfPositions[0][0] - allLocalPositions[0][0]) / allErrorMatrix[0](0, 0) +
754 square(vectorOfPositions[1][0] - allLocalPositions[1][0]) / allErrorMatrix[1](0, 0) +
755 square(vectorOfPositions[0][1] - allLocalPositions[0][1]) / allErrorMatrix[0](1, 1) +
756 square(vectorOfPositions[1][1] - allLocalPositions[1][1]) / allErrorMatrix[1](1, 1);
759 square(vectorOfPositions[1][0] - allLocalPositions[0][0]) / allErrorMatrix[0](0, 0) +
760 square(vectorOfPositions[0][0] - allLocalPositions[1][0]) / allErrorMatrix[1](0, 0) +
761 square(vectorOfPositions[1][1] - allLocalPositions[0][1]) / allErrorMatrix[0](1, 1) +
762 square(vectorOfPositions[0][1] - allLocalPositions[1][1]) / allErrorMatrix[1](1, 1);
765 " Old pix (1) x: " << vectorOfPositions[0][0] <<
" y: " << vectorOfPositions[0][1] <<
"\n"
766 <<
" Old pix (2) x: " << vectorOfPositions[1][0] <<
" y: " << vectorOfPositions[1][1] <<
"\n"
767 <<
" Pix (1) x: " << allLocalPositions[0][0] <<
" +/- " << std::sqrt(allErrorMatrix[0](0, 0))
768 <<
" y: " << allLocalPositions[0][1] <<
" +/- " << std::sqrt(allErrorMatrix[0](1, 1)) <<
"\n"
769 <<
" Pix (2) x: " << allLocalPositions[1][0] <<
" +/- " << std::sqrt(allErrorMatrix[1](0, 0))
770 <<
" y: " << allLocalPositions[1][1] <<
" +/- " << std::sqrt(allErrorMatrix[1](1, 1)) <<
"\n"
771 <<
" Old (1) new (1) dist: " << std::sqrt(distancesq1) <<
" Old (1) new (2) " << std::sqrt(distancesq2) );
774 if (distancesq1 < distancesq2) {
775 finalposition = allLocalPositions[0];
776 finalerrormatrix = allErrorMatrix[0];
778 finalposition = allLocalPositions[1];
779 finalerrormatrix = allErrorMatrix[1];
784 else if (numberOfSubclusters == 3) {
787 distances[0] =
distance(vectorOfPositions, allLocalPositions, allErrorMatrix, 0, 1, 2);
788 distances[1] =
distance(vectorOfPositions, allLocalPositions, allErrorMatrix, 0, 2, 1);
789 distances[2] =
distance(vectorOfPositions, allLocalPositions, allErrorMatrix, 1, 0, 2);
790 distances[3] =
distance(vectorOfPositions, allLocalPositions, allErrorMatrix, 1, 2, 0);
791 distances[4] =
distance(vectorOfPositions, allLocalPositions, allErrorMatrix, 2, 0, 1);
792 distances[5] =
distance(vectorOfPositions, allLocalPositions, allErrorMatrix, 2, 1, 0);
794 int smallestDistanceIndex = -10;
795 double minDistance = 1e10;
797 for (
int i = 0;
i < 6;
i++) {
800 if (distances[
i] < minDistance) {
801 minDistance = distances[
i];
802 smallestDistanceIndex =
i;
806 ATH_MSG_DEBUG(
" The minimum distance is : " << minDistance <<
" for index: " << smallestDistanceIndex);
808 if (smallestDistanceIndex == 0 || smallestDistanceIndex == 1) {
809 finalposition = allLocalPositions[0];
810 finalerrormatrix = allErrorMatrix[0];
812 if (smallestDistanceIndex == 2 || smallestDistanceIndex == 4) {
813 finalposition = allLocalPositions[1];
814 finalerrormatrix = allErrorMatrix[1];
816 if (smallestDistanceIndex == 3 || smallestDistanceIndex == 5) {
817 finalposition = allLocalPositions[2];
818 finalerrormatrix = allErrorMatrix[2];
830 std::vector<Amg::Vector2D> vectorOfPositions;
831 int numberOfSubclusters = 1;
834 numberOfSubclusters = 1 + splitClusterMap->count(pixelPrepCluster);
836 if (splitClusterMap->count(pixelPrepCluster) == 0 && splitProb.
isSplit()) {
837 numberOfSubclusters = 2;
839 if (splitClusterMap->count(pixelPrepCluster) != 0 && !splitProb.
isSplit()) {
840 numberOfSubclusters = 1;
847 ATH_MSG_WARNING(
"Parameters are not at a plane surface ! Aborting cluster "
852 std::vector<Amg::Vector2D> allLocalPositions;
853 std::vector<Amg::MatrixX> allErrorMatrix;
859 numberOfSubclusters);
861 if (allLocalPositions.empty()) {
863 " Cluster cannot be treated by NN. Giving back to default clusterization, too big: " <<
868 if (allLocalPositions.size() !=
size_t(numberOfSubclusters)) {
870 "Returned position vector size " << allLocalPositions.size() <<
" not according to expected number of subclusters: " << numberOfSubclusters <<
871 ". Abort cluster correction...");
878 if (numberOfSubclusters == 1) {
879 finalposition = allLocalPositions[0];
880 finalerrormatrix = allErrorMatrix[0];
888 if (trackPar.covariance()) {
889 localerr =
Amg::Vector2D(std::sqrt((*trackPar.covariance())(0, 0)), std::sqrt((*trackPar.covariance())(1, 1)));
892 double minDistance(1e300);
895 for (
unsigned int i(0);
i < allLocalPositions.size(); ++
i) {
897 square(localpos[0] - allLocalPositions[
i][0]) / localerr[0]
898 + square(localpos[1] - allLocalPositions[
i][1]) / localerr[1];
906 finalposition = allLocalPositions[
index];
907 finalerrormatrix = allErrorMatrix[
index];