18 #include "Identifier/Identifier.h"
33 #include "CLHEP/Geometry/Point3D.h"
37 #define AUXDATA(OBJ, TYP, NAME) \
38 static const SG::AuxElement::Accessor<TYP> acc_##NAME (#NAME); acc_##NAME(*(OBJ))
41 unsigned int makeKey(
short phi,
char eta,
char layer) {
53 m_PixelHelper(nullptr),
54 m_useSiHitsGeometryMatching(true),
55 m_firstEventWarnings(true),
61 declareProperty(
"UseSiHitsGeometryMatching", m_useSiHitsGeometryMatching=
true);
66 declareProperty(
"WriteExtendedPRDinformation", m_writeExtendedPRDinformation =
false);
69 declareProperty(
"SiClusterContainer", m_clustercontainer_key =
"PixelClusters");
72 declareProperty(
"PRD_MultiTruth", m_multiTruth_key =
"PRD_MultiTruthPixel");
74 declareProperty(
"OutputClusterContainer", m_write_xaod_key =
"PixelClusters");
77 declare(m_write_xaod_key);
78 declare(m_write_offsets);
122 return StatusCode::SUCCESS;
132 const EventContext& ctx = Gaudi::Hive::currentContext();
139 return StatusCode::FAILURE;
146 if (prdmtCollHandle.
isValid()) {
147 prdmtColl = &*prdmtCollHandle;
151 if (truthParticleLinksHandle.isValid()) {
152 truth_particle_links = truthParticleLinksHandle.cptr();
160 if (sdoCollectionHandle.
isValid()) {
161 sdoCollection = &*sdoCollectionHandle;
163 ATH_MSG_WARNING(
"SDO information requested, but SDO collection not available!");
170 if (!splitProbContainer.
isValid()) {
178 if (siHitCollectionHandle.
isValid()) {
179 for (
const SiHit& siHit: *siHitCollectionHandle) {
181 if (!siHit.isPixel())
continue;
184 siHit.getLayerDisk(),
185 siHit.getPhiModule(),
186 siHit.getEtaModule()));
189 siHits[wafer_hash].push_back(&siHit);
192 ATH_MSG_WARNING(
"SiHit information requested, but SiHit collection not available!");
199 if (!calibData_handle.
isValid()) {
202 calibData=calibData_handle.
cptr();
207 ATH_CHECK(xaod.
record(std::make_unique<xAOD::TrackMeasurementValidationContainer>(),
208 std::make_unique<xAOD::TrackMeasurementValidationAuxContainer>()));
213 unsigned int have_truth_link=0
u;
214 unsigned int missing_truth_particle=0
u;
215 unsigned int missing_parent_particle=0
u;
223 std::unordered_map< unsigned int , std::vector<unsigned int> > cluster_map;
227 (*offsets)[clusterCollection->identifyHash()] =
counter;
230 if( clusterCollection->empty() )
continue;
243 unsigned int cluster_idx = xaod->
size();
261 if(localCov.size() == 1){
264 }
else if(localCov.size() == 4){
273 std::vector< uint64_t > rdoIdentifierList;
274 rdoIdentifierList.reserve(prd->rdoList().size());
275 int rowmin=9999;
int rowmax=-9999;
276 int colmin=9999;
int colmax=-9999;
277 for(
const auto &hitIdentifier : prd->rdoList() ){
278 rdoIdentifierList.push_back( hitIdentifier.get_compact() );
282 if(rowmin >
row) rowmin =
row;
283 if(rowmax <
row) rowmax =
row;
284 if(colmin >
col) colmin =
col;
285 if(colmax <
col) colmax =
col;
295 AUXDATA(xprd,
int,phi_module) = the_phi ;
296 AUXDATA(xprd,
int,eta_module) = the_eta ;
300 cluster_map[ makeKey(the_phi, the_eta, the_layer)].push_back(cluster_idx);
305 AUXDATA(xprd,
int,nRDO) = (
int)prd->rdoList().size();
308 AUXDATA(xprd,
int,ToT) = prd->totalToT();
309 AUXDATA(xprd,
int,LVL1A) = prd->LVL1A();
312 AUXDATA(xprd,
char,gangedPixel) = (
char)prd->gangedPixel();
328 float deplVoltage = 0.0;
331 AUXDATA(xprd,
float,DepletionVoltage) = deplVoltage;
364 AUXDATA(xprd,
float,omegax) = prd->omegax();
365 AUXDATA(xprd,
float,omegay) = prd->omegay();
370 auto range{prdmtColl->equal_range(clusterId)};
371 if (truth_particle_links) {
372 std::vector<unsigned int> tp_indices;
375 if (a_truth_particle_link) {
377 if (truth_particle) {
379 tp_indices.push_back(
static_cast<int>(truth_particle->
index()));
382 ++missing_parent_particle;
387 ++missing_truth_particle;
391 AUXDATA(xprd,std::vector<unsigned int>, truth_index) = tp_indices;
401 std::vector< std::vector< int > > sdo_tracks;
423 for (
auto clusItr = xaod->
begin(); clusItr != xaod->
end(); ++clusItr ) {
424 AUXDATA(*clusItr,
char,broken) =
false;
434 for (
auto clusItr = xaod->
begin(); clusItr != xaod->
end(); ++clusItr)
441 for (
unsigned int cluster_idx : cluster_idx_list) {
442 auto pixelCluster2 = xaod->
at(cluster_idx);
443 if ( acc_layer(*pixelCluster2) !=
layer )
445 if ( acc_eta_module(*
pixelCluster) != acc_eta_module(*pixelCluster2) )
447 if ( acc_phi_module(*
pixelCluster) != acc_phi_module(*pixelCluster2) )
450 std::vector<int> barcodes2 = acc_sihit_barcode(*pixelCluster2);
453 if (
std::find(barcodes2.begin(), barcodes2.end(), bc ) == barcodes2.end())
continue;
456 acc_broken(*pixelCluster2) =
true;
466 return StatusCode::SUCCESS;
474 std::vector<int> sdo_word;
475 std::vector< std::vector< int > > sdo_depositsBarcode;
476 std::vector< std::vector< float > > sdo_depositsEnergy;
478 for(
const auto &hitIdentifier : prd->
rdoList() ){
479 auto pos = sdoCollection.find(hitIdentifier);
480 if(
pos == sdoCollection.end() )
continue;
481 sdo_word.push_back(
pos->second.word() ) ;
482 std::vector<int> sdoDepBC(
pos->second.getdeposits().size(), -1);
483 std::vector<float> sdoDepEnergy(
pos->second.getdeposits().size());
484 unsigned int nDepos{0};
485 for (
auto& deposit:
pos->second.getdeposits()) {
486 if (deposit.first) sdoDepBC[nDepos] =
HepMC::barcode(deposit.first);
488 sdoDepEnergy[nDepos] = deposit.second;
491 sdo_depositsBarcode.push_back( sdoDepBC );
492 sdo_depositsEnergy.push_back( sdoDepEnergy );
494 AUXDATA(xprd,std::vector<int>,sdo_words) = sdo_word;
495 AUXDATA(xprd,std::vector< std::vector<int> >,sdo_depositsBarcode) = sdo_depositsBarcode;
496 AUXDATA(xprd,std::vector< std::vector<float> >,sdo_depositsEnergy) = sdo_depositsEnergy;
498 return sdo_depositsBarcode;
505 const std::vector<SiHit> & matchingHits )
const
508 int numHits = matchingHits.size();
510 std::vector<float> sihit_energyDeposit(numHits,0);
511 std::vector<float> sihit_meanTime(numHits,0);
512 std::vector<int> sihit_barcode(numHits,0);
513 std::vector<int> sihit_pdgid(numHits,0);
515 std::vector<float> sihit_startPosX(numHits,0);
516 std::vector<float> sihit_startPosY(numHits,0);
517 std::vector<float> sihit_startPosZ(numHits,0);
519 std::vector<float> sihit_endPosX(numHits,0);
520 std::vector<float> sihit_endPosY(numHits,0);
521 std::vector<float> sihit_endPosZ(numHits,0);
526 for (
const auto& sihit : matchingHits ) {
527 sihit_energyDeposit[hitNumber] = sihit.energyLoss() ;
528 sihit_meanTime[hitNumber] = sihit.meanTime() ;
532 sihit_pdgid[hitNumber] = HMPL->pdg_id();
536 const HepGeom::Point3D<double>& startPos=sihit.localStartPosition();
539 sihit_startPosX[hitNumber] =
pos[0];
540 sihit_startPosY[hitNumber] =
pos[1];
541 sihit_startPosZ[hitNumber] = startPos.x();
544 const HepGeom::Point3D<double>& endPos=sihit.localEndPosition();
546 sihit_endPosX[hitNumber] =
pos[0];
547 sihit_endPosY[hitNumber] =
pos[1];
548 sihit_endPosZ[hitNumber] = endPos.x();
553 AUXDATA(xprd,std::vector<float>,sihit_energyDeposit) = sihit_energyDeposit;
554 AUXDATA(xprd,std::vector<float>,sihit_meanTime) = sihit_meanTime;
555 AUXDATA(xprd,std::vector<int>,sihit_barcode) = sihit_barcode;
556 AUXDATA(xprd,std::vector<int>,sihit_pdgid) = sihit_pdgid;
558 AUXDATA(xprd,std::vector<float>,sihit_startPosX) = sihit_startPosX;
559 AUXDATA(xprd,std::vector<float>,sihit_startPosY) = sihit_startPosY;
560 AUXDATA(xprd,std::vector<float>,sihit_startPosZ) = sihit_startPosZ;
562 AUXDATA(xprd,std::vector<float>,sihit_endPosX) = sihit_endPosX;
563 AUXDATA(xprd,std::vector<float>,sihit_endPosY) = sihit_endPosY;
564 AUXDATA(xprd,std::vector<float>,sihit_endPosZ) = sihit_endPosZ;
575 const std::vector<const SiHit*>* sihits,
576 std::vector< std::vector< int > > & trkBCs )
const
578 ATH_MSG_VERBOSE(
"Got " << sihits->size() <<
" SiHits to look through" );
579 std::vector<SiHit> matchingHits;
586 std::vector<const SiHit* > multiMatchingHits;
588 for (
const SiHit* siHit : *sihits) {
594 HepGeom::Point3D<double> averagePosition = siHit->localStartPosition() + siHit->localEndPosition();
595 averagePosition *= 0.5;
599 for(
const auto &hitIdentifier : prd->
rdoList() ){
605 multiMatchingHits.push_back(siHit);
613 for (
const auto& barcodeSDOColl : trkBCs ) {
614 if (
std::find(barcodeSDOColl.begin(),barcodeSDOColl.end(),bc) == barcodeSDOColl.end() )
continue;
615 multiMatchingHits.push_back(siHit);
623 ATH_MSG_DEBUG(
"Found " << multiMatchingHits.size() <<
" SiHit " );
624 for ( ; siHitIter != multiMatchingHits.end(); ++siHitIter) {
625 const SiHit* lowestXPos = *siHitIter;
626 const SiHit* highestXPos = *siHitIter;
630 std::vector<const SiHit* > ajoiningHits;
631 ajoiningHits.push_back( *siHitIter );
633 siHitIter2 = siHitIter+1;
635 while ( siHitIter2 != multiMatchingHits.end() ) {
644 if (std::abs((highestXPos->
localEndPosition().x()-(*siHitIter2)->localStartPosition().x()))<0.00005 &&
645 std::abs((highestXPos->
localEndPosition().y()-(*siHitIter2)->localStartPosition().y()))<0.00005 &&
646 std::abs((highestXPos->
localEndPosition().z()-(*siHitIter2)->localStartPosition().z()))<0.00005 )
648 highestXPos = *siHitIter2;
649 ajoiningHits.push_back( *siHitIter2 );
652 siHitIter2 = multiMatchingHits.erase( siHitIter2 );
653 }
else if (std::abs((lowestXPos->
localStartPosition().x()-(*siHitIter2)->localEndPosition().x()))<0.00005 &&
654 std::abs((lowestXPos->
localStartPosition().y()-(*siHitIter2)->localEndPosition().y()))<0.00005 &&
655 std::abs((lowestXPos->
localStartPosition().z()-(*siHitIter2)->localEndPosition().z()))<0.00005)
657 lowestXPos = *siHitIter2;
658 ajoiningHits.push_back( *siHitIter2 );
661 siHitIter2 = multiMatchingHits.erase( siHitIter2 );
667 if( ajoiningHits.size() == 0){
670 }
else if(ajoiningHits.size() == 1){
672 matchingHits.push_back( *ajoiningHits[0] );
676 ATH_MSG_DEBUG(
"Merging " << ajoiningHits.size() <<
" SiHits together." );
681 for(
const auto& siHit : ajoiningHits){
682 energyDep += siHit->energyLoss();
683 time += siHit->meanTime();
693 (*siHitIter)->getBarrelEndcap(),
694 (*siHitIter)->getLayerDisk(),
695 (*siHitIter)->getEtaModule(),
696 (*siHitIter)->getPhiModule(),
697 (*siHitIter)->getSide() );
698 ATH_MSG_DEBUG(
"Finished Merging " << ajoiningHits.size() <<
" SiHits together." );
715 const std::vector<Identifier>& rdos =
pixelCluster->rdoList();
718 const std::vector<int> &totList =
pixelCluster->totList();
722 std::vector<int> etaIndexList;
723 std::vector<int> phiIndexList;
724 std::vector<float> CTerm;
725 std::vector<float> ATerm;
726 std::vector<float> ETerm;
731 std::vector<Identifier>::const_iterator rdosBegin = rdos.begin();
732 std::vector<Identifier>::const_iterator rdosEnd = rdos.end();
734 ATH_MSG_VERBOSE(
" Putting together the n. " << rdos.size() <<
" rdos into a matrix.");
736 phiIndexList.reserve( rdos.size());
737 etaIndexList.reserve( rdos.size());
738 CTerm.reserve( rdos.size());
739 ATerm.reserve( rdos.size());
740 ETerm.reserve( rdos.size());
741 for (; rdosBegin!= rdosEnd; ++rdosBegin)
760 AUXDATA(xprd, std::vector<int>,rdo_phi_pixel_index) = phiIndexList;
761 AUXDATA(xprd, std::vector<int>,rdo_eta_pixel_index) = etaIndexList;
763 AUXDATA(xprd, std::vector<int>,rdo_tot) = totList;
765 AUXDATA(xprd, std::vector<float>,rdo_Cterm) = CTerm;
766 AUXDATA(xprd, std::vector<float>,rdo_Aterm) = ATerm;
767 AUXDATA(xprd, std::vector<float>,rdo_Eterm) = ETerm;
774 const unsigned int sizeX,
const unsigned int sizeY )
const
787 ATH_MSG_WARNING(
"PixelModuleDesign was not retrieved in function 'addNNInformation'");
790 const std::vector<Identifier>& rdos =
pixelCluster->rdoList();
793 const std::vector<int>& totList =
pixelCluster->totList();
801 int phiPixelIndexMin, phiPixelIndexMax, etaPixelIndexMin, etaPixelIndexMax;
804 if (!cellIdWeightedPosition.
isValid())
809 int etaPixelIndexWeightedPosition=cellIdWeightedPosition.
etaIndex();
810 int phiPixelIndexWeightedPosition=cellIdWeightedPosition.
phiIndex();
813 ATH_MSG_DEBUG(
" weighted pos phiPixelIndex: " << phiPixelIndexWeightedPosition <<
" etaPixelIndex: " << etaPixelIndexWeightedPosition );
831 double localEtaPixelIndexWeightedPosition =
w.xEta();
832 double localPhiPixelIndexWeightedPosition =
w.xPhi();
834 int centralIndexX=(sizeX-1)/2;
835 int centralIndexY=(sizeY-1)/2;
841 if (abs(phiPixelIndexWeightedPosition-phiPixelIndexMin)>centralIndexX ||
842 abs(phiPixelIndexWeightedPosition-phiPixelIndexMax)>centralIndexX)
844 ATH_MSG_DEBUG(
" Cluster too large phiPixelIndexMin " << phiPixelIndexMin <<
" phiPixelIndexMax " << phiPixelIndexMax <<
" centralX " << centralIndexX);
848 if (abs(etaPixelIndexWeightedPosition-etaPixelIndexMin)>centralIndexY ||
849 abs(etaPixelIndexWeightedPosition-etaPixelIndexMax)>centralIndexY)
851 ATH_MSG_DEBUG(
" Cluster too large etaPixelIndexMin" << etaPixelIndexMin <<
" etaPixelIndexMax " << etaPixelIndexMax <<
" centralY " << centralIndexY);
855 std::vector< std::vector<float> > matrixOfToT (sizeX, std::vector<float>(sizeY,0) );
856 std::vector< std::vector<float> > matrixOfCharge(sizeX, std::vector<float>(sizeY,0));
857 std::vector<float> vectorOfPitchesY(sizeY,0.4);
861 std::vector<Identifier>::const_iterator rdosBegin = rdos.begin();
862 std::vector<Identifier>::const_iterator rdosEnd = rdos.end();
864 auto tot = totList.begin();
866 ATH_MSG_VERBOSE(
" Putting together the n. " << rdos.size() <<
" rdos into a matrix.");
868 for (; rdosBegin!= rdosEnd; ++rdosBegin)
877 if (absphiPixelIndex <0 || absphiPixelIndex >= (
int)sizeX)
879 ATH_MSG_DEBUG(
" problem with index: " << absphiPixelIndex <<
" min: " << 0 <<
" max: " << sizeX);
883 if (absetaPixelIndex <0 || absetaPixelIndex >= (
int)sizeY)
885 ATH_MSG_DEBUG(
" problem with index: " << absetaPixelIndex <<
" min: " << 0 <<
" max: " << sizeY);
891 float pitchY = diodeParameters.
width().
xEta();
893 if ( (not totList.empty()) &&
tot != totList.end()) {
894 matrixOfToT[absphiPixelIndex][absetaPixelIndex] =*
tot;
896 }
else matrixOfToT[absphiPixelIndex][absetaPixelIndex] = -1;
899 matrixOfCharge[absphiPixelIndex][absetaPixelIndex]=*
charge;
901 }
else matrixOfCharge[absphiPixelIndex][absetaPixelIndex] = -1;
905 vectorOfPitchesY[absetaPixelIndex]=pitchY;
919 trackDir.normalize();
926 float trkphicomp = trackDir.dot(module_phiax);
927 float trketacomp = trackDir.dot(module_etaax);
928 float trknormcomp = trackDir.dot(module_normal);
929 double bowphi = atan2(trkphicomp,trknormcomp);
930 double boweta = atan2(trketacomp,trknormcomp);
932 if(bowphi > TMath::Pi()/2) bowphi -= TMath::Pi();
933 if(bowphi < -TMath::Pi()/2) bowphi += TMath::Pi();
934 int readoutside = design->readoutSide();
940 if (boweta>TMath::Pi()/2.) boweta-=TMath::Pi();
941 if (boweta<-TMath::Pi()/2.) boweta+=TMath::Pi();
945 ATH_MSG_VERBOSE(
" PhiPixelIndexWeightedPosition: " << phiPixelIndexWeightedPosition <<
" EtaPixelIndexWeightedPosition: " << etaPixelIndexWeightedPosition );
948 std::vector<float> vectorOfCharge(sizeX*sizeY,0);
949 std::vector<float> vectorOfToT(sizeX*sizeY,0);
951 for (
unsigned int u=0;
u<sizeX;
u++)
953 for (
unsigned int s=0;
s<sizeY;
s++)
955 vectorOfToT[
counter] = matrixOfToT[
u][
s];
956 vectorOfCharge[
counter] = matrixOfCharge[
u][
s];
965 AUXDATA(xprd,
int, NN_sizeX) = sizeX;
966 AUXDATA(xprd,
int, NN_sizeY) = sizeY;
969 AUXDATA(xprd,
float, NN_thetaBS) = boweta;
971 AUXDATA(xprd, std::vector<float>, NN_matrixOfToT) = vectorOfToT;
972 AUXDATA(xprd, std::vector<float>, NN_matrixOfCharge) = vectorOfCharge;
973 AUXDATA(xprd, std::vector<float>, NN_vectorOfPitchesY) = vectorOfPitchesY;
976 AUXDATA(xprd,
int, NN_etaPixelIndexWeightedPosition) = etaPixelIndexWeightedPosition;
977 AUXDATA(xprd,
int, NN_phiPixelIndexWeightedPosition) = phiPixelIndexWeightedPosition;
979 AUXDATA(xprd,
float, NN_localEtaPixelIndexWeightedPosition) = localEtaPixelIndexWeightedPosition;
980 AUXDATA(xprd,
float, NN_localPhiPixelIndexWeightedPosition) = localPhiPixelIndexWeightedPosition;
987 const std::vector<SiHit> & matchingHits )
const
1025 ATH_MSG_WARNING(
"PixelModuleDesign was not retrieved in function 'addNNTruthInfo'");
1030 unsigned hitNumber(0);
1031 for(
const auto& siHit : matchingHits ){
1033 HepGeom::Point3D<double> averagePosition = (siHit.localStartPosition() + siHit.localEndPosition()) * 0.5;
1035 ATH_MSG_VERBOSE(
"Truth Part X: " << averagePosition.y() <<
" shift " << shift <<
" Y: " << averagePosition.z() );
1038 float YposC = averagePosition.y()-shift;
1040 if (std::abs(YposC)>design->width()/2 &&
1041 std::abs(averagePosition.y())<design->width()/2)
1043 if (YposC>design->width()/2)
1045 YposC=design->width()/2-1
e-6;
1046 }
else if (YposC<-design->
width()/2)
1048 YposC=-design->width()/2+1
e-6;
1052 positionsX[hitNumber] = YposC;
1053 positionsY[hitNumber] = averagePosition.z();
1055 HepGeom::Point3D<double> deltaPosition = siHit.localEndPosition() - siHit.localStartPosition();
1057 pathlengthX[hitNumber] = deltaPosition.y();
1058 pathlengthY[hitNumber] = deltaPosition.z();
1059 pathlengthZ[hitNumber] = deltaPosition.x();
1064 InDetDD::SiCellId cellIdOfTruthPosition = design->cellIdOfPosition(siLocalTruthPosition);
1070 int truthEtaIndex = cellIdOfTruthPosition.
etaIndex();
1071 int truthPhiIndex = cellIdOfTruthPosition.
phiIndex();
1074 double pitchY = diodeParameters.
width().
xEta();
1075 double pitchX = diodeParameters.
width().
xPhi();
1092 double pixelCenterY = siLocalPositionCenter.
xEta();
1093 double pixelCenterX = siLocalPositionCenter.
xPhi();
1099 double truthIndexY = truthEtaIndex + (siLocalTruthPosition[
Trk::distEta] - pixelCenterY)/pitchY;
1100 double truthIndexX = truthPhiIndex + (siLocalTruthPosition[
Trk::distPhi] - pixelCenterX)/pitchX;
1103 positions_indexX[hitNumber] = truthIndexX - cellIdWeightedPosition.
phiIndex();
1104 positions_indexY[hitNumber] = truthIndexY - cellIdWeightedPosition.
etaIndex();
1106 HepGeom::Point3D<double> diffPositions = (siHit.localEndPosition() - siHit.localStartPosition());
1107 double bowphi = std::atan2( diffPositions.y(), diffPositions.x() );
1111 theta[hitNumber] = std::atan2(diffPositions.z() ,diffPositions.x());
1115 int readoutside = design->readoutSide();
1121 pdgid[hitNumber] =
particle->pdg_id();
1128 const auto& mother_of_particle=
vertex->particles_in().front();
1130 motherPdgid[hitNumber] = mother_of_particle->pdg_id();
1134 if(
vertex->particles_in_const_begin() !=
vertex->particles_in_const_end() ){
1136 motherPdgid[hitNumber] = (*
vertex->particles_in_const_begin())->pdg_id();
1141 chargeDep[hitNumber] = siHit.energyLoss() ;
1147 AUXDATA(xprd, std::vector<float>, NN_positionsX) = positionsX;
1148 AUXDATA(xprd, std::vector<float>, NN_positionsY) = positionsY;
1150 AUXDATA(xprd, std::vector<float>, NN_positions_indexX) = positions_indexX;
1151 AUXDATA(xprd, std::vector<float>, NN_positions_indexY) = positions_indexY;
1154 AUXDATA(xprd, std::vector<float>, NN_phi) =
phi;
1157 AUXDATA(xprd, std::vector<int>, NN_pdgid) = pdgid;
1158 AUXDATA(xprd, std::vector<float>, NN_energyDep) = chargeDep;
1159 AUXDATA(xprd, std::vector<float>, NN_trueP) = truep;
1161 AUXDATA(xprd, std::vector<int>, NN_motherBarcode) = motherBarcode;
1162 AUXDATA(xprd, std::vector<int>, NN_motherPdgid) = motherPdgid;
1165 AUXDATA(xprd, std::vector<float>, NN_pathlengthX) = pathlengthX;
1166 AUXDATA(xprd, std::vector<float>, NN_pathlengthY) = pathlengthY;
1167 AUXDATA(xprd, std::vector<float>, NN_pathlengthZ) = pathlengthZ;
1176 int *rphiPixelIndexMin,
1177 int *rphiPixelIndexMax,
1178 int *retaPixelIndexMin,
1179 int *retaPixelIndexMax )
const
1190 ATH_MSG_WARNING(
"PixelModuleDesign was not retrieved in function 'getCellIdWeightedPosition'");
1193 const std::vector<Identifier>& rdos =
pixelCluster->rdoList();
1199 std::vector<Identifier>::const_iterator rdosBegin = rdos.begin();
1200 std::vector<Identifier>::const_iterator rdosEnd = rdos.end();
1205 double sumOfCharge=0;
1207 int phiPixelIndexMin = 99999;
1208 int phiPixelIndexMax = -99999;
1209 int etaPixelIndexMin = 99999;
1210 int etaPixelIndexMax = -99999;
1212 for (; rdosBegin!= rdosEnd; ++rdosBegin, ++
charge)
1219 ATH_MSG_VERBOSE(
" Adding pixel phiPixelIndex: " << phiPixelIndex <<
" etaPixelIndex: " << etaPixelIndex <<
" charge: " << *
charge );
1237 sumOfWeightedPositions += (*charge)*siLocalPosition;
1238 sumOfCharge += (*charge);
1240 if (phiPixelIndex < phiPixelIndexMin)
1241 phiPixelIndexMin = phiPixelIndex;
1243 if (phiPixelIndex > phiPixelIndexMax)
1244 phiPixelIndexMax = phiPixelIndex;
1246 if (etaPixelIndex < etaPixelIndexMin)
1247 etaPixelIndexMin = etaPixelIndex;
1249 if (etaPixelIndex > etaPixelIndexMax)
1250 etaPixelIndexMax = etaPixelIndex;
1253 sumOfWeightedPositions /= sumOfCharge;
1255 ATH_MSG_VERBOSE (
"Wighted position: Row = " << sumOfWeightedPositions.
xRow() <<
", Col = " << sumOfWeightedPositions.
xColumn() );
1257 if(rphiPixelIndexMin) *rphiPixelIndexMin = phiPixelIndexMin;
1258 if(rphiPixelIndexMax) *rphiPixelIndexMax = phiPixelIndexMax;
1259 if(retaPixelIndexMin) *retaPixelIndexMin = etaPixelIndexMin;
1260 if(retaPixelIndexMax) *retaPixelIndexMax = etaPixelIndexMax;
1265 InDetDD::SiCellId cellIdWeightedPosition=design->cellIdOfPosition(sumOfWeightedPositions);
1268 return cellIdWeightedPosition;
1285 return StatusCode::SUCCESS;