00001 #include "GaudiKernel/MsgStream.h"
00002 #include "GaudiKernel/SmartDataPtr.h"
00003
00004 #include "EventModel/EventHeader.h"
00005 #include "EvtRecEvent/EvtRecEvent.h"
00006 #include "EvtRecEvent/EvtRecTrack.h"
00007 #include "EvtRecEvent/EvtRecPi0.h"
00008 #include "VertexFit/KinematicFit.h"
00009 #include "Pi0RecAlg/Pi0RecAlg.h"
00010 #include "Pi0RecAlg/MakeGroupList.h"
00011 const double xmpi0 = 0.13497;
00012
00013
00014 Pi0RecAlg::Pi0RecAlg(const std::string& name, ISvcLocator* pSvcLocator) :
00015 Algorithm(name, pSvcLocator)
00016 {
00017
00018 declareProperty("EMinOfBarrel", _pi0_cut.MinEnergyOfBarrelPhoton = 0.025);
00019 declareProperty("EMinOfEndcap", _pi0_cut.MinEnergyOfEndcapPhoton = 0.045);
00020 declareProperty("Angle", _pi0_cut.MinAngle = 20);
00021 declareProperty("TimeLimit", _pi0_cut.TimeLimit = 20);
00022
00023 declareProperty("Pi0MinMass", _pi0_cut.MinMass = 0.10);
00024 declareProperty("Pi0MaxMass", _pi0_cut.MaxMass = 0.18);
00025 declareProperty("ChisqCut", _pi0_cut.Chisq = 50.0);
00026 }
00027
00028
00029 StatusCode Pi0RecAlg::initialize() {
00030
00031 MsgStream log(msgSvc(), name());
00032 log << MSG::INFO << "in initialize()" <<endreq;
00033
00034 return StatusCode::SUCCESS;
00035 }
00036
00037
00038
00039 StatusCode Pi0RecAlg::execute() {
00040
00041 MsgStream log(msgSvc(), name());
00042 log << MSG::INFO << "in execute()" <<endreq;
00043
00044
00045 SmartDataPtr<Event::EventHeader> eventHeader(eventSvc(),"/Event/EventHeader");
00046 SmartDataPtr<EvtRecEvent> recEvt(eventSvc(), EventModel::EvtRec::EvtRecEvent);
00047 SmartDataPtr<EvtRecTrackCol> recTrkCol(eventSvc(), EventModel::EvtRec::EvtRecTrackCol);
00048
00049 bool save2TDS = false;
00050 SmartDataPtr<EvtRecPi0Col> recPi0Col(eventSvc(), EventModel::EvtRec::EvtRecPi0Col);
00051 if ( !recPi0Col ) {
00052 recPi0Col = new EvtRecPi0Col;
00053 save2TDS = true;
00054 }
00055 using namespace Pi0;
00056
00057 static Criteria cri_inv(_pi0_cut.MinMass, _pi0_cut.MaxMass);
00058
00059
00060
00061 UserPi0Cut::SetForTrack(recEvt, recTrkCol);
00062
00063
00064 Pi0::make_gamma_list(_pi0_cut);
00065 Pi0::make_pi0_list(Pi0::GetDefaultGammaList());
00066
00067 Pi0::apply_criteria(cri_inv);
00068 Pi0::GetCandidatePi0List().sort(high_momentum());
00069
00070 Pi0::Pi0ListToTDS(Pi0::GetCandidatePi0List(), recPi0Col);
00071
00072 if ( save2TDS ) {
00073 StatusCode sc = eventSvc()->registerObject(EventModel::EvtRec::EvtRecPi0Col, recPi0Col);
00074 if ( sc.isFailure() ) {
00075 log << MSG::ERROR << "could not register EvtRecPi0Col in TDS" <<endreq;
00076 return StatusCode::FAILURE;
00077 }
00078 }
00079 return StatusCode::SUCCESS;
00080 }
00081
00082
00083
00084 StatusCode Pi0RecAlg::finalize() {
00085
00086 MsgStream log(msgSvc(), name());
00087 log << MSG::INFO << "in finalize()" << endreq;
00088
00089 return StatusCode::SUCCESS;
00090 }
00091
00092 HepLorentzVector Pi0RecAlg::getP4(RecEmcShower* gTrk) {
00093
00094 double eraw = gTrk->energy();
00095 double phi = gTrk->phi();
00096 double the = gTrk->theta();
00097
00098 return HepLorentzVector( eraw * sin(the) * cos(phi),
00099 eraw * sin(the) * sin(phi),
00100 eraw * cos(the),
00101 eraw );
00102 }
00103
00104
00105
00106
00107
00108
00109
00110
00111
00112
00113
00114
00115
00116
00117
00118
00119
00120
00121
00122
00123
00124
00125
00126
00127
00128
00129
00130
00131
00132
00133
00134
00135
00136
00137
00138
00139
00140
00141
00142
00143
00144
00145
00146
00147
00148
00149
00150
00151
00152
00153
00154
00155
00156
00157
00158
00159
00160
00161
00162
00163
00164
00165
00166
00167
00168
00169
00170
00171
00172
00173
00174
00175
00176
00177
00178
00179
00180
00181
00182