ATLAS Offline Software
Loading...
Searching...
No Matches
InDet::TRT_TrackSegmentsMaker_ECcosmics Class Reference

AlgTool that creates TrackSegments out of TRT Driftcircles in the special case of cosmic runs in SR1. More...

#include <TRT_TrackSegmentsMaker_ECcosmics.h>

Inheritance diagram for InDet::TRT_TrackSegmentsMaker_ECcosmics:
Collaboration diagram for InDet::TRT_TrackSegmentsMaker_ECcosmics:

Classes

class  EventData

Public Member Functions

 TRT_TrackSegmentsMaker_ECcosmics (const std::string &, const std::string &, const IInterface *)
 Constructor with parameters.
virtual ~TRT_TrackSegmentsMaker_ECcosmics ()
virtual StatusCode initialize () override
virtual StatusCode finalize () override
virtual std::unique_ptr< InDet::ITRT_TrackSegmentsMaker::IEventDatanewEvent (const EventContext &ctx) const override
virtual std::unique_ptr< InDet::ITRT_TrackSegmentsMaker::IEventDatanewRegion (const EventContext &ctx, const std::vector< IdentifierHash > &) const override
void endEvent (InDet::ITRT_TrackSegmentsMaker::IEventData &event_data) const override
virtual void find (const EventContext &ctx, InDet::ITRT_TrackSegmentsMaker::IEventData &event_data, InDet::TRT_DetElementLink_xk::TRT_DetElemUsedMap &used) const override
virtual Trk::TrackSegmentnext (InDet::ITRT_TrackSegmentsMaker::IEventData &event_data) const override
virtual MsgStream & dump (MsgStream &out) const override
virtual std::ostream & dump (std::ostream &out) const override
ServiceHandle< StoreGateSvc > & evtStore ()
 The standard StoreGateSvc (event store) Returns (kind of) a pointer to the StoreGateSvc.
const ServiceHandle< StoreGateSvc > & detStore () const
 The standard StoreGateSvc/DetectorStore Returns (kind of) a pointer to the StoreGateSvc.
virtual StatusCode sysInitialize () override
 Perform system initialization for an algorithm.
virtual StatusCode sysStart () override
 Handle START transition.
virtual std::vector< Gaudi::DataHandle * > inputHandles () const override
 Return this algorithm's input handles.
virtual std::vector< Gaudi::DataHandle * > outputHandles () const override
 Return this algorithm's output handles.
Gaudi::Details::PropertyBase & declareProperty (Gaudi::Property< T, V, H > &t)
void updateVHKA (Gaudi::Details::PropertyBase &)
MsgStream & msg () const
bool msgLvl (const MSG::Level lvl) const

Static Public Member Functions

static const InterfaceID & interfaceID ()

Protected Member Functions

void retrieveHits (TRT_TrackSegmentsMaker_ECcosmics::EventData &event_data) const
 sort hits into good/noise lists
bool find_seed (int endcap, int zslice, int sector, TRT_TrackSegmentsMaker_ECcosmics::EventData &event_data) const
 Find seed in the given sector/zslice/endcap.
int evaluate_seed (int endcap, int zslice, int sector, const double *p, TRT_TrackSegmentsMaker_ECcosmics::EventData &event_data) const
 Evaluate how many dc match this seed.
void create_segment (std::vector< const InDet::TRT_DriftCircle * > *seed, TRT_TrackSegmentsMaker_ECcosmics::EventData &event_data) const
 Create segment out of a seed.
void setFitFunctions (TRT_TrackSegmentsMaker_ECcosmics::EventData &event_data) const
TF1 * perform_fit (int count, TRT_TrackSegmentsMaker_ECcosmics::EventData &event_data) const
 Perform the fit and return a function that provides the fitted phi information.
bool is_suspicious (const InDet::TRT_DriftCircle *dc, std::vector< const InDet::TRT_DriftCircle * > *seed) const
 checks if a hit that matches the segment looks suspicious (i.e.
bool accepted (const std::list< const InDet::TRT_DriftCircle * >::iterator compareIt, std::list< const InDet::TRT_DriftCircle * > &container, double phiLimit, double dzLimit) const
 is the hit accepted?
void renounceArray (SG::VarHandleKeyArray &handlesArray)
 remove all handles from I/O resolution
std::enable_if_t< std::is_void_v< std::result_of_t< decltype(&T::renounce)(T)> > &&!std::is_base_of_v< SG::VarHandleKeyArray, T > &&std::is_base_of_v< Gaudi::DataHandle, T >, void > renounce (T &h)
void extraDeps_update_handler (Gaudi::Details::PropertyBase &ExtraDeps)
 Add StoreName to extra input/output deps as needed.

Static Protected Member Functions

static double phidiff (double a, double b)
 provide the proper subtraction of two phi values

Protected Attributes

BooleanProperty m_phaseMode
const TRT_IDm_trtid {}
SG::ReadHandleKey< InDet::TRT_DriftCircleContainer > m_trtname {this,"TRT_ClustersContainer","TRT_DriftCircles","RHK to retrieve TRT_DriftCircles"}
 TRTs container.
SG::ReadHandleKey< Trk::PRDtoTrackMapm_prdToTrackMap {this,"PRDtoTrackMap",""}
ToolHandle< Trk::IRIO_OnTrackCreatorm_riomakerD
ToolHandle< Trk::IRIO_OnTrackCreatorm_riomakerN
BooleanProperty m_useDriftTime
DoubleProperty m_scaleTube
DoubleProperty m_scaleFactorDrift
DoubleProperty m_scaleTubeNoise
DoubleProperty m_cutToTLoose
DoubleProperty m_cutToTTight
DoubleProperty m_cutToTUpper
IntegerProperty m_minDCSeed
IntegerProperty m_hitLimit

Static Protected Attributes

static std::mutex s_fitMutex

Private Types

typedef ServiceHandle< StoreGateSvcStoreGateSvc_t

Private Member Functions

Gaudi::Details::PropertyBase & declareGaudiProperty (Gaudi::Property< T, V, H > &hndl, const SG::VarHandleKeyType &)
 specialization for handling Gaudi::Property<SG::VarHandleKey>

Private Attributes

StoreGateSvc_t m_evtStore
 Pointer to StoreGate (event store by default)
StoreGateSvc_t m_detStore
 Pointer to StoreGate (detector store by default)
std::vector< SG::VarHandleKeyArray * > m_vhka
bool m_varHandleArraysDeclared

Detailed Description

AlgTool that creates TrackSegments out of TRT Driftcircles in the special case of cosmic runs in SR1.

Author
Chris.nosp@m.tian.nosp@m..Schm.nosp@m.itt@.nosp@m.cern..nosp@m.ch

Definition at line 64 of file TRT_TrackSegmentsMaker_ECcosmics.h.

Member Typedef Documentation

◆ StoreGateSvc_t

typedef ServiceHandle<StoreGateSvc> AthCommonDataStore< AthCommonMsg< AlgTool > >::StoreGateSvc_t
privateinherited

Definition at line 388 of file AthCommonDataStore.h.

Constructor & Destructor Documentation

◆ TRT_TrackSegmentsMaker_ECcosmics()

InDet::TRT_TrackSegmentsMaker_ECcosmics::TRT_TrackSegmentsMaker_ECcosmics ( const std::string & t,
const std::string & n,
const IInterface * p )

Constructor with parameters.

Definition at line 35 of file TRT_TrackSegmentsMaker_ECcosmics.cxx.

37 : AthAlgTool(t,n,p)
38{
39 declareInterface<ITRT_TrackSegmentsMaker>(this);
40}
AthAlgTool()
Default constructor:

◆ ~TRT_TrackSegmentsMaker_ECcosmics()

InDet::TRT_TrackSegmentsMaker_ECcosmics::~TRT_TrackSegmentsMaker_ECcosmics ( )
virtualdefault

Member Function Documentation

◆ accepted()

bool InDet::TRT_TrackSegmentsMaker_ECcosmics::accepted ( const std::list< const InDet::TRT_DriftCircle * >::iterator compareIt,
std::list< const InDet::TRT_DriftCircle * > & container,
double phiLimit,
double dzLimit ) const
protected

is the hit accepted?

Definition at line 1693 of file TRT_TrackSegmentsMaker_ECcosmics.cxx.

1696 {
1697 //
1698 const auto id = (*it)->identify();
1699 const auto ns = m_trtid->straw(id);
1700 const auto endcap= m_trtid->barrel_ec(id);
1701 const auto wheel = std::abs(m_trtid->layer_or_wheel(id));
1702 const Amg::Vector3D& sc = (*it)->detectorElement()->strawCenter(ns);
1703 //
1704 bool accept{false};
1705 for(auto it2=container.begin();it2!=container.end();++it2){
1706 if(it==it2) continue;
1707 Identifier id2 = (*it2)->identify();
1708 int endcap2=m_trtid->barrel_ec(id2);
1709 int wheel2=std::abs(m_trtid->layer_or_wheel(id2));
1710 //only look inside the same wheel
1711 if(endcap!=endcap2 || wheel2!=wheel) continue;
1712 int ns2 = m_trtid->straw(id2);
1713 const Amg::Vector3D& sc2 = (*it2)->detectorElement()->strawCenter(ns2);
1714 double dz=sc.z()-sc2.z();
1715 double dphi=phidiff(sc.phi(),sc2.phi());
1716 dphi*=300.0; // scale factor for dphi to get equal weight
1717 dz=std::abs(dz);
1718 dphi=std::abs(dphi/300.0);
1719 if(dphi<dPhiLimit and dz<dzLimit){
1720 accept=true;
1721 //break was not in original code...
1722 break; //...but no need to continue once accept has been set
1723 }
1724 } //end inner loop
1725 return accept;
1726}
static Double_t sc
HWIdentifier id2
static double phidiff(double a, double b)
provide the proper subtraction of two phi values
Eigen::Matrix< double, 3, 1 > Vector3D
StatusCode accept(const xAOD::Muon *mu)

◆ create_segment()

void InDet::TRT_TrackSegmentsMaker_ECcosmics::create_segment ( std::vector< const InDet::TRT_DriftCircle * > * seed,
TRT_TrackSegmentsMaker_ECcosmics::EventData & event_data ) const
protected

Create segment out of a seed.

Definition at line 799 of file TRT_TrackSegmentsMaker_ECcosmics.cxx.

801{
802
803 ATH_MSG_VERBOSE("Number of hits in initial seed: "<<seed->size());
804
805 double shift=0.;
806 double phimin=10.;
807 double phimax=-10.;
808
809 double zmin=1e5;
810 double zmax=-1e5;
811
812 std::vector<const InDet::TRT_DriftCircle*>::iterator vit,vitE;
813
814 vit=seed->begin();
815 vitE=seed->end();
816
817 int count=0;
818 for(;vit!=vitE;++vit){
819 int straw = m_trtid->straw((*vit)->identify());
820 const Amg::Vector3D& sc = (*vit)->detectorElement()->strawCenter(straw);
821
822 if(sc.z()<zmin) zmin=sc.z();
823 if(sc.z()>zmax) zmax=sc.z();
824
825 if(count>500){
826 return; // a seed cannot contain that many hits
827 }
828
829 event_data.m_a_z[count]=sc.z();
830 event_data.m_a_phi[count]=sc.phi();
831 event_data.m_a_tan[count]=std::tan(sc.phi());
832 const auto cosPhi{std::cos(sc.phi())};
833 event_data.m_a_tan_err[count]=(1.15/700.)/(cosPhi*cosPhi);
834 event_data.m_a_z_err[count]=1.15; // 4/sqrt(12)
835 event_data.m_a_phi_err[count]=1.15/700.; //take uncertainty at a radius of 700 mm
836
837 if(sc.phi()<phimin)
838 phimin=sc.phi();
839
840 if(sc.phi()>phimax)
841 phimax=sc.phi();
842
843 count++;
844 }
845
846 shift=(phimax+phimin)/2.-M_PI_4;
847 ATH_MSG_VERBOSE("Phimin="<<phimin<<" Phimax="<<phimax<<" shift="<<shift);
848
849 //protection for boundary at pi,-pi
850 if(phimin<-M_PI_2 && phimax>M_PI_2){
851 shift=3.*M_PI_4;
852 ATH_MSG_VERBOSE("Boundary between pi and -pi!!! use the following shift: "<<shift);
853 }
854
855
856 //correct the phi values
857
858 for(int i=0;i<count;i++){
859 double orig=event_data.m_a_phi[i];
860
861 double corr=orig-shift;
862 if (corr > M_PI) corr = std::fmod(corr+M_PI,2.*M_PI)-M_PI;
863 else if(corr <-M_PI) corr = std::fmod(corr-M_PI,2.*M_PI)+M_PI;
864
865 event_data.m_a_phi[i]=corr;
866 event_data.m_a_tan[i]=std::tan(corr);
867 const auto cosCorr{std::cos(corr)};
868 event_data.m_a_tan_err[i]=(1.15/700.)/(cosCorr*cosCorr);
869 }
870
871 TF1 *fit=perform_fit(count, event_data);
872 if(!fit) return;
873
874 std::list<const InDet::TRT_DriftCircle*>::iterator lit,litE;
875
876 //add the ones from the good list (but skip the ones that we already have)
877 lit=event_data.m_goodHits.begin();
878 litE=event_data.m_goodHits.end();
879
880 std::vector<const InDet::TRT_DriftCircle*> tobeadded;
881
882
883 for(;lit!=litE;++lit){
884 int straw = m_trtid->straw((*lit)->identify());
885 const Amg::Vector3D& sc = (*lit)->detectorElement()->strawCenter(straw);
886
887 double fitted_phi=fit->Eval(sc.z(),0.,0.);
888 fitted_phi+=shift;
889
890 if (fitted_phi > M_PI) fitted_phi = std::fmod(fitted_phi+M_PI,2.*M_PI)-M_PI;
891 else if(fitted_phi <-M_PI) fitted_phi = std::fmod(fitted_phi-M_PI,2.*M_PI)+M_PI;
892
893 if(std::abs(phidiff(fitted_phi,sc.phi()))<m_scaleTube*4./800.){
894 //verify that this hit is not already included
895
896 vit=seed->begin();
897 vitE=seed->end();
898 bool new_hit=true;
899 for(;vit!=vitE && new_hit;++vit){
900 if((*vit)==(*lit)){
901 new_hit=false;
902 }
903 }
904 if(new_hit){
905 ATH_MSG_VERBOSE("\t\t good hit matched! fitted_phi="<<fitted_phi<<" z=" <<sc.z());
906
907 if(sc.z()>zmin-200 && sc.z()<zmax+200){
908 if(sc.z()<zmin) zmin=sc.z();
909 if(sc.z()>zmax) zmax=sc.z();
910
911 seed->push_back(*lit);
912 }else{
913 //not sure if this is really a good hit yet ...
914 tobeadded.push_back(*lit);
915 }
916 }
917 }
918 }
919
920 std::sort(tobeadded.begin(),tobeadded.end(), CompareDCz);
921
922 int index1=-1,index2=-1;
923
924 double z1=-1e5;
925
926 count=0;
927 ATH_MSG_VERBOSE("zmin="<<zmin<<" zmax="<<zmax);
928 for(vit=tobeadded.begin();vit!=tobeadded.end();++vit){
929 int straw = m_trtid->straw((*vit)->identify());
930 const Amg::Vector3D& sc = (*vit)->detectorElement()->strawCenter(straw);
931 ATH_MSG_VERBOSE("Shall we add hit at z="<<sc.z());
932
933 if(sc.z()<zmin && sc.z()>z1) {
934 z1=sc.z();
936 }
937 if(sc.z()>zmax && index2==-1) {
939 }
940 count++;
941 }
942
943 ATH_MSG_VERBOSE("index1="<<index1<<" index2="<<index2);
944
945 if(index2>=0){
946 for(int i=index2;i<count;i++){
947 int straw = m_trtid->straw(tobeadded[i]->identify());
948 const Amg::Vector3D& sc = tobeadded[i]->detectorElement()->strawCenter(straw);
949 if(sc.z()<zmax+200){
950 zmax=sc.z();
951 seed->push_back(tobeadded[i]);
952 ATH_MSG_VERBOSE(" added goot hit at z="<<sc.z());
953 }else{
954 ATH_MSG_VERBOSE(" rejected goot hit at z="<<sc.z());
955 i=count;
956 }
957 }
958 }
959 if(index1>=0){
960 for(int i=index1;i>=0;i--){
961 int straw = m_trtid->straw(tobeadded[i]->identify());
962 const Amg::Vector3D& sc = tobeadded[i]->detectorElement()->strawCenter(straw);
963 if(sc.z()>zmin-200){
964 zmin=sc.z();
965 seed->push_back(tobeadded[i]);
966 ATH_MSG_VERBOSE(" added goot hit at z="<<sc.z());
967 }else{
968 ATH_MSG_VERBOSE(" rejected goot hit at z="<<sc.z());
969 i=-1;
970 }
971 }
972 }
973 tobeadded.clear();
974
975
976
977 //fit again
978 std::sort(seed->begin(),seed->end(), CompareDCz);
979
980 vit=seed->begin();
981 vitE=seed->end();
982
983 count=0;
984 for(;vit!=vitE;++vit){
985 int straw = m_trtid->straw((*vit)->identify());
986 const Amg::Vector3D& sc = (*vit)->detectorElement()->strawCenter(straw);
987
988 if(count>500){
989 return; // a seed cannot contain that many hits
990 }
991
992 event_data.m_a_z[count]=sc.z();
993 event_data.m_a_phi[count]=sc.phi();
994 event_data.m_a_tan[count]=std::tan(sc.phi());
995 const auto cosPhi{std::cos(sc.phi())};
996 event_data.m_a_tan_err[count]=(1.15/700.)/(cosPhi*cosPhi);
997 event_data.m_a_z_err[count]=1.15; // 4/sqrt(12)
998 event_data.m_a_phi_err[count]=1.15/700.; //take uncertainty at a radius of 700 mm
999
1000 if(sc.phi()<phimin)
1001 phimin=sc.phi();
1002
1003 if(sc.phi()>phimax)
1004 phimax=sc.phi();
1005
1006 count++;
1007 }
1008
1009 shift=(phimax+phimin)/2.-M_PI_4;
1010 ATH_MSG_VERBOSE("Phimin="<<phimin<<" Phimax="<<phimax<<" shift="<<shift);
1011
1012 //protection for boundary at pi,-pi
1013 if(phimin<-M_PI_2 && phimax>M_PI_2){
1014 shift=3.*M_PI_4;
1015 ATH_MSG_VERBOSE("Boundary between pi and -pi!!! use the following shift: "<<shift);
1016 }
1017
1018 //correct the phi values
1019
1020 for(int i=0;i<count;i++){
1021 double orig=event_data.m_a_phi[i];
1022 double corr=orig-shift;
1023
1024 if (corr > M_PI) corr = std::fmod(corr+M_PI,2.*M_PI)-M_PI;
1025 else if(corr <-M_PI) corr = std::fmod(corr-M_PI,2.*M_PI)+M_PI;
1026
1027 event_data.m_a_phi[i]=corr;
1028 event_data.m_a_tan[i]=std::tan(corr);
1029 const auto cosCorr{std::cos(corr)};
1030 event_data.m_a_tan_err[i]=(1.15/700.)/(cosCorr*cosCorr);
1031 }
1032
1033
1034 fit=perform_fit(count, event_data);
1035 if(!fit) return;
1036
1037 //add the ones from the noise list
1038 lit=event_data.m_noiseHits.begin();
1039 litE=event_data.m_noiseHits.end();
1040
1041 for(;lit!=litE;++lit){
1042
1043 //don't add hits that really look like noise
1044
1045 int straw = m_trtid->straw((*lit)->identify());
1046 const Amg::Vector3D& sc = (*lit)->detectorElement()->strawCenter(straw);
1047
1048 double z=sc.z();
1049 double phi=sc.phi();
1050 double fitted_phi=fit->Eval(sc.z(),0.,0.);
1051 fitted_phi+=shift;
1052
1053 if (fitted_phi > M_PI) fitted_phi = std::fmod(fitted_phi+M_PI,2.*M_PI)-M_PI;
1054 else if(fitted_phi <-M_PI) fitted_phi = std::fmod(fitted_phi-M_PI,2.*M_PI)+M_PI;
1055
1056 //only look into the same endcap for noise hits
1057 if(z*event_data.m_a_z[0]<0) continue;
1058
1059 if(std::abs(phidiff(fitted_phi,phi))<m_scaleTubeNoise*4./800.){
1060 vit=seed->begin();
1061 vitE=seed->end();
1062 bool new_hit=true;
1063 for(;vit!=vitE && new_hit;++vit){
1064 if((*vit)==(*lit)){
1065 new_hit=false;
1066 }
1067 }
1068 if(new_hit){
1069 if(sc.z()>zmin-50 && sc.z()<zmax+50){
1070 ATH_MSG_VERBOSE("\t\t noise hit matched! fitted_phi="<<fitted_phi<<" z=" <<sc.z());
1071 seed->push_back(*lit);
1072 }else{
1073 ATH_MSG_VERBOSE("\t\t noise hit matched but too far away! fitted_phi="<<fitted_phi<<" z=" <<sc.z());
1074 }
1075 }
1076 }
1077 }
1078
1079 //first: sort along y
1080 std::sort(seed->begin(),seed->end(), CompareDCy);
1081
1082 //check z-coordinate of first and last hit
1083 Amg::Vector3D firststraw=((*seed)[0]->detectorElement()->strawCenter( m_trtid->straw((*seed)[0]->identify()) ));
1084 Amg::Vector3D laststraw=((*seed)[seed->size()-1]->detectorElement()->strawCenter( m_trtid->straw((*seed)[seed->size()-1]->identify()) ));
1085 z1=firststraw.z();
1086 double z2=laststraw.z();
1087 double yfirst=firststraw.y();
1088 double ylast=laststraw.y();
1089
1090 ATH_MSG_VERBOSE("zfirst="<<z1<<" zlast="<<z2);
1091
1092 double Theta=0;
1093 int state=0;
1094
1095 //sort again, but this time with the correct ordering
1096 if (std::abs(yfirst-ylast)>50){
1097 if(z1>z2) {
1098 std::sort(seed->begin(),seed->end(), CompareDCzreverse);
1099 }else{
1100 std::sort(seed->begin(),seed->end(), CompareDCz);
1101 }
1102 }
1103 else {
1104 if ((yfirst>0 && z1<0) || (yfirst<0 && z1>0)) std::sort(seed->begin(),seed->end(), CompareDCz);
1105 else std::sort(seed->begin(),seed->end(), CompareDCzreverse);
1106 }
1107 firststraw=(*seed)[0]->detectorElement()->strawCenter( m_trtid->straw((*seed)[0]->identify()) );
1108 laststraw=(*seed)[seed->size()-1]->detectorElement()->strawCenter( m_trtid->straw((*seed)[seed->size()-1]->identify()) );
1109 z1=firststraw.z();
1110 z2=laststraw.z();
1111
1112
1113 if(z1>z2) {
1114 Theta=M_PI-std::atan2(400.,std::abs(z1-z2));
1115 state=2;
1116 }else{
1117 Theta=std::atan2(400.,std::abs(z1-z2));
1118 state=1;
1119 }
1120
1121 double globaly=((*seed)[0]->detectorElement()->strawCenter( m_trtid->straw((*seed)[0]->identify()) )).y();
1122 double firstphi=((*seed)[0]->detectorElement()->strawCenter( m_trtid->straw((*seed)[0]->identify()) )).phi();
1123
1124
1125 ATH_MSG_VERBOSE("Number of tube hits matching straight line: "<<seed->size());
1126
1127
1128 //estimate Segment parameters from first and last hit
1129
1130 const Amg::Vector3D& sc_first = ((*seed)[0]->detectorElement()->strawCenter( m_trtid->straw((*seed)[0]->identify()) ));
1131 const Amg::Vector3D& sc_last = ((*seed)[seed->size()-1]->detectorElement()->strawCenter( m_trtid->straw((*seed)[seed->size()-1]->identify()) ));
1132
1133
1134 double tx1,tx2,ty1,ty2;
1135
1136 tx1=1050*std::cos(sc_first.phi());
1137 ty1=1050*std::sin(sc_first.phi());
1138
1139 tx2=650*std::cos(sc_last.phi());
1140 ty2=650*std::sin(sc_last.phi());
1141
1142 double tphi=std::atan2(ty2-ty1,tx2-tx1);
1143
1144
1145 ATH_MSG_VERBOSE("First hit: "<<sc_first<<" \nLast hit: "<<sc_last);
1146
1147
1148 ATH_MSG_VERBOSE(" (x1,y1) = ("<<tx1<<","<<ty1<<")");
1149 ATH_MSG_VERBOSE(" (x2,y2) = ("<<tx2<<","<<ty2<<")");
1150
1151
1152
1153 //correct phi to point downwards
1154 if(tphi>0)
1155 tphi-=M_PI;
1156
1157 //make sure that Theta is positive and less than Pi!
1158
1159 while(Theta<0)
1160 Theta+=M_PI;
1161
1162 while(Theta>M_PI)
1163 Theta-=M_PI;
1164
1165 ATH_MSG_VERBOSE(" tphi="<<tphi);
1166
1167 ATH_MSG_VERBOSE("globaly = "<<globaly<<" phi="<<firstphi<<" theta="<<Theta<<" --> state = "<<state);
1168
1169
1170 vit=seed->begin();
1171 vitE=seed->end();
1172
1173 count=0;
1174 for(;vit!=vitE;++vit){
1175 count++;
1176 int straw = m_trtid->straw((*vit)->identify());
1177 const Amg::Vector3D& sc = (*vit)->detectorElement()->strawCenter(straw);
1178 ATH_MSG_VERBOSE("Hit "<<count<<" z="<<sc.z()<<" phi="<<sc.phi()<<" y="<<sc.y());
1179 }
1180
1181
1182 if(m_useDriftTime){
1183
1184 ATH_MSG_VERBOSE("drifttime should be used!!!");
1185
1186 //refit the seed but this time including the drift time information
1187
1188 //take the left-rigfht assignment from the tube fit and then fit again
1189 count=0;
1190
1191 vit=seed->begin();
1192 vitE=seed->end();
1193
1194 for(;vit!=vitE;++vit){
1195 int straw = m_trtid->straw((*vit)->identify());
1196 const Amg::Vector3D& sc = (*vit)->detectorElement()->strawCenter(straw);
1197
1198 double fitted_phi=fit->Eval(sc.z(),0.,0.);
1199 fitted_phi+=shift;
1200
1201 if (fitted_phi > M_PI) fitted_phi = std::fmod(fitted_phi+M_PI,2.*M_PI)-M_PI;
1202 else if(fitted_phi <-M_PI) fitted_phi = std::fmod(fitted_phi-M_PI,2.*M_PI)+M_PI;
1203
1204 int sign=1;
1205 if(fitted_phi<sc.phi()) sign=-1;
1206
1207
1208 double driftr = (*vit)->localPosition()(Trk::driftRadius);
1209 double drifte = 2*std::sqrt((*vit)->localCovariance()(0,0));
1210
1211 double straw_r=640+400./std::abs(z2-z1)*std::abs(sc.z()-z1);
1212
1213 double dphi=sign*driftr/straw_r;
1214
1215
1216 event_data.m_a_z[count]=sc.z();
1217 event_data.m_a_phi[count]=sc.phi()+dphi;
1218 event_data.m_a_phi_err[count]=drifte/straw_r;
1219
1220 event_data.m_a_tan[count]=std::tan(event_data.m_a_phi[count]-shift);
1221 const auto cosPhiMinShift{std::cos(event_data.m_a_phi[count]-shift)};
1222 event_data.m_a_tan_err[count]=(event_data.m_a_phi_err[count])/( cosPhiMinShift
1223 * cosPhiMinShift);
1224
1225 event_data.m_a_z_err[count]=0.;
1226
1227 ATH_MSG_VERBOSE(count<<" z="<<event_data.m_a_z[count]<<" phi="<<event_data.m_a_phi[count]<<" +- "<< event_data.m_a_phi_err[count]);
1228 count++;
1229 }
1230
1231 ATH_MSG_VERBOSE("Number of hits in first fit (driftradius): "<<count);
1232
1233 fit=perform_fit(count,event_data);
1234 if(!fit) return;
1235
1236
1237 //check if hits are matching
1238 vit=seed->begin();
1239 vitE=seed->end();
1240
1241 count=0;
1242 for(;vit!=vitE;++vit){
1243 int straw = m_trtid->straw((*vit)->identify());
1244 const Amg::Vector3D& sc = (*vit)->detectorElement()->strawCenter(straw);
1245 double fitted_phi=fit->Eval(sc.z(),0.,0.);
1246 fitted_phi+=shift;
1247
1248 double driftr = (*vit)->localPosition()(Trk::driftRadius);
1249 double drifte = 2*std::sqrt((*vit)->localCovariance()(0,0));
1250 int sign=1;
1251
1252 if(fitted_phi<sc.phi()) sign=-1;
1253
1254 double straw_r=640+400./std::abs(z2-z1)*std::abs(sc.z()-z1);
1255 double dphi=sign*driftr/straw_r;
1256
1257 if(count>500){
1258 return; // a seed cannot contain that many hits
1259 }
1260
1261 if (msgLvl(MSG::VERBOSE)) msg()<<"z="<<sc.z()<<" phi="<<sc.phi()<<" fitted_phi = "<<fitted_phi<<endmsg;
1262 if(std::abs(phidiff(fitted_phi,sc.phi()-dphi)) < m_scaleFactorDrift* drifte/straw_r){
1263 if (msgLvl(MSG::VERBOSE)) msg()<<"\t\t matched!!!"<<endmsg;
1264
1265 event_data.m_a_z[count]=sc.z();
1266 event_data.m_a_phi[count]=sc.phi()+dphi;
1267 event_data.m_a_phi_err[count]=drifte/straw_r;
1268
1269 event_data.m_a_tan[count]=std::tan(event_data.m_a_phi[count]-shift);
1270 const auto cosPhiMinShift{ std::cos(event_data.m_a_phi[count]-shift)};
1271 event_data.m_a_tan_err[count]=(event_data.m_a_phi_err[count])/( cosPhiMinShift
1272 *cosPhiMinShift);
1273
1274 event_data.m_a_z_err[count]=0.;
1275 count++;
1276
1277 }
1278 }
1279
1280 //refit
1281
1282 ATH_MSG_VERBOSE("Number of hits in second fit (driftradius): "<<count);
1283
1284 if(count>4){
1285 fit=perform_fit(count,event_data);
1286 if(!fit) return;
1287 }
1288
1289 vit=seed->begin();
1290 vitE=seed->end();
1291
1292 count=0;
1293 for(;vit!=vitE;++vit){
1294 int straw = m_trtid->straw((*vit)->identify());
1295 const Amg::Vector3D& sc = (*vit)->detectorElement()->strawCenter(straw);
1296
1297 double fitted_phi=fit->Eval(sc.z(),0.,0.);
1298 fitted_phi+=shift;
1299
1300 if (fitted_phi > M_PI) fitted_phi = std::fmod(fitted_phi+M_PI,2.*M_PI)-M_PI;
1301 else if(fitted_phi <-M_PI) fitted_phi = std::fmod(fitted_phi-M_PI,2.*M_PI)+M_PI;
1302
1303 double driftr = (*vit)->localPosition()(Trk::driftRadius);
1304 double drifte = 2*std::sqrt((*vit)->localCovariance()(0,0));
1305 int sign=1;
1306
1307 if(fitted_phi<sc.phi()) sign=-1;
1308
1309 double straw_r=640+400./std::abs(z2-z1)*std::abs(sc.z()-z1);
1310 double dphi=sign*driftr/straw_r;
1311
1312
1313 ATH_MSG_VERBOSE("z="<<sc.z()<<" phi="<<sc.phi()<<" fitted_phi = "<<fitted_phi);
1314 if(std::abs(phidiff(fitted_phi,sc.phi()-dphi)) < m_scaleFactorDrift* drifte/straw_r){
1315 ATH_MSG_VERBOSE("\t\t matched!!!");
1316 count++;
1317 }
1318 }
1319
1320 ATH_MSG_VERBOSE("Number of hits matched after second fit: "<<count);
1321
1322 }
1323
1324
1325 //create the segment
1326
1327
1328
1329 if(count<m_minDCSeed) return; //too short segments don't make sense
1330
1331 if(count>=400) return; //too long segments don't make sense
1332
1333
1334
1335 // RIOonTrack production
1336 //
1337 auto rio = DataVector<const Trk::MeasurementBase>{};
1338
1339 vit=seed->begin();
1340 vitE=seed->end();
1341
1342 double chi2=0.;
1343
1344 double initial_r=-10.;
1345 double initial_locz=0.;
1346
1347 count=0;
1348 for(;vit!=vitE;++vit){
1349 const Trk::StraightLineSurface* line = static_cast<const Trk::StraightLineSurface*>
1350 (&((*vit)->detectorElement())->surface( (*vit)->identify()));
1351
1352
1353 int straw = m_trtid->straw((*vit)->identify());
1354 const Amg::Vector3D& sc = (*vit)->detectorElement()->strawCenter(straw);
1355
1356 double z=sc.z();
1357 // double phi=sc.phi();
1358
1359 double locz=0.;
1360
1361 if(std::abs(z2-z1)>0.000001){
1362 if(state==1){
1363 locz=400./std::abs(z2-z1)*std::abs(z-z1)-200.;
1364 }else{
1365 locz=200.-400./std::abs(z2-z1)*std::abs(z-z1);
1366 }
1367 }
1368 double fitted_phi=fit->Eval(z,0.,0.);
1369 fitted_phi+=shift;
1370
1371 if (fitted_phi > M_PI) fitted_phi = std::fmod(fitted_phi+M_PI,2.*M_PI)-M_PI;
1372 else if(fitted_phi <-M_PI) fitted_phi = std::fmod(fitted_phi-M_PI,2.*M_PI)+M_PI;
1373
1374 double fitted_r=(phidiff(fitted_phi,sc.phi()))*(840+locz);
1375 double temp_phi=fitted_phi;
1376
1377 if(initial_r==-10. && std::abs(fitted_r)<4.0){
1378 initial_r=fitted_r;
1379 initial_locz=locz;
1380 }
1381
1382 if (temp_phi > M_PI) temp_phi = std::fmod(temp_phi+M_PI,2.*M_PI)-M_PI;
1383 else if(temp_phi <-M_PI) temp_phi = std::fmod(temp_phi-M_PI,2.*M_PI)+M_PI;
1384
1385 ATH_MSG_VERBOSE(count<<" fitted_r = "<<fitted_r<<" locz="<<locz<<" temp_phi="<<temp_phi);
1386
1387
1388 fitted_r*=-1.0;
1389
1390
1391 temp_phi=firstphi;
1392
1393 const Trk::AtaStraightLine Tp(fitted_r,locz,tphi,Theta,.000000001,*line);
1394
1395 //decide if drifttime should be used
1396
1397
1398
1399
1400 bool useDrift=false;
1401
1402 if(m_useDriftTime){
1403 double driftr = (*vit)->localPosition()(Trk::driftRadius);
1404 double drifte = 2*std::sqrt((*vit)->localCovariance()(0,0));
1405 int sign=1;
1406
1407 if(fitted_phi<sc.phi()) sign=-1;
1408
1409 double straw_r=640+400./std::abs(z2-z1)*std::abs(sc.z()-z1);
1410 double dphi=sign*driftr/straw_r;
1411
1412 if(std::abs(phidiff(fitted_phi,sc.phi()-dphi)) < 3* drifte/straw_r){
1413 useDrift=true;
1414 chi2+=(std::abs(fitted_r)-std::abs(driftr))*(std::abs(fitted_r)-std::abs(driftr))/(drifte*drifte);
1415 }
1416
1417 }
1418
1419
1420 if(useDrift){
1421 ATH_MSG_VERBOSE("RIO using drift time!");
1422 rio.push_back(m_riomakerD->correct(*(*vit),Tp, Gaudi::Hive::currentContext()));
1423 }else{
1424 ATH_MSG_VERBOSE("RIO without using drift time!");
1425 chi2+=(fitted_r/1.15)*(fitted_r/1.15); // no drift time used
1426 ATH_MSG_VERBOSE(count<<"\t\t chi2 contribution: "<<(fitted_r/1.15)*(fitted_r/1.15));
1427 rio.push_back(m_riomakerN->correct(*(*vit),Tp, Gaudi::Hive::currentContext()));
1428 }
1429
1430 count++;
1431
1432 //add a pseudomeasurement for the first and the last hit constraining loc_z
1433 if(count==1 || (size_t)count==seed->size()){
1435 Trk::LocalParameters par(dp);
1436 Amg::MatrixX cov(1,1);
1437 cov<<1.;
1438
1439 //create new surface
1440 Amg::Vector3D C = line->transform().translation();
1441
1442 if(state==2 || state==3){
1443 C[2]-=1.;
1444 }else{
1445 C[2]+=1.;
1446 }
1447 Amg::Transform3D T = Amg::Transform3D(line->transform().rotation()*Amg::Translation3D(C));
1448 Trk::StraightLineSurface* surface = new Trk::StraightLineSurface(T);
1449
1450
1451 Trk::PseudoMeasurementOnTrack *pseudo=new Trk::PseudoMeasurementOnTrack(
1452 std::move(par),
1453 std::move(cov),
1454 *surface);
1455
1456 delete surface;
1457
1458 rio.push_back(pseudo);
1459
1460 }
1461 }
1462
1463
1464 //Track Segment production
1465
1466 Trk::FitQuality * fqu = new Trk::FitQuality(chi2,count-2);
1467 const Trk::Surface * sur = &((*seed)[0]->detectorElement())->surface((*seed)[0]->identify());
1468
1469
1470 //phi should always be negative!
1471
1472 if(firstphi>0) firstphi=-firstphi;
1473
1474
1475 Trk::DefinedParameter dp0(initial_r,Trk::locR);
1476 Trk::DefinedParameter dp1(initial_locz,Trk::locZ);
1480
1481 std::array<Trk::DefinedParameter,5> defPar = {dp0,dp1,dp2,dp3,dp4};
1482
1483 Trk::LocalParameters par(defPar);
1484
1485
1486 ATH_MSG_VERBOSE("Track parameters from segment:"<<par);
1487
1488 double universal=1.15;
1489
1490
1491 Amg::MatrixX cov(5,5);
1492 cov<<
1493 universal*universal, 0., 0., 0., 0.,
1494 0. , 160000./12., 0., 0., 0.,
1495 0. , 0., 0.1, 0., 0.,
1496 0. , 0., 0., 1., 0.,
1497 0. , 0., 0., 0., 1.;
1498 Trk::TrackSegment* segment = new Trk::TrackSegment(
1499 std::move(par), std::move(cov), sur, std::move(rio), fqu, Trk::Segment::TRT_SegmentMaker);
1500
1501 //add segment to list of segments
1502
1503 ATH_MSG_VERBOSE("Add segment to list");
1504 event_data.m_segments.push_back(segment);
1505
1506
1507}
#define M_PI
Scalar phi() const
phi method
#define endmsg
#define ATH_MSG_VERBOSE(x)
static Double_t Tp(Double_t *t, Double_t *par)
int sign(int a)
#define z
bool msgLvl(const MSG::Level lvl) const
MsgStream & msg() const
ToolHandle< Trk::IRIO_OnTrackCreator > m_riomakerN
ToolHandle< Trk::IRIO_OnTrackCreator > m_riomakerD
TF1 * perform_fit(int count, TRT_TrackSegmentsMaker_ECcosmics::EventData &event_data) const
Perform the fit and return a function that provides the fitted phi information.
double chi2(TH1 *h0, TH1 *h1)
int count(std::string s, const std::string &regx)
count how many occurances of a regx are in a string
Definition hcg.cxx:146
struct color C
Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic > MatrixX
Dynamic Matrix - dynamic allocation.
Eigen::Affine3d Transform3D
Eigen::Translation< double, 3 > Translation3D
bool CompareDCy(const InDet::TRT_DriftCircle *x, const InDet::TRT_DriftCircle *y)
bool CompareDCzreverse(const InDet::TRT_DriftCircle *x, const InDet::TRT_DriftCircle *y)
bool CompareDCz(const InDet::TRT_DriftCircle *x, const InDet::TRT_DriftCircle *y)
@ straw
Definition HitInfo.h:82
unsigned long long T
ParametersT< TrackParametersDim, Charged, StraightLineSurface > AtaStraightLine
@ driftRadius
trt, straws
Definition ParamDefs.h:53
@ locR
Definition ParamDefs.h:44
@ theta
Definition ParamDefs.h:66
@ qOverP
perigee
Definition ParamDefs.h:67
@ phi
Definition ParamDefs.h:75
@ locZ
local cylindrical
Definition ParamDefs.h:42
std::pair< double, ParamDefs > DefinedParameter
Typedef to of a std::pair<double, ParamDefs> to identify a passed-through double as a specific type o...
void sort(typename DataModel_detail::iterator< DVL > beg, typename DataModel_detail::iterator< DVL > end)
Specialization of sort for DataVector/List.
const Identifier & identify(const UncalibratedMeasurement *meas)
Returns the associated identifier from the muon measurement.

◆ declareGaudiProperty()

Gaudi::Details::PropertyBase & AthCommonDataStore< AthCommonMsg< AlgTool > >::declareGaudiProperty ( Gaudi::Property< T, V, H > & hndl,
const SG::VarHandleKeyType &  )
inlineprivateinherited

specialization for handling Gaudi::Property<SG::VarHandleKey>

Definition at line 156 of file AthCommonDataStore.h.

158 {
160 hndl.value(),
161 hndl.documentation());
162
163 }
Gaudi::Details::PropertyBase & declareProperty(Gaudi::Property< T, V, H > &t)

◆ declareProperty()

Gaudi::Details::PropertyBase & AthCommonDataStore< AthCommonMsg< AlgTool > >::declareProperty ( Gaudi::Property< T, V, H > & t)
inlineinherited

Definition at line 145 of file AthCommonDataStore.h.

145 {
146 typedef typename SG::HandleClassifier<T>::type htype;
148 }
Gaudi::Details::PropertyBase & declareGaudiProperty(Gaudi::Property< T, V, H > &hndl, const SG::VarHandleKeyType &)
specialization for handling Gaudi::Property<SG::VarHandleKey>

◆ detStore()

const ServiceHandle< StoreGateSvc > & AthCommonDataStore< AthCommonMsg< AlgTool > >::detStore ( ) const
inlineinherited

The standard StoreGateSvc/DetectorStore Returns (kind of) a pointer to the StoreGateSvc.

Definition at line 95 of file AthCommonDataStore.h.

◆ dump() [1/2]

MsgStream & InDet::TRT_TrackSegmentsMaker_ECcosmics::dump ( MsgStream & out) const
overridevirtual

Implements InDet::ITRT_TrackSegmentsMaker.

Definition at line 1528 of file TRT_TrackSegmentsMaker_ECcosmics.cxx.

1529{
1530 return out;
1531}

◆ dump() [2/2]

std::ostream & InDet::TRT_TrackSegmentsMaker_ECcosmics::dump ( std::ostream & out) const
overridevirtual

Implements InDet::ITRT_TrackSegmentsMaker.

Definition at line 1538 of file TRT_TrackSegmentsMaker_ECcosmics.cxx.

1539{
1540 return out;
1541}

◆ endEvent()

void InDet::TRT_TrackSegmentsMaker_ECcosmics::endEvent ( InDet::ITRT_TrackSegmentsMaker::IEventData & event_data) const
overridevirtual

Implements InDet::ITRT_TrackSegmentsMaker.

Definition at line 248 of file TRT_TrackSegmentsMaker_ECcosmics.cxx.

249{
250}

◆ evaluate_seed()

int InDet::TRT_TrackSegmentsMaker_ECcosmics::evaluate_seed ( int endcap,
int zslice,
int sector,
const double * p,
TRT_TrackSegmentsMaker_ECcosmics::EventData & event_data ) const
protected

Evaluate how many dc match this seed.

Definition at line 744 of file TRT_TrackSegmentsMaker_ECcosmics.cxx.

749{
750
751 std::list<const InDet::TRT_DriftCircle*>::iterator it,itE;
752
753 int count=0;
754
755 //look at all nearest neighbours ...
756
757 for(int zsl=zslice-1;zsl<=zslice+1;zsl++){
758 for(int ps=sector-1;ps<=sector+1;ps++){
759
760 //make sure that the slices exist ...
761 if(zsl<0 || zsl>=20) continue;
762
763 int ps2=ps;
764 if(ps2<0) ps2+=16;
765 else if(ps2>15) ps2-=16;
766
767
768 itE=event_data.m_sectors[endcap][zsl][ps2].end();
769 for(it=event_data.m_sectors[endcap][zsl][ps2].begin();it!=itE;++it){
770 int straw = m_trtid->straw((*it)->identify());
771 const Amg::Vector3D& sc = (*it)->detectorElement()->strawCenter(straw);
772
773 double pred_phi=sc.z()*p[0]+p[1];
774 double corr=sc.phi()-p[2];
775
776 if(p[0]==1e10){
777 if(sc.z()==p[1]) {
778 pred_phi=corr; //z=const ...
779 }else{
780 pred_phi=corr+0.5;
781 }
782 }
783 if (corr > M_PI) corr = std::fmod(corr+M_PI,2.*M_PI)-M_PI;
784 else if(corr <-M_PI) corr = std::fmod(corr-M_PI,2.*M_PI)+M_PI;
785
786
787 if(std::abs(phidiff(pred_phi,corr))<=4./650.){
788 count=count+1;
789 }
790 }
791 }
792 }
793
794
795 return count;
796}

◆ evtStore()

ServiceHandle< StoreGateSvc > & AthCommonDataStore< AthCommonMsg< AlgTool > >::evtStore ( )
inlineinherited

The standard StoreGateSvc (event store) Returns (kind of) a pointer to the StoreGateSvc.

Definition at line 85 of file AthCommonDataStore.h.

◆ extraDeps_update_handler()

void AthCommonDataStore< AthCommonMsg< AlgTool > >::extraDeps_update_handler ( Gaudi::Details::PropertyBase & ExtraDeps)
protectedinherited

Add StoreName to extra input/output deps as needed.

use the logic of the VarHandleKey to parse the DataObjID keys supplied via the ExtraInputs and ExtraOuputs Properties to add the StoreName if it's not explicitly given

◆ finalize()

StatusCode InDet::TRT_TrackSegmentsMaker_ECcosmics::finalize ( )
overridevirtual

Implements InDet::ITRT_TrackSegmentsMaker.

Definition at line 145 of file TRT_TrackSegmentsMaker_ECcosmics.cxx.

146{
147 StatusCode sc = AthAlgTool::finalize();
148 return sc;
149}
::StatusCode StatusCode
StatusCode definition for legacy code.

◆ find()

void InDet::TRT_TrackSegmentsMaker_ECcosmics::find ( const EventContext & ctx,
InDet::ITRT_TrackSegmentsMaker::IEventData & event_data,
InDet::TRT_DetElementLink_xk::TRT_DetElemUsedMap & used ) const
overridevirtual

Implements InDet::ITRT_TrackSegmentsMaker.

Definition at line 302 of file TRT_TrackSegmentsMaker_ECcosmics.cxx.

305{
306 TRT_TrackSegmentsMaker_ECcosmics::EventData &
308
309
310
311 //loop over both endcaps
312 for(int endcap=0;endcap<2;endcap++){
313
314 // 1st step: find sector with the maximal number of hits until no one is left
315 // with at least 5 hits
316
317 int max=m_minDCSeed;
318
319 while(max>=m_minDCSeed){
320 max=-1;
321 int sector=-1;
322 int zslice=-1;
323
324 for(int j=0;j<20;j++){
325 for(int i=0;i<16;i++){
326 int count=event_data.m_sectors[endcap][j][i].size();
327
328 if(count>0)
329 ATH_MSG_VERBOSE("Endcap "<<endcap<<" zslice "<<j<<" sector "<<i<<" : "<<count);
330
331 if(count>max){
332 max=count;
333 sector=i;
334 zslice=j;
335 }
336 }
337 }
338
339 if(sector<0) continue; // no hits in this endcap ...
340
341 max=event_data.m_sectors[endcap][zslice][sector].size();
342
343 //check if enough driftcircles present
344 if(max>=3){
345 ATH_MSG_VERBOSE("Number of DCs in seed sector: "<<max);
346
347 //2nd step: find seed among those driftcircles
348
349 bool found=find_seed(endcap,zslice,sector,event_data);
350
351 //if no seed is found we clear the sector
352 if(!found){
353 event_data.m_sectors[endcap][zslice][sector].clear();
354 }
355 }
356
357 }
358 }
359
360 ATH_MSG_DEBUG("Found "<<event_data.m_seeds.size()<<" initial seeds !");
361
362
363 //now loop over all seeds and try to create segments out of them
364
365 std::vector<std::vector<const InDet::TRT_DriftCircle*> *>::iterator sit,sitE;
366
367 sit=event_data.m_seeds.begin();
368 sitE=event_data.m_seeds.end();
369
370 for(;sit!=sitE;++sit){
371 create_segment(*sit, event_data);
372 }
373
374
375 ATH_MSG_DEBUG("Found "<<event_data.m_segments.size()<<" TRT Segments");
376
377 event_data.m_segiterator = event_data.m_segments.begin();
378}
#define ATH_MSG_DEBUG(x)
#define max(a, b)
Definition cfImp.cxx:41
void create_segment(std::vector< const InDet::TRT_DriftCircle * > *seed, TRT_TrackSegmentsMaker_ECcosmics::EventData &event_data) const
Create segment out of a seed.
bool find_seed(int endcap, int zslice, int sector, TRT_TrackSegmentsMaker_ECcosmics::EventData &event_data) const
Find seed in the given sector/zslice/endcap.
static EventData & getPrivateEventData(InDet::ITRT_TrackSegmentsMaker::IEventData &virt_event_data)

◆ find_seed()

bool InDet::TRT_TrackSegmentsMaker_ECcosmics::find_seed ( int endcap,
int zslice,
int sector,
TRT_TrackSegmentsMaker_ECcosmics::EventData & event_data ) const
protected

Find seed in the given sector/zslice/endcap.

Definition at line 386 of file TRT_TrackSegmentsMaker_ECcosmics.cxx.

390{
391 double shift=0.;
392 double phimin=10.;
393 double phimax=-10.;
394
395 std::vector<const InDet::TRT_DriftCircle*> dc;
396 std::list<const InDet::TRT_DriftCircle*>::iterator it,itE;
397
398 dc.reserve(event_data.m_sectors[endcap][zslice][sector].size());
399 itE=event_data.m_sectors[endcap][zslice][sector].end();
400 for(it=event_data.m_sectors[endcap][zslice][sector].begin();it!=itE;++it){
401 dc.push_back(*it);
402 }
403 //sort them along z
404
405 std::sort(dc.begin(),dc.end(), CompareDCz);
406
407 std::vector<const InDet::TRT_DriftCircle*>::iterator vit,vitE;
408 vit=dc.begin();
409 vitE=dc.end();
410
411 ATH_MSG_VERBOSE("Sorted driftcircles: "<<dc.size());
412
413 for(;vit!=vitE;++vit){
414 int straw = m_trtid->straw((*vit)->identify());
415 const Amg::Vector3D& sc = (*vit)->detectorElement()->strawCenter(straw);
416
417 ATH_MSG_VERBOSE("z="<<sc.z()<<" phi="<<sc.phi());
418
419 }
420
421
422 //sanity check: if more than 40 good hits per seed -> not a pure cosmic event!
423
424 if(dc.size()>=40) return false;
425
426 //fill arrays for the fit
427
428 vit=dc.begin();
429 vitE=dc.end();
430
431
432 for(;vit!=vitE;++vit){
433 int straw = m_trtid->straw((*vit)->identify());
434 const Amg::Vector3D& sc = (*vit)->detectorElement()->strawCenter(straw);
435
436 if(sc.phi()<phimin)
437 phimin=sc.phi();
438
439 if(sc.phi()>phimax)
440 phimax=sc.phi();
441 }
442
443 shift=(phimax+phimin)/2.-M_PI_4;
444 ATH_MSG_VERBOSE("Phimin="<<phimin<<" Phimax="<<phimax<<" -> shift = "<<shift);
445
446 //protection for boundary at pi,-pi
447 if(phimin<-M_PI_2 && phimax>M_PI_2){
448 shift=3.*M_PI_4;
449 ATH_MSG_VERBOSE("Boundary between pi and -pi!!! use the following shift: "<<shift);
450 }
451
452
453 int count=0;
454 for(vit=dc.begin();vit!=vitE;++vit){
455 int straw = m_trtid->straw((*vit)->identify());
456 const Amg::Vector3D& sc = (*vit)->detectorElement()->strawCenter(straw);
457
458 double orig=sc.phi();
459
460 double corr=orig-shift;
461
462 if (corr > M_PI) corr = std::fmod(corr+M_PI,2.*M_PI)-M_PI;
463 else if(corr <-M_PI) corr = std::fmod(corr-M_PI,2.*M_PI)+M_PI;
464
465 if(count>500){
466 return false; // a seed cannot contain that many hits
467 }
468
469 event_data.m_a_z[count]=sc.z();
470 event_data.m_a_phi[count]=corr;
471 count++;
472 }
473
474 int maxdc=0;
475 double p_best[3]={0,0,0};
476 p_best[2]=shift;
477
478 for(int i=0;i<count;i++){
479 for(int j=i+1;j<count;j++){
480
481 double p[3];
482
483 p[2]=shift;
484
485 if(event_data.m_a_z[i]==event_data.m_a_z[j]){
486 p[0]=1e10;
487 p[1]=event_data.m_a_z[i];
488 }else{
489 p[0]=(event_data.m_a_phi[i]-event_data.m_a_phi[j])/(event_data.m_a_z[i]-event_data.m_a_z[j]);
490 p[1]=event_data.m_a_phi[j]-p[0]*event_data.m_a_z[j];
491 }
492
493
494
495 int match=evaluate_seed(endcap,zslice,sector,p,event_data);
496
497 if(match>maxdc){
498 maxdc=match;
499 p_best[0]=p[0];
500 p_best[1]=p[1];
501 }
502 }
503 }
504
505 ATH_MSG_VERBOSE("Number of matching dc from best seed: "<<maxdc);
506
507
508 if(maxdc>=4){
509 //ok, we've found a seed
510
511 double zmin=1e4;
512 double zmax=-1e4;
513
514
515 int muster[3][3];
516 for(auto & i : muster)
517 for(int & j : i)
518 j=0;
519
520 std::vector<const InDet::TRT_DriftCircle*> *seed=new std::vector<const InDet::TRT_DriftCircle*>;
521
522 for(int zsl=zslice-1;zsl<=zslice+1;zsl++){
523 for(int ps=sector-1;ps<=sector+1;ps++){
524
525 //make sure that the slices exist ...
526 if(zsl<0 || zsl>=20) continue;
527
528 int ps2=ps;
529 if(ps2<0) ps2+=16;
530 else if(ps2>15) ps2-=16;
531
532 for(it=event_data.m_sectors[endcap][zsl][ps2].begin();it!=event_data.m_sectors[endcap][zsl][ps2].end();++it){
533 int straw = m_trtid->straw((*it)->identify());
534 const Amg::Vector3D& sc = (*it)->detectorElement()->strawCenter(straw);
535
536 double pred_phi=sc.z()*p_best[0]+p_best[1];
537 double corr=sc.phi()-p_best[2];
538
539 if(p_best[0]==1e10){
540 if(sc.z()==p_best[1]) {
541 pred_phi=corr; //z=const ...
542 }else{
543 pred_phi=corr+0.5;
544 }
545 }
546 if (corr > M_PI) corr = std::fmod(corr+M_PI,2.*M_PI)-M_PI;
547 else if(corr <-M_PI) corr = std::fmod(corr-M_PI,2.*M_PI)+M_PI;
548
549
550 if(std::abs(phidiff(pred_phi,corr))<=4./650.){
551 seed->push_back(*it);
552 if(sc.z()<zmin) zmin=sc.z();
553 if(sc.z()>zmax) zmax=sc.z();
554 muster[zsl-zslice+1][ps-sector+1]=1;
555 }
556 }
557 }
558 }
559
560 std::vector<const InDet::TRT_DriftCircle*>::iterator vit,vitE;
561
562
563 //do we need to extend the range ??
564
565 ATH_MSG_VERBOSE("Pattern:");
566 ATH_MSG_VERBOSE("\t"<<muster[0][0]<<" "<<muster[1][0]<<" "<<muster[2][0]);
567 ATH_MSG_VERBOSE("\t"<<muster[0][1]<<" "<<muster[1][1]<<" "<<muster[2][1]);
568 ATH_MSG_VERBOSE("\t"<<muster[0][2]<<" "<<muster[1][2]<<" "<<muster[2][2]);
569
570 int zsl=-1;
571 if(muster[2][0] || muster[2][1] || muster[2][2]){
572 //extend towards higher z-values
573 zsl=zslice+2;
574
575 for(int ps=sector-2;ps<=sector+2;ps++){
576
577 //make sure that the slices exist ...
578 if(zsl<0 || zsl>=20) continue;
579
580 int ps2=ps;
581 if(ps2<0) ps2+=16;
582 else if(ps2>15) ps2-=16;
583
584 for(it=event_data.m_sectors[endcap][zsl][ps2].begin();it!=event_data.m_sectors[endcap][zsl][ps2].end();++it){
585 int straw = m_trtid->straw((*it)->identify());
586 const Amg::Vector3D& sc = (*it)->detectorElement()->strawCenter(straw);
587
588 double pred_phi=sc.z()*p_best[0]+p_best[1];
589 double corr=sc.phi()-p_best[2];
590
591 if(p_best[0]==1e10){
592 if(sc.z()==p_best[1]) {
593 pred_phi=corr; //z=const ...
594 }else{
595 pred_phi=corr+0.5;
596 }
597 }
598 if (corr > M_PI) corr = std::fmod(corr+M_PI,2.*M_PI)-M_PI;
599 else if(corr <-M_PI) corr = std::fmod(corr-M_PI,2.*M_PI)+M_PI;
600
601
602 if(std::abs(phidiff(pred_phi,corr))<=4./650.){
603 if (msgLvl(MSG::VERBOSE)) msg()<<"extended hit at z= "<<sc.z()<<endmsg;
604 if(sc.z()<zmin) zmin=sc.z();
605 if(sc.z()>zmax) zmax=sc.z();
606 seed->push_back(*it);
607 }
608 }
609 }
610 }
611
612 if(muster[0][0] || muster[0][1] || muster[0][2]){
613 //extend towards lower z-values
614 zsl=zslice-2;
615
616 for(int ps=sector-2;ps<=sector+2;ps++){
617
618 //make sure that the slices exist ...
619 if(zsl<0 || zsl>=20) continue;
620
621 int ps2=ps;
622 if(ps2<0) ps2+=16;
623 else if(ps2>15) ps2-=16;
624
625 for(it=event_data.m_sectors[endcap][zsl][ps2].begin();it!=event_data.m_sectors[endcap][zsl][ps2].end();++it){
626 int straw = m_trtid->straw((*it)->identify());
627 const Amg::Vector3D& sc = (*it)->detectorElement()->strawCenter(straw);
628
629 double pred_phi=sc.z()*p_best[0]+p_best[1];
630 double corr=sc.phi()-p_best[2];
631
632 if(p_best[0]==1e10){
633 if(sc.z()==p_best[1]) {
634 pred_phi=corr; //z=const ...
635 }else{
636 pred_phi=corr+0.5;
637 }
638 }
639 if (corr > M_PI) corr = std::fmod(corr+M_PI,2.*M_PI)-M_PI;
640 else if(corr <-M_PI) corr = std::fmod(corr-M_PI,2.*M_PI)+M_PI;
641
642
643 if(std::abs(phidiff(pred_phi,corr))<=4./650.){
644 ATH_MSG_VERBOSE("extended hit at z= "<<sc.z());
645 if(sc.z()<zmin) zmin=sc.z();
646 if(sc.z()>zmax) zmax=sc.z();
647 seed->push_back(*it);
648 }
649 }
650 }
651 }
652
653
654 if(seed->size()<(size_t)m_minDCSeed){
655 ATH_MSG_VERBOSE("Seed has too few hits: "<<seed->size());
656 seed->clear();
657 delete seed;
658 return false;
659 }
660
661
662 vit=seed->begin();
663 vitE=seed->end();
664
665 for(;vit!=vitE;++vit){
666 int straw = m_trtid->straw((*vit)->identify());
667 const Amg::Vector3D& sc = (*vit)->detectorElement()->strawCenter(straw);
668
669 int phi = int((std::atan2(sc.y(),sc.x())+M_PI)*8./M_PI);
670
671 if ( phi<0 ) phi = 15;
672 else if( phi>15 ) phi = 0 ;
673
674 int z=0;
675 if(sc.z()<0) z=1;
676 int zslic=int((std::abs(sc.z())-800.)/100.);
677
678 event_data.m_sectors[z][zslic][phi].remove(*vit);
679
680 }
681
682 //just to be sure: check also the noise hits in these 9 slices
683
684 std::list<const InDet::TRT_DriftCircle*>::iterator lit,litE;
685 lit=event_data.m_noiseHits.begin();
686 litE=event_data.m_noiseHits.end();
687
688 for(;lit!=litE;++lit){
689
690 //don't add hits that really look like noise
691 if((*lit)->timeOverThreshold()<m_cutToTLoose) continue;
692
693 int straw = m_trtid->straw((*lit)->identify());
694 const Amg::Vector3D& sc = (*lit)->detectorElement()->strawCenter(straw);
695
696 int phi = int((std::atan2(sc.y(),sc.x())+M_PI)*8./M_PI);
697
698 if ( phi<0 ) phi = 15;
699 else if( phi>15 ) phi = 0 ;
700
701
702 int zsl=int((std::abs(sc.z())-800.)/100.);
703
704 if(std::abs(zsl-zslice)<2){
705 double pred_phi=sc.z()*p_best[0]+p_best[1];
706 double corr=sc.phi()-p_best[2];
707
708 if(p_best[0]==1e10){
709 if(sc.z()==p_best[1]) {
710 pred_phi=corr; //z=const ...
711 }else{
712 pred_phi=corr+0.5;
713 }
714 }
715 if (corr > M_PI) corr = std::fmod(corr+M_PI,2.*M_PI)-M_PI;
716 else if(corr <-M_PI) corr = std::fmod(corr-M_PI,2.*M_PI)+M_PI;
717
718
719 if(std::abs(phidiff(pred_phi,corr))<=4./650.){
720 if(sc.z()>zmin && sc.z()<zmax){
721 seed->push_back(*lit);
722 ATH_MSG_VERBOSE("\t\t noise hit matched seed! z="<< sc.z());
723 }
724 }
725
726 }
727
728
729 }
730
731
732
733
734 event_data.m_seeds.push_back(seed);
735
736 return true;
737 }
738
739
740 return false;
741
742}
int evaluate_seed(int endcap, int zslice, int sector, const double *p, TRT_TrackSegmentsMaker_ECcosmics::EventData &event_data) const
Evaluate how many dc match this seed.
bool match(std::string s1, std::string s2)
match the individual directories of two strings
Definition hcg.cxx:357

◆ initialize()

StatusCode InDet::TRT_TrackSegmentsMaker_ECcosmics::initialize ( )
overridevirtual

Implements InDet::ITRT_TrackSegmentsMaker.

Definition at line 76 of file TRT_TrackSegmentsMaker_ECcosmics.cxx.

77{
78 StatusCode sc = AthAlgTool::initialize();
79
80 // Initialize ReadHandle
81 ATH_CHECK(m_trtname.initialize());
82
83 // TRT
84 if (detStore()->retrieve(m_trtid, "TRT_ID").isFailure()) {
85 ATH_MSG_FATAL("Could not get TRT ID helper");
86 return StatusCode::FAILURE;
87 }
88
89 // Get RIO_OnTrack creator with drift time information
90 //
91 if( m_riomakerD.retrieve().isFailure()) {
92 ATH_MSG_FATAL("Failed to retrieve tool "<< m_riomakerD);
93 return StatusCode::FAILURE;
94 }
95 ATH_MSG_DEBUG("Retrieved tool " << m_riomakerD);
96
97 // Get RIO_OnTrack creator without drift time information
98 //
99
100 if( m_riomakerN.retrieve().isFailure()) {
101 ATH_MSG_FATAL("Failed to retrieve tool "<< m_riomakerN);
102 return StatusCode::FAILURE;
103 }
104 ATH_MSG_DEBUG("Retrieved tool " << m_riomakerN);
105
106 // optional PRD to track association map
107 //
108 ATH_CHECK( m_prdToTrackMap.initialize( !m_prdToTrackMap.key().empty() ) );
109
110
111 return sc;
112}
#define ATH_CHECK
Evaluate an expression and check for errors.
#define ATH_MSG_FATAL(x)
const ServiceHandle< StoreGateSvc > & detStore() const
SG::ReadHandleKey< InDet::TRT_DriftCircleContainer > m_trtname
TRTs container.
SG::ReadHandleKey< Trk::PRDtoTrackMap > m_prdToTrackMap
retrieve(aClass, aKey=None)
Definition PyKernel.py:110

◆ inputHandles()

virtual std::vector< Gaudi::DataHandle * > AthCommonDataStore< AthCommonMsg< AlgTool > >::inputHandles ( ) const
overridevirtualinherited

Return this algorithm's input handles.

We override this to include handle instances from key arrays if they have not yet been declared. See comments on updateVHKA.

◆ interfaceID()

const InterfaceID & InDet::ITRT_TrackSegmentsMaker::interfaceID ( )
inlinestaticinherited

Definition at line 110 of file ITRT_TrackSegmentsMaker.h.

111 {
113 }
static const InterfaceID IID_ITRT_TrackSegmentsMaker("InDet::ITRT_TrackSegmentsMaker", 1, 0)

◆ is_suspicious()

bool InDet::TRT_TrackSegmentsMaker_ECcosmics::is_suspicious ( const InDet::TRT_DriftCircle * dc,
std::vector< const InDet::TRT_DriftCircle * > * seed ) const
protected

checks if a hit that matches the segment looks suspicious (i.e.

isolated)

Definition at line 1820 of file TRT_TrackSegmentsMaker_ECcosmics.cxx.

1822{
1823
1824 double mean=0,var=0;
1825 double meanz=0,varz=0;
1826
1827 std::vector<const InDet::TRT_DriftCircle*>::iterator vit,vitE;
1828
1829 vit=seed->begin();
1830 vitE=seed->end();
1831
1832 double dcz[500];
1833 int count=0;
1834 int index=-1;
1835 for(;vit!=vitE;++vit){
1836 int straw = m_trtid->straw((*vit)->identify());
1837 const Amg::Vector3D& sc = (*vit)->detectorElement()->strawCenter(straw);
1838
1839 if(dc==*vit) index=count;
1840 dcz[count++]=sc.z();
1841
1842 meanz+=sc.z();
1843 varz+=sc.z()*sc.z();
1844 }
1845
1846 if(count < 2) return true;
1847
1848 varz=(count*varz-meanz*meanz);
1849 varz/=(double)(count*(count-1));
1850
1851 meanz/=(double)count;
1852
1853 varz=std::sqrt(varz);
1854
1855 double zmin_sel=100000;
1856
1857 for(int i=0;i<count;i++){
1858 double zmin=10000;
1859
1860 if(i-1>=0)
1861 if(std::abs(dcz[i]-dcz[i-1])<zmin) zmin=std::abs(dcz[i]-dcz[i-1]);
1862 if(i+1<count)
1863 if(std::abs(dcz[i]-dcz[i+1])<zmin) zmin=std::abs(dcz[i]-dcz[i+1]);
1864
1865 if(i==index)
1866 zmin_sel=zmin;
1867
1868 mean+=zmin;
1869 var+=zmin*zmin;
1870 }
1871
1872 var=(count*var-mean*mean);
1873 var/=(double)(count*(count-1));
1874
1875 mean/=(double)count;
1876
1877 var=std::sqrt(var);
1878
1879 if(var<6) var=6;
1880
1881 if(std::abs(zmin_sel-mean)>2*var){// || std::abs(dcz[index]-meanz)>2*varz){
1882 if(index>=0) {
1883 ATH_MSG_VERBOSE("Hit is suspicious, remove it! "<<std::abs(zmin_sel-mean)<<" , "<<2*var<<" -> "<<zmin_sel<< " : "<<dcz[index]<<" <-> "<<meanz<<" ("<<varz<<")");
1884 }
1885 return true;
1886 }else{
1887 if(index>=0) {
1888 ATH_MSG_VERBOSE("Hit seems to be ok: "<<std::abs(zmin_sel-mean)<<" < "<<2*var<<" && "<<std::abs(dcz[index]-meanz)<<" < "<<1.5*varz);
1889 }
1890 }
1891 return false;
1892}
void mean(std::vector< double > &bins, std::vector< double > &values, const std::vector< std::string > &files, const std::string &histname, const std::string &tplotname, const std::string &label="")
str index
Definition DeMoScan.py:362

◆ msg()

MsgStream & AthCommonMsg< AlgTool >::msg ( ) const
inlineinherited

Definition at line 24 of file AthCommonMsg.h.

24 {
25 return this->msgStream();
26 }

◆ msgLvl()

bool AthCommonMsg< AlgTool >::msgLvl ( const MSG::Level lvl) const
inlineinherited

Definition at line 30 of file AthCommonMsg.h.

30 {
31 return this->msgLevel(lvl);
32 }

◆ newEvent()

std::unique_ptr< InDet::ITRT_TrackSegmentsMaker::IEventData > InDet::TRT_TrackSegmentsMaker_ECcosmics::newEvent ( const EventContext & ctx) const
overridevirtual

Implements InDet::ITRT_TrackSegmentsMaker.

Definition at line 156 of file TRT_TrackSegmentsMaker_ECcosmics.cxx.

157{
158
159 std::unique_ptr<TRT_TrackSegmentsMaker_ECcosmics::EventData>
160 event_data = std::make_unique<TRT_TrackSegmentsMaker_ECcosmics::EventData>();
161
162 setFitFunctions(*event_data);
163 //sort into good/noise hits
164 retrieveHits(*event_data);
165
166 std::list<const InDet::TRT_DriftCircle*>::iterator it,itE;
167 it=event_data->m_goodHits.begin();
168 itE=event_data->m_goodHits.end();
169
170 //loop over good hits and sort them into sectors
171
172 for(;it!=itE;++it){
173
174 int straw = m_trtid->straw((*it)->identify());
175 const Amg::Vector3D& sc = (*it)->detectorElement()->strawCenter(straw);
176
177 int phi = int((std::atan2(sc.y(),sc.x())+M_PI)*8./M_PI);
178
179 if ( phi<0 ) phi = 15;
180 else if( phi>15 ) phi = 0 ;
181
182 int z=0;
183 if(sc.z()<0) z=1;
184
185 int zslice=int((std::abs(sc.z())-800.)/100.);
186
187 event_data->m_sectors[z][zslice][phi].push_back(*it);
188
189 }
190
191 ATH_MSG_DEBUG("Number of good driftcircles: "<<event_data->m_goodHits.size());
192 ATH_MSG_DEBUG("Number of noise driftcircles: "<<event_data->m_noiseHits.size());
193
194 return std::unique_ptr<InDet::ITRT_TrackSegmentsMaker::IEventData>(event_data.release());
195}
void setFitFunctions(TRT_TrackSegmentsMaker_ECcosmics::EventData &event_data) const
void retrieveHits(TRT_TrackSegmentsMaker_ECcosmics::EventData &event_data) const
sort hits into good/noise lists

◆ newRegion()

std::unique_ptr< InDet::ITRT_TrackSegmentsMaker::IEventData > InDet::TRT_TrackSegmentsMaker_ECcosmics::newRegion ( const EventContext & ctx,
const std::vector< IdentifierHash > & vTRT ) const
overridevirtual

Implements InDet::ITRT_TrackSegmentsMaker.

Definition at line 202 of file TRT_TrackSegmentsMaker_ECcosmics.cxx.

204{
205 (void) vTRT;
206 std::unique_ptr<TRT_TrackSegmentsMaker_ECcosmics::EventData>
207 event_data = std::make_unique<TRT_TrackSegmentsMaker_ECcosmics::EventData>();
208
209 setFitFunctions(*event_data);
210
211 // @TODO here the retrieve_hits method should be called with the list of the idhash
212
213
214 std::list<const InDet::TRT_DriftCircle*>::iterator it,itE;
215 it=event_data->m_goodHits.begin();
216 itE=event_data->m_goodHits.end();
217
218 //loop over good hits and sort them into sectors
219
220 for(;it!=itE;++it){
221
222 int straw = m_trtid->straw((*it)->identify());
223 const Amg::Vector3D& sc = (*it)->detectorElement()->strawCenter(straw);
224
225 int phi = int((std::atan2(sc.y(),sc.x())+M_PI)*8./M_PI);
226
227 if ( phi<0 ) phi = 15;
228 else if( phi>15 ) phi = 0 ;
229
230 int z=0;
231 if(sc.z()<0) z=1;
232 int zslice=int((std::abs(sc.z())-800.)/100.);
233
234 event_data->m_sectors[z][zslice][phi].push_back(*it);
235
236 }
237
238 ATH_MSG_DEBUG("Number of good driftcircles: "<<event_data->m_goodHits.size());
239 ATH_MSG_DEBUG("Number of noise driftcircles: "<<event_data->m_noiseHits.size());
240
241 return std::unique_ptr<InDet::ITRT_TrackSegmentsMaker::IEventData>(event_data.release());
242}

◆ next()

Trk::TrackSegment * InDet::TRT_TrackSegmentsMaker_ECcosmics::next ( InDet::ITRT_TrackSegmentsMaker::IEventData & event_data) const
overridevirtual

Implements InDet::ITRT_TrackSegmentsMaker.

Definition at line 1515 of file TRT_TrackSegmentsMaker_ECcosmics.cxx.

1516{
1517 TRT_TrackSegmentsMaker_ECcosmics::EventData &
1519
1520 if(event_data.m_segiterator!=event_data.m_segments.end()) return (*event_data.m_segiterator++);
1521 return nullptr;
1522}

◆ outputHandles()

virtual std::vector< Gaudi::DataHandle * > AthCommonDataStore< AthCommonMsg< AlgTool > >::outputHandles ( ) const
overridevirtualinherited

Return this algorithm's output handles.

We override this to include handle instances from key arrays if they have not yet been declared. See comments on updateVHKA.

◆ perform_fit()

TF1 * InDet::TRT_TrackSegmentsMaker_ECcosmics::perform_fit ( int count,
TRT_TrackSegmentsMaker_ECcosmics::EventData & event_data ) const
protected

Perform the fit and return a function that provides the fitted phi information.

Definition at line 1734 of file TRT_TrackSegmentsMaker_ECcosmics.cxx.

1736{
1737 // @TODO not thread-safe enough! need a globak root lock !
1738 std::lock_guard<std::mutex> lock(s_fitMutex);
1739 TVirtualFitter::SetMaxIterations(100000);
1740
1741 RootUtils::WithRootErrorHandler hand ([] (int, bool, const char*, const char*) { return false; });
1742
1743
1744 event_data.m_fitf_ztanphi->SetParameter(0,0.);
1745 event_data.m_fitf_ztanphi->SetParameter(1,0.);
1746 event_data.m_fitf_ztanphi->SetParameter(2,0.);
1747
1748 TGraphErrors *gre=new TGraphErrors(count,event_data.m_a_z,event_data.m_a_tan,event_data.m_a_z_err,event_data.m_a_tan_err);
1749 int fitresult1=gre->Fit(event_data.m_fitf_ztanphi,"Q");
1750
1751 ATH_MSG_VERBOSE("Fit result (ztanphi): a0="<<event_data.m_fitf_ztanphi->GetParameter(0)
1752 <<" a1="<<event_data.m_fitf_ztanphi->GetParameter(1)<<" a2="
1753 <<event_data.m_fitf_ztanphi->GetParameter(2)<<" fitresult="<<fitresult1);
1754
1755 //copy over the parameters
1756 event_data.m_fitf_zphi->SetParameter(0,event_data.m_fitf_ztanphi->GetParameter(0));
1757 event_data.m_fitf_zphi->SetParameter(1,event_data.m_fitf_ztanphi->GetParameter(1));
1758 event_data.m_fitf_zphi->SetParameter(2,event_data.m_fitf_ztanphi->GetParameter(2));
1759
1760 delete gre;
1761
1762 gre=new TGraphErrors(count,event_data.m_a_z,event_data.m_a_phi,event_data.m_a_z_err,event_data.m_a_phi_err);
1763 int fitresult2=gre->Fit(event_data.m_fitf_zphi_approx,"Q");
1764
1765
1766 ATH_MSG_VERBOSE( "Fit result (zphi): a0="<< event_data.m_fitf_zphi_approx->GetParameter(0)
1767 <<" a1="<< event_data.m_fitf_zphi_approx->GetParameter(1)<<" a2="
1768 << event_data.m_fitf_zphi_approx->GetParameter(2)<<" fitresult="<<fitresult2);
1769
1770 double p1=event_data.m_fitf_ztanphi->GetProb();
1771 double p2=event_data.m_fitf_zphi_approx->GetProb();
1772
1773 ATH_MSG_VERBOSE("tanphi prob: "<<p1<<" pol2 prob : "<<p2);
1774
1775 delete gre;
1776
1777
1778 if(fitresult1!=0 && fitresult2!=0) return nullptr;
1779 if(fitresult1!=0 && fitresult2==0) return event_data.m_fitf_zphi_approx;
1780 if(fitresult1==0 && fitresult2!=0) return event_data.m_fitf_zphi;
1781
1782
1783 int match_tan=0;
1784 int match_phi=0;
1785 for(int i=0;i<count;i++){
1786 double phi1=event_data.m_fitf_zphi->Eval(event_data.m_a_z[i]);
1787 double phi2=event_data.m_fitf_zphi_approx->Eval(event_data.m_a_z[i]);
1788
1789 if(std::abs(phi1-event_data.m_a_phi[i])<2*event_data.m_a_phi_err[i])
1790 match_tan++;
1791
1792 if(std::abs(phi2-event_data.m_a_phi[i])<2*event_data.m_a_phi_err[i])
1793 match_phi++;
1794
1795 }
1796
1797
1798 ATH_MSG_VERBOSE("Number of hits matching this fit: tan="<<match_tan<<" pol2="<<match_phi);
1799
1800 if(match_tan>match_phi){
1801 return event_data.m_fitf_zphi;
1802 }else if(match_tan<match_phi){
1803 return event_data.m_fitf_zphi_approx;
1804 }
1805
1806 if(p1>p2){
1807 return event_data.m_fitf_zphi;
1808 }else{
1809 return event_data.m_fitf_zphi_approx;
1810 }
1811
1812}

◆ phidiff()

double InDet::TRT_TrackSegmentsMaker_ECcosmics::phidiff ( double a,
double b )
inlinestaticprotected

provide the proper subtraction of two phi values

Definition at line 250 of file TRT_TrackSegmentsMaker_ECcosmics.h.

251{
252 double c=a-b;
253
254 if(std::abs(c)>100)
255 c/=int(c/(2*M_PI));
256 while(c>=M_PI) c-=2*M_PI;
257 while(c<-M_PI) c+=2*M_PI;
258 return c;
259}
static Double_t a

◆ renounce()

std::enable_if_t< std::is_void_v< std::result_of_t< decltype(&T::renounce)(T)> > &&!std::is_base_of_v< SG::VarHandleKeyArray, T > &&std::is_base_of_v< Gaudi::DataHandle, T >, void > AthCommonDataStore< AthCommonMsg< AlgTool > >::renounce ( T & h)
inlineprotectedinherited

Definition at line 380 of file AthCommonDataStore.h.

381 {
382 h.renounce();
384 }
std::enable_if_t< std::is_void_v< std::result_of_t< decltype(&T::renounce)(T)> > &&!std::is_base_of_v< SG::VarHandleKeyArray, T > &&std::is_base_of_v< Gaudi::DataHandle, T >, void > renounce(T &h)

◆ renounceArray()

void AthCommonDataStore< AthCommonMsg< AlgTool > >::renounceArray ( SG::VarHandleKeyArray & handlesArray)
inlineprotectedinherited

remove all handles from I/O resolution

Definition at line 364 of file AthCommonDataStore.h.

364 {
366 }

◆ retrieveHits()

void InDet::TRT_TrackSegmentsMaker_ECcosmics::retrieveHits ( TRT_TrackSegmentsMaker_ECcosmics::EventData & event_data) const
protected

sort hits into good/noise lists

Definition at line 1548 of file TRT_TrackSegmentsMaker_ECcosmics.cxx.

1549{
1550
1551 event_data.m_noiseHits.clear();
1552 event_data.m_goodHits.clear();
1553
1554
1555 SG::ReadHandle<InDet::TRT_DriftCircleContainer> trtcontainer(m_trtname);
1556 if (not trtcontainer.isValid()) {
1557 ATH_MSG_DEBUG("Could not get TRT_DriftCircleContainer");
1558 return;
1559 }
1560
1561 SG::ReadHandle<Trk::PRDtoTrackMap> prd_to_track_map;
1562 const Trk::PRDtoTrackMap *prd_to_track_map_cptr=nullptr;
1563 if (!m_prdToTrackMap.key().empty()) {
1564 prd_to_track_map = SG::ReadHandle<Trk::PRDtoTrackMap>(m_prdToTrackMap);
1565 if (!prd_to_track_map.isValid()) {
1566 ATH_MSG_ERROR("Failed to read PRD to track association map: " << m_prdToTrackMap );
1567 }
1568 prd_to_track_map_cptr = prd_to_track_map.cptr();
1569 }
1570
1571 const InDet::TRT_DriftCircleContainer* mjo_trtcontainer = trtcontainer.get();
1572 if (!mjo_trtcontainer) return;
1573
1574 InDet::TRT_DriftCircleContainer::const_iterator
1575 w = trtcontainer->begin(),we = trtcontainer->end();
1576 for(; w!=we; ++w) {
1577
1578 Identifier ID = (*w)->identify();
1579 if(std::abs(m_trtid->barrel_ec(ID))!=2) continue; // only endcaps please!!!
1580
1581 InDet::TRT_DriftCircleCollection::const_iterator
1582 c = (*w)->begin(), ce = (*w)->end();
1583
1584 for(; c!=ce; ++c) {
1585
1586 if(prd_to_track_map_cptr && prd_to_track_map_cptr->isUsed(*(*c))) continue; //don't use hits that have already been used
1587
1588 if( (*c)->isNoise() || (*c)->timeOverThreshold()<m_cutToTLoose) {
1589 event_data.m_noiseHits.push_back(*c);
1590 }else{
1591 event_data.m_goodHits.push_back(*c);
1592 }
1593
1594 }
1595 }
1596
1597 ATH_MSG_DEBUG("good hits in endcap: "<<event_data.m_goodHits.size());
1598 ATH_MSG_DEBUG("noise hits in endcap: "<<event_data.m_noiseHits.size());
1599
1600
1601 if(event_data.m_goodHits.size()>(size_t)m_hitLimit){
1602 // this event is too busy ...
1603 ATH_MSG_DEBUG("This event has more than "<<m_hitLimit<<" good hits. Aborting segment finding ...");
1604 event_data.m_goodHits.clear();
1605 event_data.m_noiseHits.clear();
1606 return;
1607 }
1608
1609 //move isolated hits from the good hits vector to the noise hits vector
1610 //
1611 auto it=event_data.m_goodHits.begin();
1612 auto itE=event_data.m_goodHits.end();
1613 auto it_remove(it);
1614 bool remove=false;
1615 //
1616 const double dPhiLimit{0.03};
1617 const double dzLimit{31.};
1618 for(;it!=itE;++it){
1619 if(remove){
1620 event_data.m_goodHits.erase(it_remove);
1621 remove=false;
1622 }
1623 const bool accept=accepted(it, event_data.m_goodHits, dPhiLimit, dzLimit);
1624 if(!accept){
1625 event_data.m_noiseHits.push_back(*it);
1626 it_remove=it;
1627 remove=true;
1628 }
1629 }
1630 if(remove){
1631 event_data.m_goodHits.erase(it_remove);
1632 remove=false;
1633 }
1634
1635 //2nd pass through good hits with increased ToT cut
1636
1637 it=event_data.m_goodHits.begin();
1638 itE=event_data.m_goodHits.end();
1639
1640 remove=false;
1641
1642 for(;it!=itE;++it){
1643
1644 if(remove){
1645 event_data.m_goodHits.erase(it_remove);
1646 remove=false;
1647 }
1648
1649 if((*it)->timeOverThreshold()<m_cutToTTight || (*it)->timeOverThreshold()>m_cutToTUpper){
1650 event_data.m_noiseHits.push_back(*it);
1651 it_remove=it;
1652 remove=true;
1653 }
1654 }
1655 if(remove){
1656 event_data.m_goodHits.erase(it_remove);
1657 remove=false;
1658 }
1659
1660
1661 //2nd pass through good hits with loosened isolation criteria
1662
1663 it=event_data.m_goodHits.begin();
1664 itE=event_data.m_goodHits.end();
1665
1666 remove=false;
1667
1668 for(;it!=itE;++it){
1669 if(remove){
1670 event_data.m_goodHits.erase(it_remove);
1671 remove=false;
1672 }
1673 const bool accept=accepted(it, event_data.m_goodHits, 0.25, 200.);
1674 if(!accept){
1675 ATH_MSG_DEBUG("Removing hit in 2nd pass isolation");
1676 event_data.m_noiseHits.push_back(*it);
1677 it_remove=it;
1678 remove=true;
1679 }
1680 }
1681
1682 if(remove){
1683 event_data.m_goodHits.erase(it_remove);
1684 remove=false;
1685 }
1686
1687 ATH_MSG_DEBUG("good hits in endcap: "<<event_data.m_goodHits.size());
1688 ATH_MSG_DEBUG("noise hits in endcap: "<<event_data.m_noiseHits.size());
1689
1690}
#define ATH_MSG_ERROR(x)
std::vector< Identifier > ID
bool accepted(const std::list< const InDet::TRT_DriftCircle * >::iterator compareIt, std::list< const InDet::TRT_DriftCircle * > &container, double phiLimit, double dzLimit) const
is the hit accepted?
Trk::PrepRawDataContainer< TRT_DriftCircleCollection > TRT_DriftCircleContainer
virtual bool isValid() override final
Can the handle be successfully dereferenced?
const_pointer_type cptr()
Dereference the pointer.
bool isUsed(const PrepRawData &prd) const
does this PRD belong to at least one track?

◆ setFitFunctions()

void InDet::TRT_TrackSegmentsMaker_ECcosmics::setFitFunctions ( TRT_TrackSegmentsMaker_ECcosmics::EventData & event_data) const
protected

Definition at line 115 of file TRT_TrackSegmentsMaker_ECcosmics.cxx.

115 {
116 std::string fit_name1="ztanphi";
117 std::string fit_name2="zphi";
118 std::string fit_name3="zphi_ap";
119
120 if(m_phaseMode){
121 fit_name1+="_phase";
122 fit_name2+="_phase";
123 fit_name3+="_phase";
124 }
125
126 if(!m_prdToTrackMap.key().empty()){
127 fit_name1+="_asso";
128 fit_name2+="_asso";
129 fit_name3+="_asso";
130 }
131
132 fit_name1+=name();
133 fit_name2+=name();
134 fit_name3+=name();
135
136 event_data.m_fitf_ztanphi = new TF1(fit_name1.c_str(),fitf_ztanphi,-3000,3000,3);
137 event_data.m_fitf_zphi = new TF1(fit_name2.c_str(),fitf_zphi,-3000,3000,3);
138 event_data.m_fitf_zphi_approx = new TF1(fit_name3.c_str(),fitf_zphi_approx,-3000,3000,3);
139}
Double_t fitf_ztanphi(const Double_t *x, const Double_t *par)
Double_t fitf_zphi_approx(const Double_t *x, const Double_t *par)
Double_t fitf_zphi(const Double_t *x, const Double_t *par)

◆ sysInitialize()

virtual StatusCode AthCommonDataStore< AthCommonMsg< AlgTool > >::sysInitialize ( )
overridevirtualinherited

Perform system initialization for an algorithm.

We override this to declare all the elements of handle key arrays at the end of initialization. See comments on updateVHKA.

Reimplemented in asg::AsgMetadataTool, AthCheckedComponent< AthAlgTool >, AthCheckedComponent<::AthAlgTool >, and DerivationFramework::CfAthAlgTool.

◆ sysStart()

virtual StatusCode AthCommonDataStore< AthCommonMsg< AlgTool > >::sysStart ( )
overridevirtualinherited

Handle START transition.

We override this in order to make sure that conditions handle keys can cache a pointer to the conditions container.

◆ updateVHKA()

void AthCommonDataStore< AthCommonMsg< AlgTool > >::updateVHKA ( Gaudi::Details::PropertyBase & )
inlineinherited

Definition at line 308 of file AthCommonDataStore.h.

308 {
309 // debug() << "updateVHKA for property " << p.name() << " " << p.toString()
310 // << " size: " << m_vhka.size() << endmsg;
311 for (auto &a : m_vhka) {
313 for (auto k : keys) {
314 k->setOwner(this);
315 }
316 }
317 }
std::vector< SG::VarHandleKeyArray * > m_vhka

Member Data Documentation

◆ m_cutToTLoose

DoubleProperty InDet::TRT_TrackSegmentsMaker_ECcosmics::m_cutToTLoose
protected
Initial value:
{this, "ToTCutLoose", 7.,
"Loose cut on ToT (preselection)"}

Definition at line 191 of file TRT_TrackSegmentsMaker_ECcosmics.h.

191 {this, "ToTCutLoose", 7.,
192 "Loose cut on ToT (preselection)"};

◆ m_cutToTTight

DoubleProperty InDet::TRT_TrackSegmentsMaker_ECcosmics::m_cutToTTight
protected
Initial value:
{this, "ToTCutTight", 18.,
"Hard cut on ToT (preselection)"}

Definition at line 193 of file TRT_TrackSegmentsMaker_ECcosmics.h.

193 {this, "ToTCutTight", 18.,
194 "Hard cut on ToT (preselection)"};

◆ m_cutToTUpper

DoubleProperty InDet::TRT_TrackSegmentsMaker_ECcosmics::m_cutToTUpper
protected
Initial value:
{this, "ToTCutUpper", 32.,
"Upper cut on ToT (preselection)"}

Definition at line 195 of file TRT_TrackSegmentsMaker_ECcosmics.h.

195 {this, "ToTCutUpper", 32.,
196 "Upper cut on ToT (preselection)"};

◆ m_detStore

StoreGateSvc_t AthCommonDataStore< AthCommonMsg< AlgTool > >::m_detStore
privateinherited

Pointer to StoreGate (detector store by default)

Definition at line 393 of file AthCommonDataStore.h.

◆ m_evtStore

StoreGateSvc_t AthCommonDataStore< AthCommonMsg< AlgTool > >::m_evtStore
privateinherited

Pointer to StoreGate (event store by default)

Definition at line 390 of file AthCommonDataStore.h.

◆ m_hitLimit

IntegerProperty InDet::TRT_TrackSegmentsMaker_ECcosmics::m_hitLimit
protected
Initial value:
{this, "HitLimit", 2000,
"Maximum number of good hits (i.e. after noise cut) in endcap"}

Definition at line 199 of file TRT_TrackSegmentsMaker_ECcosmics.h.

199 {this, "HitLimit", 2000,
200 "Maximum number of good hits (i.e. after noise cut) in endcap"};

◆ m_minDCSeed

IntegerProperty InDet::TRT_TrackSegmentsMaker_ECcosmics::m_minDCSeed
protected
Initial value:
{this, "MinDCperSeed", 7,
"Minimum number of driftcircles to form a seed"}

Definition at line 197 of file TRT_TrackSegmentsMaker_ECcosmics.h.

197 {this, "MinDCperSeed", 7,
198 "Minimum number of driftcircles to form a seed"};

◆ m_phaseMode

BooleanProperty InDet::TRT_TrackSegmentsMaker_ECcosmics::m_phaseMode
protected
Initial value:
{this, "Phase", false,
"Switch to destinguish between phase calculation and full reco"}

Definition at line 166 of file TRT_TrackSegmentsMaker_ECcosmics.h.

166 {this, "Phase", false,
167 "Switch to destinguish between phase calculation and full reco"};

◆ m_prdToTrackMap

SG::ReadHandleKey<Trk::PRDtoTrackMap> InDet::TRT_TrackSegmentsMaker_ECcosmics::m_prdToTrackMap {this,"PRDtoTrackMap",""}
protected

Definition at line 171 of file TRT_TrackSegmentsMaker_ECcosmics.h.

172{this,"PRDtoTrackMap",""};

◆ m_riomakerD

ToolHandle<Trk::IRIO_OnTrackCreator> InDet::TRT_TrackSegmentsMaker_ECcosmics::m_riomakerD
protected
Initial value:
{this, "RIOonTrackToolYesDr",
"InDet::TRT_DriftCircleOnTrackTool/TRT_DriftCircleOnTrackTool",
"RI0_onTrack creator with drift information"}

Definition at line 174 of file TRT_TrackSegmentsMaker_ECcosmics.h.

175 {this, "RIOonTrackToolYesDr",
176 "InDet::TRT_DriftCircleOnTrackTool/TRT_DriftCircleOnTrackTool",
177 "RI0_onTrack creator with drift information"};

◆ m_riomakerN

ToolHandle<Trk::IRIO_OnTrackCreator> InDet::TRT_TrackSegmentsMaker_ECcosmics::m_riomakerN
protected
Initial value:
{this, "RIOonTrackToolNoDr",
"InDet::TRT_DriftCircleOnTrackNoDriftTimeTool/TRT_DriftCircleOnTrackNoDriftTimeTool",
"RI0_onTrack creator without drift information"}

Definition at line 178 of file TRT_TrackSegmentsMaker_ECcosmics.h.

179 {this, "RIOonTrackToolNoDr",
180 "InDet::TRT_DriftCircleOnTrackNoDriftTimeTool/TRT_DriftCircleOnTrackNoDriftTimeTool",
181 "RI0_onTrack creator without drift information"};

◆ m_scaleFactorDrift

DoubleProperty InDet::TRT_TrackSegmentsMaker_ECcosmics::m_scaleFactorDrift
protected
Initial value:
{this, "ScaleFactorDrift", 1.,
"Scalefactor for uncertainty of drifttime hits"}

Definition at line 187 of file TRT_TrackSegmentsMaker_ECcosmics.h.

187 {this, "ScaleFactorDrift", 1.,
188 "Scalefactor for uncertainty of drifttime hits"};

◆ m_scaleTube

DoubleProperty InDet::TRT_TrackSegmentsMaker_ECcosmics::m_scaleTube
protected
Initial value:
{this, "ScaleFactorTube", 2.,
"Scalefactor for uncertainty of tube hits"}

Definition at line 185 of file TRT_TrackSegmentsMaker_ECcosmics.h.

185 {this, "ScaleFactorTube", 2.,
186 "Scalefactor for uncertainty of tube hits"};

◆ m_scaleTubeNoise

DoubleProperty InDet::TRT_TrackSegmentsMaker_ECcosmics::m_scaleTubeNoise
protected
Initial value:
{this, "ScaleFactorTubeNoise", 1.,
"Scalefactor for uncertainty of tube hits flagged as noise"}

Definition at line 189 of file TRT_TrackSegmentsMaker_ECcosmics.h.

189 {this, "ScaleFactorTubeNoise", 1.,
190 "Scalefactor for uncertainty of tube hits flagged as noise"};

◆ m_trtid

const TRT_ID* InDet::TRT_TrackSegmentsMaker_ECcosmics::m_trtid {}
protected

Definition at line 168 of file TRT_TrackSegmentsMaker_ECcosmics.h.

168{} ;

◆ m_trtname

SG::ReadHandleKey<InDet::TRT_DriftCircleContainer> InDet::TRT_TrackSegmentsMaker_ECcosmics::m_trtname {this,"TRT_ClustersContainer","TRT_DriftCircles","RHK to retrieve TRT_DriftCircles"}
protected

TRTs container.

Definition at line 170 of file TRT_TrackSegmentsMaker_ECcosmics.h.

170{this,"TRT_ClustersContainer","TRT_DriftCircles","RHK to retrieve TRT_DriftCircles"};

◆ m_useDriftTime

BooleanProperty InDet::TRT_TrackSegmentsMaker_ECcosmics::m_useDriftTime
protected
Initial value:
{this, "UseDriftTime", false,
"Shall the drifttime be used or only tube hits?"}

Definition at line 183 of file TRT_TrackSegmentsMaker_ECcosmics.h.

183 {this, "UseDriftTime", false,
184 "Shall the drifttime be used or only tube hits?"};

◆ m_varHandleArraysDeclared

bool AthCommonDataStore< AthCommonMsg< AlgTool > >::m_varHandleArraysDeclared
privateinherited

Definition at line 399 of file AthCommonDataStore.h.

◆ m_vhka

std::vector<SG::VarHandleKeyArray*> AthCommonDataStore< AthCommonMsg< AlgTool > >::m_vhka
privateinherited

Definition at line 398 of file AthCommonDataStore.h.

◆ s_fitMutex

std::mutex InDet::TRT_TrackSegmentsMaker_ECcosmics::s_fitMutex
staticprotected

Definition at line 202 of file TRT_TrackSegmentsMaker_ECcosmics.h.


The documentation for this class was generated from the following files: