37 square(
const double x) {
42 distance(
const std::vector<Amg::Vector2D> &vectorOfPositions,
43 const std::vector<Amg::Vector2D> &allLocalPositions,
44 const std::vector<Amg::MatrixX> &allErrorMatrix,
45 const int i,
const int j,
const int k) {
47 square(vectorOfPositions[
i][0] - allLocalPositions[0][0]) / allErrorMatrix[0](0, 0) +
48 square(vectorOfPositions[j][0] - allLocalPositions[1][0]) / allErrorMatrix[1](0, 0) +
49 square(vectorOfPositions[
k][0] - allLocalPositions[2][0]) / allErrorMatrix[2](0, 0) +
50 square(vectorOfPositions[
i][1] - allLocalPositions[0][1]) / allErrorMatrix[0](1, 1) +
51 square(vectorOfPositions[j][1] - allLocalPositions[1][1]) / allErrorMatrix[1](1, 1) +
52 square(vectorOfPositions[
k][1] - allLocalPositions[2][1]) / allErrorMatrix[2](1, 1);
60 (
const std::string &
t,
const std::string &
n,
const IInterface *
p) :
63 m_NnClusterizationFactory(
"InDet::NnClusterizationFactory/NnClusterizationFactory", this),
64 m_doNotRecalibrateNN(false),
65 m_noNNandBroadErrors(false),
66 m_usingTIDE_Ambi(false),
67 m_splitClusterMapKey(
"")
69 declareInterface<IRIO_OnTrackCreator>(
this);
109 return StatusCode::SUCCESS;
161 const double TOPHAT_SIGMA = 1. / std::sqrt(12.);
188 if (
pix->rdoList().empty()) {
189 ATH_MSG_WARNING(
"Pixel RDO-list size is 0, check integrity of pixel clusters! stop ROT creation.");
200 float trkphicomp = my_track.dot(my_phiax);
201 float trketacomp = my_track.dot(my_etaax);
202 float trknormcomp = my_track.dot(my_normal);
203 double bowphi = std::atan2(trkphicomp, trknormcomp);
204 double boweta = std::atan2(trketacomp, trknormcomp);
210 if (bowphi >
M_PI *0.5) {
213 if (bowphi < -
M_PI *0.5) {
222 double thetaloc = -999.;
223 if (boweta > -0.5 *
M_PI && boweta <
M_PI / 2.) {
224 thetaloc = M_PI_2 - boweta;
225 }
else if (boweta > M_PI_2 && boweta <
M_PI) {
226 thetaloc = 1.5 *
M_PI - boweta;
228 thetaloc = -M_PI_2 - boweta;
230 double etaloc = -1 *
log(
tan(thetaloc * 0.5));
236 double PixTrkPt = trackPar.
pT();
237 double PixTrkEta = trackPar.
eta();
238 ATH_MSG_VERBOSE(
"tanl = " << tanl <<
" readout side is " << readoutside <<
239 " module " << PixEtaModule <<
" " << PixPhiModule <<
240 " track pt, eta = " << PixTrkPt <<
" " << PixTrkEta <<
241 " track momentum phi, norm = " << trkphicomp <<
" " <<
242 trknormcomp <<
" bowphi = " << bowphi <<
" angle = " <<
angle);
244 float omegaphi =
pix->omegax();
245 float omegaeta =
pix->omegay();
246 double localphi = -9999.;
247 double localeta = -9999.;
249 const std::vector<Identifier> & rdos =
pix->rdoList();
255 for (
const auto & rId:rdos) {
264 meanpos = meanpos / rdos.size();
276 int nrows = rowmax - rowmin + 1;
277 int ncol = colmax - colmin + 1;
287 double delta_phi = nrows != 1 ? delta.first : 0.;
288 double delta_eta = ncol != 1 ? delta.second : 0.;
290 localeta += delta_eta*(omegaeta-0.5);
294 localphi = meanpos.
xPhi() + shift;
295 localeta = meanpos.
xEta();
307 if (std::abs(
angle) > 1) {
311 ATH_MSG_VERBOSE(
"Shallow track with tanl = " << tanl <<
" bowphi = " <<
312 bowphi <<
" angle = " <<
angle <<
" width.z = " <<
width.z() <<
313 " errphi = " << errphi <<
" erreta = " << erreta);
315 errphi =
width.phiR() * TOPHAT_SIGMA;
316 erreta =
width.z() * TOPHAT_SIGMA;
318 errphi = (
width.phiR() / nrows) * TOPHAT_SIGMA;
319 erreta = (
width.z() / ncol) * TOPHAT_SIGMA;
322 errphi = nrows != 1 ? delta_err.first : (
width.phiR()/nrows)*TOPHAT_SIGMA;
323 erreta = ncol != 1 ? delta_err.second : (
width.z()/ncol)*TOPHAT_SIGMA;
338 cov(0, 0) = errphi * errphi;
341 cov(1, 1) = erreta * erreta;
349 cov = Trk::ErrorScalingCast<PixelRIO_OnTrackErrorScaling>(*error_scaling)
356 iH, glob,
pix->gangedPixel(), isbroad);
364 int initial_errorStrategy;
372 newROT =
correct(rio, trackPar);
392 if (pixelPrepCluster ==
nullptr) {
399 ATH_MSG_WARNING(
"Cannot access detector element. Aborting cluster correction...");
413 std::move(
cov), iH, glob,
423 if (!
getErrorsTIDE_Ambi(pixelPrepCluster, trackPar, finalposition, finalerrormatrix)) {
435 <<
" +/- " << std::sqrt(pixelPrepCluster->
localCovariance()(1, 1)) <<
"\n"
436 <<
" Final position x: " << finalposition[0]
437 <<
" +/- " << std::sqrt(finalerrormatrix(0, 0))
438 <<
" y: " << finalposition[1] <<
" +/- "
439 <<std::sqrt(finalerrormatrix(1, 1)) );
445 cov = Trk::ErrorScalingCast<PixelRIO_OnTrackErrorScaling>(*error_scaling)
471 std::vector<Amg::Vector2D> vectorOfPositions;
472 int numberOfSubclusters = 1;
473 vectorOfPositions.push_back(pixelPrepCluster->
localPosition());
477 InDet::PixelGangedClusterAmbiguities::const_iterator mapBegin = splitClusterMap->begin();
478 InDet::PixelGangedClusterAmbiguities::const_iterator mapEnd = splitClusterMap->end();
479 for (InDet::PixelGangedClusterAmbiguities::const_iterator mapIter = mapBegin; mapIter != mapEnd; ++mapIter) {
482 if (
first == pixelPrepCluster &&
second != pixelPrepCluster) {
483 ATH_MSG_DEBUG(
"Found additional split cluster in ambiguity map (+=1).");
484 numberOfSubclusters += 1;
490 if (pixelAddCluster ==
nullptr) {
491 ATH_MSG_WARNING(
"Pixel ambiguity map has empty pixel cluster. Please DEBUG!");
494 vectorOfPositions.push_back(pixelAddCluster->
localPosition());
507 "Parameters are not at a plane ! Aborting cluster correction... ");
511 std::vector<Amg::Vector2D> allLocalPositions;
512 std::vector<Amg::MatrixX> allErrorMatrix;
518 numberOfSubclusters);
520 if (allLocalPositions.empty()) {
521 ATH_MSG_DEBUG(
" Cluster cannot be treated by NN. Giving back to default clusterization " );
526 if (allLocalPositions.size() !=
size_t(numberOfSubclusters)) {
527 ATH_MSG_WARNING(
"Returned position vector size " << allLocalPositions.size() <<
528 " not according to expected number of subclusters: " << numberOfSubclusters <<
". Abort cluster correction..." );
536 if (numberOfSubclusters == 1) {
537 finalposition = allLocalPositions[0];
538 finalerrormatrix = allErrorMatrix[0];
541 else if (numberOfSubclusters == 2) {
543 square(vectorOfPositions[0][0] - allLocalPositions[0][0]) / allErrorMatrix[0](0, 0) +
544 square(vectorOfPositions[1][0] - allLocalPositions[1][0]) / allErrorMatrix[1](0, 0) +
545 square(vectorOfPositions[0][1] - allLocalPositions[0][1]) / allErrorMatrix[0](1, 1) +
546 square(vectorOfPositions[1][1] - allLocalPositions[1][1]) / allErrorMatrix[1](1, 1);
549 square(vectorOfPositions[1][0] - allLocalPositions[0][0]) / allErrorMatrix[0](0, 0) +
550 square(vectorOfPositions[0][0] - allLocalPositions[1][0]) / allErrorMatrix[1](0, 0) +
551 square(vectorOfPositions[1][1] - allLocalPositions[0][1]) / allErrorMatrix[0](1, 1) +
552 square(vectorOfPositions[0][1] - allLocalPositions[1][1]) / allErrorMatrix[1](1, 1);
555 " Old pix (1) x: " << vectorOfPositions[0][0] <<
" y: " << vectorOfPositions[0][1] <<
"\n"
556 <<
" Old pix (2) x: " << vectorOfPositions[1][0] <<
" y: " << vectorOfPositions[1][1] <<
"\n"
557 <<
" Pix (1) x: " << allLocalPositions[0][0] <<
" +/- " << std::sqrt(allErrorMatrix[0](0, 0))
558 <<
" y: " << allLocalPositions[0][1] <<
" +/- " << std::sqrt(allErrorMatrix[0](1, 1)) <<
"\n"
559 <<
" Pix (2) x: " << allLocalPositions[1][0] <<
" +/- " << std::sqrt(allErrorMatrix[1](0, 0))
560 <<
" y: " << allLocalPositions[1][1] <<
" +/- " << std::sqrt(allErrorMatrix[1](1, 1)) <<
"\n"
561 <<
" Old (1) new (1) dist: " << std::sqrt(distancesq1) <<
" Old (1) new (2) " << std::sqrt(distancesq2) );
564 if (distancesq1 < distancesq2) {
565 finalposition = allLocalPositions[0];
566 finalerrormatrix = allErrorMatrix[0];
568 finalposition = allLocalPositions[1];
569 finalerrormatrix = allErrorMatrix[1];
574 else if (numberOfSubclusters == 3) {
577 distances[0] =
distance(vectorOfPositions, allLocalPositions, allErrorMatrix, 0, 1, 2);
578 distances[1] =
distance(vectorOfPositions, allLocalPositions, allErrorMatrix, 0, 2, 1);
579 distances[2] =
distance(vectorOfPositions, allLocalPositions, allErrorMatrix, 1, 0, 2);
580 distances[3] =
distance(vectorOfPositions, allLocalPositions, allErrorMatrix, 1, 2, 0);
581 distances[4] =
distance(vectorOfPositions, allLocalPositions, allErrorMatrix, 2, 0, 1);
582 distances[5] =
distance(vectorOfPositions, allLocalPositions, allErrorMatrix, 2, 1, 0);
584 int smallestDistanceIndex = -10;
585 double minDistance = 1e10;
587 for (
int i = 0;
i < 6;
i++) {
590 if (distances[
i] < minDistance) {
591 minDistance = distances[
i];
592 smallestDistanceIndex =
i;
596 ATH_MSG_DEBUG(
" The minimum distance is : " << minDistance <<
" for index: " << smallestDistanceIndex);
598 if (smallestDistanceIndex == 0 || smallestDistanceIndex == 1) {
599 finalposition = allLocalPositions[0];
600 finalerrormatrix = allErrorMatrix[0];
602 if (smallestDistanceIndex == 2 || smallestDistanceIndex == 4) {
603 finalposition = allLocalPositions[1];
604 finalerrormatrix = allErrorMatrix[1];
606 if (smallestDistanceIndex == 3 || smallestDistanceIndex == 5) {
607 finalposition = allLocalPositions[2];
608 finalerrormatrix = allErrorMatrix[2];
620 std::vector<Amg::Vector2D> vectorOfPositions;
621 int numberOfSubclusters = 1;
624 numberOfSubclusters = 1 + splitClusterMap->count(pixelPrepCluster);
626 if (splitClusterMap->count(pixelPrepCluster) == 0 && splitProb.
isSplit()) {
627 numberOfSubclusters = 2;
629 if (splitClusterMap->count(pixelPrepCluster) != 0 && !splitProb.
isSplit()) {
630 numberOfSubclusters = 1;
637 ATH_MSG_WARNING(
"Parameters are not at a plane surface ! Aborting cluster "
642 std::vector<Amg::Vector2D> allLocalPositions;
643 std::vector<Amg::MatrixX> allErrorMatrix;
649 numberOfSubclusters);
651 if (allLocalPositions.empty()) {
653 " Cluster cannot be treated by NN. Giving back to default clusterization, too big: " <<
658 if (allLocalPositions.size() !=
size_t(numberOfSubclusters)) {
660 "Returned position vector size " << allLocalPositions.size() <<
" not according to expected number of subclusters: " << numberOfSubclusters <<
661 ". Abort cluster correction...");
668 if (numberOfSubclusters == 1) {
669 finalposition = allLocalPositions[0];
670 finalerrormatrix = allErrorMatrix[0];
678 if (trackPar.covariance()) {
679 localerr =
Amg::Vector2D(std::sqrt((*trackPar.covariance())(0, 0)), std::sqrt((*trackPar.covariance())(1, 1)));
682 double minDistance(1e300);
685 for (
unsigned int i(0);
i < allLocalPositions.size(); ++
i) {
687 square(localpos[0] - allLocalPositions[
i][0]) / localerr[0]
688 + square(localpos[1] - allLocalPositions[
i][1]) / localerr[1];
696 finalposition = allLocalPositions[
index];
697 finalerrormatrix = allErrorMatrix[
index];