User will overwrite this function. Histogram booking is no longer done in C++. This function is called in execute once the filters are all passed.
242 {
244
245 const std::map<Identifier, unsigned int>* badChips{nullptr};
246 if (m_vetoBadChips) {
247 badChips = m_configConditions->badChips(ctx);
248 }
249
250 double timecor{-20.};
251 if (m_useTRTPhase or m_isCosmic) {
254 timecor = theComTime->getTime();
256 } else {
257 timecor = -18.;
259 }
260 }
261
262 const EventIDBase& pEvent{ctx.eventID()};
263 unsigned BCID{pEvent.bunch_crossing_id()};
264 SG::ReadCondHandle<BunchCrossingCondData> bcidHdl(m_bunchCrossingKey,ctx);
265 if (!bcidHdl.isValid()) {
266 ATH_MSG_ERROR(
"Unable to retrieve BunchCrossing conditions object" );
267 return StatusCode::FAILURE;
268 }
269 const BunchCrossingCondData* bcData{*bcidHdl};
271
273 const AtlasFieldCacheCondObj* fieldCondObj{*fieldHandle};
274 if (fieldCondObj==nullptr) {
275 ATH_MSG_ERROR(
"AtlasFieldCacheCondObj cannot be retrieved.");
276 return StatusCode::RECOVERABLE;
277 }
278 MagField::AtlasFieldCache fieldCache;
281
282
283 SG::ReadHandle<TrackCollection>tracks{
m_TrackName, ctx};
286 return StatusCode::SUCCESS;
287 } else {
289 }
290
291 SG::ReadHandle<InDet::SCT_ClusterContainer> p_sctclcontainer{
m_sctContainerName, ctx};
292 if (not p_sctclcontainer.
isValid()) {
294 return StatusCode::SUCCESS;
295 }
296
297
299 if (
failCut(tracks->size() <= m_maxTracks,
"# of tracks cut")) {
300 return StatusCode::SUCCESS;
301 }
302 }
303
304 SG::ReadCondHandle<InDetDD::SiDetectorElementCollection> sctDetEle{
m_SCTDetEleCollKey, ctx};
305 const InDetDD::SiDetectorElementCollection* elements{sctDetEle.
retrieve()};
306 if (elements==nullptr) {
308 return StatusCode::FAILURE;
309 }
310
311
312
313 for (const Trk::Track* pthisTrack: *tracks) {
314 if (pthisTrack==nullptr) {
315 continue;
316 }
317 if (
failCut(pthisTrack and pthisTrack->trackParameters() and pthisTrack->trackParameters()->size(),
318 "track cut: presence")) {
319 continue;
320 }
321
323 "track cut: inside-out only")) {
324 continue;
325 }
326 if (pthisTrack->perigeeParameters() == nullptr) {
327 continue;
328 }
329 const Trk::Perigee* perigee{pthisTrack->perigeeParameters()};
334
335 if (solenoidOn and
failCut(perigee->
pT() >= m_minPt,
"track cut: Min Pt")) {
336 continue;
337 }
338 if (not m_isCosmic and
failCut(std::abs(d0) <= m_maxD0,
"track cut: max D0")) {
339 continue;
340 }
341 if (m_maxZ0sinTheta and
failCut(std::abs(z0 *
sin(perigeeTheta)) <= m_maxZ0sinTheta,
"track cut: Max Z0sinTheta")) {
342 continue;
343 }
344 }
345
346
347 for (const Trk::Track* pthisTrack: *tracks) {
348
349
351 if (pthisTrack==nullptr) {
352 continue;
353 }
354 if (
failCut(pthisTrack and pthisTrack->trackParameters() and pthisTrack->trackParameters()->size(),
355 "track cut: presence")) {
356 continue;
357 }
358
360 "track cut: inside-out only")) {
361 continue;
362 }
363 if (pthisTrack->perigeeParameters() == nullptr) {
364 continue;
365 }
366 const Trk::Perigee* perigee{pthisTrack->perigeeParameters()};
371
372 if (
failCut(perigee->
pT() >= m_minPt,
"track cut: Min Pt")) {
373 continue;
374 }
375 if (not m_isCosmic and
failCut(std::abs(d0) <= m_maxD0,
"track cut: max D0")) {
376 continue;
377 }
378 if (m_maxZ0sinTheta and
failCut(std::abs(z0 *
sin(perigeeTheta)) <= m_maxZ0sinTheta,
"track cut: Max Z0sinTheta")) {
379 continue;
380 }
381
382 const Trk::TrackSummary*
summary{pthisTrack->trackSummary()};
383
385 continue;
386 }
387
388 std::unique_ptr<const Trk::Track> trackWithHoles(
m_holeSearchTool->getTrackWithHoles(*pthisTrack));
389 if (not trackWithHoles) {
391 continue;
392 }
393 ATH_MSG_VERBOSE(
"Found " << trackWithHoles->trackStateOnSurfaces()->size() <<
" states on track");
394
396 0, 0, 0
397 };
398 int pixelNHits{0};
399 int pixelNHoles{0};
400 int trtNHits{0};
401
404
405
406
407
408
409
410 std::map < Identifier, double > mapOfTrackHitResiduals;
411 double zmin = std::numeric_limits<float>::max();
412 double zmax = -std::numeric_limits<float>::max();
413 double zpos{0.};
414 float layerSide{-1};
415 float min_layerSide{999.};
416 float max_layerSide{-1.};
417 Identifier surfaceID;
418
419
420 for (const Trk::TrackStateOnSurface* tsos: *(trackWithHoles->trackStateOnSurfaces())) {
421 surfaceID = surfaceOnTrackIdentifier(tsos);
422
424 continue;
425 }
426
427
428 int waferIndex = -1;
429
434 }
435
438 pixelNHits++;
439 }
441 trtNHits++;
442 }
445 mapOfTrackHitResiduals[surfaceID] =
getResidual(surfaceID, tsos->trackParameters(), &*p_sctclcontainer);
446 sctNHitsPerRegion[waferIndex]++;
447 }
450 pixelNHoles++;
452 sctNHolesPerRegion[waferIndex]++;
453 }
454 }
455
457
458 if (m_isCosmic) {
459 if (tsos->trackParameters()) {
460 zpos = tsos->trackParameters()->position().z();
461 zmax = std::max(zpos, zmax);
462 zmin = std::min(zpos, zmin);
463 } else {
464 ATH_MSG_WARNING(
"No track parameter found. Zmin and Zmax not recalculated.");
465 }
466 } else {
470 min_layerSide = std::min(min_layerSide, layerSide);
471 max_layerSide = std::max(max_layerSide, layerSide);
473 min_layerSide = -1;
476 }
477 }
478 }
479 }
480
482 std::vector<bool> layersCrossedByTrack[
N_REGIONS];
483 std::vector<int> nHolesOnLayer[
N_REGIONS];
484 std::vector<int> nHolesDistOnLayer[
N_REGIONS];
486 nHolesDistOnLayer[
i].resize(
n_layers[i] * 2, 0);
487 nHolesOnLayer[
i].resize(
n_layers[i] * 2, 0);
488 layersCrossedByTrack[
i].resize(
n_layers[i] * 2,
false);
489 }
490
491
492 for (const Trk::TrackStateOnSurface* tsos: *(trackWithHoles->trackStateOnSurfaces())) {
494 surfaceID = surfaceOnTrackIdentifier(tsos);
495
497 continue;
498 }
499
500
506
508
509 Int_t sctNHitsExceptThisWafer{0};
510 Int_t sctNHolesExceptThisWafer{0};
511
513 if (i != waferIndex) {
514 sctNHitsExceptThisWafer += sctNHitsPerRegion[
i];
515 sctNHolesExceptThisWafer += sctNHolesPerRegion[
i];
516 }
517 }
518
519
520
521
522
523 if ((unsigned int)(sctNHitsExceptThisWafer + pixelNHits) < m_minSiHits) {
524 ATH_MSG_VERBOSE(
"This track is rejected due to the number of hits: " << sctNHitsExceptThisWafer * pixelNHits);
525 continue;
526 }
527 if ((unsigned int)(sctNHolesExceptThisWafer + pixelNHoles) > m_maxSiHoles) {
528 ATH_MSG_VERBOSE(
"This track is rejected due to the number of holes: " << sctNHolesExceptThisWafer * pixelNHoles);
529 continue;
530 }
531
532 std::string etaPhiSuffix = "_" + std::to_string(layer) + "_" + std::to_string(side);
534 if (detIndex == -1) {
535 ATH_MSG_WARNING(
"The detector region (barrel, endcap A, endcap C) could not be determined");
536 continue;
537 }
541 float layerPlusHalfSide{
static_cast<float>(
layer) +
static_cast<float>(side) * 0.5f};
542 float dedicated_layerPlusHalfSide{
static_cast<float>(
layer) +
static_cast<float>((side + 1) % 2) * 0.5f};
544 double trackHitResidual{
getResidual(surfaceID, trkParamOnSurface, &*p_sctclcontainer)};
545
547
552 }
553
554 bool otherFaceFound{false};
555 IdentifierHash otherSideHash;
556 Identifier otherSideSurfaceID;
559 m_sctId->
get_id(otherSideHash, otherSideSurfaceID, &context);
560 otherFaceFound = mapOfTrackHitResiduals.find(otherSideSurfaceID) != mapOfTrackHitResiduals.end();
561
562 int nOther{sctNHits};
564 --nOther;
565 }
566
567
568 double phiUp{90.};
572 }
573
574 if (m_useSCTorTRT) {
575 if (
failCut(trtNHits >= m_minTRTHits or
576 sctNHits >= m_minSCTHits, "track cut: min TRT or SCT hits")) {
577 continue;
578 }
579 } else {
580 if (
failCut(trtNHits >= m_minTRTHits,
"track cut: min TRT hits")) {
581 continue;
582 }
583 if (
failCut(sctNHits >= m_minSCTHits,
"track cut: min SCT hits")) {
584 continue;
585 }
586 if (
failCut(pixelNHits >= m_minPixelHits,
"track cut: min Pixel hits")) {
587 continue;
588 }
589 }
590
591 if (
failCut(nOther >= m_minOtherHits,
"track cut: minOtherHits")) {
592 continue;
593 }
594
595 ATH_MSG_DEBUG(
"Use TRT phase " << m_useTRTPhase <<
" is cosmic? " << m_isCosmic <<
" timecor " << timecor);
596 if (m_useTRTPhase or m_isCosmic) {
597 if (timecor == 0) {
598 continue;
599 }
600 static const double tmin{-15.};
601 static const double tmax{10.};
602 if (
failCut((timecor >= tmin) and (timecor <= tmax),
"track cut: timing cut")) {
603 continue;
604 }
606 }
607
608 bool enclosingHits{true};
609 if (m_isCosmic) {
610 if (tsos->trackParameters()) {
611 zpos = tsos->trackParameters()->position().z();
612 enclosingHits = ((zpos >
zmin) and (zpos < zmax));
613 } else {
614 ATH_MSG_WARNING(
"No track parameters found. Cannot determine whether it is an enclosed hit.");
615 enclosingHits = false;
616 }
617 } else {
620 + (
static_cast<float>(
m_sctId->
side(surfaceID)) == 0) * 0.5;
621 enclosingHits = ((layerSide > min_layerSide) and (layerSide < max_layerSide));
622 }
623
624 if (m_requireEnclosingHits and
625 (not (layerPlusHalfSide == 0.5)) and
626 (not ((isub == 1) and (layerPlusHalfSide == 3))) and
627 (not (layerPlusHalfSide == 8))) {
628 if (
failCut(enclosingHits,
"hit cut: enclosing hits")) {
629 continue;
630 }
631 }
632
633
634 double chi2{trackWithHoles->fitQuality()->chiSquared()};
635 int ndf{trackWithHoles->fitQuality()->numberDoF()};
636 double chi2_div_ndf{
ndf > 0. ?
chi2 /
ndf : -1.};
637
638 if (
failCut(std::abs(phiUp) <= m_maxPhiAngle,
"hit cut: incidence angle")) {
639 continue;
640 }
641
642 if (
failCut((ndf > 0) and (chi2_div_ndf <= m_maxChi2),
"track cut: chi2 cut")) {
643 continue;
644 }
645
646 if (m_requireOtherFace and
failCut(otherFaceFound,
"hit cut: other face found")) {
647 continue;
648 }
649
650 if (not trkParamOnSurface) continue;
651 double xl{trkParamOnSurface->localPosition()[0]};
652 double yl{trkParamOnSurface->localPosition()[1]};
653
654
655 bool insideGuardRing{true};
658 static const float yGuard{3.};
659 if (xl < -30.7 + xGuard) {
660 insideGuardRing = false;
661 }
662 if (xl > 30.7 - xGuard) {
663 insideGuardRing = false;
664 }
665
666 static const double yend{63.960 + 0.03 - 1.};
667 static const double ydead{2.06 / 2.};
668 if (yl > yend - yGuard) {
669 insideGuardRing = false;
670 }
671 if (yl < -yend + yGuard) {
672 insideGuardRing = false;
673 }
674 if ((yl < ydead + yGuard) and (yl > -ydead - yGuard)) {
675 insideGuardRing = false;
676 }
677 } else {
678
679 insideGuardRing = true;
680 }
681
682 if (m_requireGuardRing and
failCut(insideGuardRing,
"hit cut: inside guard ring")) {
683 continue;
684 }
685
686
687 if (m_vetoBadChips) {
688 bool nearBadChip{false};
694 std::map<Identifier, unsigned int>::const_iterator badChip{badChips->find(module_id)};
695 if (badChip != badChips->end()) {
696 status = (*badChip).second;
697
698 const bool nearBadChipDead{(
status & (1 << chipPos)) != 0};
699 const bool nextBadChipDead{(
status & (1 << (chipPos + 1))) != 0};
700 const bool isNotEndChip{(chipPos != 5) and (chipPos != 11)};
701
702
703
704
705 nearBadChip = nearBadChipDead or (isNotEndChip and nextBadChipDead);
706 }
707 if (
failCut(not nearBadChip,
"hit cut: not near bad chip")) {
708 continue;
709 }
710 }
712
715
716 auto effAcc{Monitored::Scalar<float>("eff", eff)};
717 auto ineffAcc{Monitored::Scalar<float>("ineff", (testOffline ? 1. : 1.-eff))};
718 auto ietaAcc{Monitored::Scalar<int>("ieta"+etaPhiSuffix, ieta)};
719 auto iphiAcc{Monitored::Scalar<int>("iphi"+etaPhiSuffix, iphi)};
720 auto layerAcc{Monitored::Scalar<float>("layerPlusHalfSide", dedicated_layerPlusHalfSide)};
721 auto lumiAcc{Monitored::Scalar<int>("LumiBlock", ctx.eventID().lumi_block())};
722 auto isubAcc{Monitored::Scalar<int>("isub", isub)};
723 auto sideHashAcc{Monitored::Scalar<int>("sideHash", sideHash)};
724 auto isFirstBCIDAcc{Monitored::Scalar<bool>("isFirstBCID", (BCIDpos <= 0))};
725
726
727 fill(regionNames[isub].
data(), effAcc, ineffAcc, ietaAcc, iphiAcc, layerAcc, lumiAcc, isFirstBCIDAcc);
728 fill(
"SCTHitEffMonitor", effAcc, lumiAcc, isubAcc, sideHashAcc, isFirstBCIDAcc);
729
730 if (testOffline) {
731 ATH_MSG_INFO(
"Filling " << detIndex <<
", " << side <<
" eta " << ieta <<
" phi " << iphi);
732 }
733 }
734 }
736
737 return StatusCode::SUCCESS;
738}
Scalar theta() const
theta method
void swap(DataVector< T > &a, DataVector< T > &b)
See DataVector<T, BASE>::swap().
char data[hepevt_bytes_allocation_ATLAS]
Environment_t environment() const
Accessor functions for the environment.
std::string print_to_string(Identifier id, const IdContext *context=0) const
or provide the printout in string form
bool is_sct(Identifier id) const
bool is_pixel(Identifier id) const
bool is_trt(Identifier id) const
void getInitializedCache(MagField::AtlasFieldCache &cache) const
get B field cache for evaluation as a function of 2-d or 3-d position.
int distanceFromFront(const bcid_type bcid, const BunchDistanceType type=NanoSec) const
The distance of the specific bunch crossing from the front of the train.
@ BunchCrossings
Distance in units of 25 nanoseconds.
bool is_valid() const
Check if id is in a valid state.
const SiDetectorElement * getDetectorElement(const IdentifierHash &hash) const
bool swapPhiReadoutDirection() const
Determine if readout direction between online and offline needs swapping.
bool solenoidOn() const
status of the magnets
SG::ReadHandleKey< TrackCollection > m_TrackName
int becIdxLayer2Index(const int becIdx, const int layer) const
ToolHandle< Trk::ITrackHoleSearchTool > m_holeSearchTool
int getWaferIndex(const int barrel_bc, const int layer_disk, const int side) const
FloatProperty m_effdistcut
const PixelID * m_pixelId
SG::ReadCondHandleKey< AtlasFieldCacheCondObj > m_fieldCondObjInputKey
SG::ReadHandleKey< ComTime > m_comTimeName
int previousChip(double xl, int side, bool swap) const
SG::ReadHandleKey< InDet::SCT_ClusterContainer > m_sctContainerName
double getResidual(const Identifier &surfaceID, const Trk::TrackParameters *trkParam, const InDet::SCT_ClusterContainer *p_sctclcontainer) const
StatusCode findAnglesToWaferSurface(const Amg::Vector3D &mom, const Identifier id, const InDetDD::SiDetectorElementCollection *elements, double &theta, double &phi) const
SG::ReadCondHandleKey< InDetDD::SiDetectorElementCollection > m_SCTDetEleCollKey
StatusCode failCut(bool value, const std::string &name) const
int layer_disk(const Identifier &id) const
int side(const Identifier &id) const
IdentifierHash wafer_hash(const Identifier &wafer_id) const
wafer hash from id - optimized
virtual int get_id(const IdentifierHash &hash_id, Identifier &id, const IdContext *context=0) const override final
Create compact id from hash id (return == 0 for OK)
int phi_module(const Identifier &id) const
Identifier module_id(int barrel_ec, int layer_disk, int phi_module, int eta_module) const
For a single crystal.
int barrel_ec(const Identifier &id) const
Values of different levels (failure returns 0)
int eta_module(const Identifier &id) const
IdContext wafer_context() const
int get_other_side(const IdentifierHash &id, IdentifierHash &other) const
Wafer hash on other side.
const_pointer_type retrieve()
virtual bool isValid() override final
Can the handle be successfully dereferenced?
const std::string & key() const
Return the StoreGate ID for the referenced object.
const Amg::Vector3D & momentum() const
Access method for the momentum.
double pT() const
Access method for transverse momentum.
@ SiSPSeededFinder
Tracks from SiSPSeedFinder.
@ Measurement
This is a measurement, and will at least contain a Trk::MeasurementBase.
@ Outlier
This TSoS contains an outlier, that is, it contains a MeasurementBase/RIO_OnTrack which was not used ...
@ Hole
A hole on the track - this is defined in the following way.
double chi2(TH1 *h0, TH1 *h1)
unsigned int bec2Index(const int becVal)
Conversion bec->index.
static const std::vector< int > n_layers
ParametersT< TrackParametersDim, Charged, PerigeeSurface > Perigee
ParametersBase< TrackParametersDim, Charged > TrackParameters
@ numberOfSCTHits
number of SCT holes
void fill(H5::Group &out_file, size_t iterations)
1.14.0