00001 #include "GaudiKernel/MsgStream.h"
00002
00003
00004
00005 #include "GaudiKernel/Bootstrap.h"
00006 #include "GaudiKernel/SmartDataPtr.h"
00007
00008 #include "EventModel/Event.h"
00009 #include "EventModel/EventHeader.h"
00010 #include "MdcRecEvent/RecMdcTrack.h"
00011 #include "HltAlgorithms/EFChargedTrack.h"
00012
00013 using namespace Event;
00014
00015 EFChargedTrack::EFChargedTrack(const std::string& name, ISvcLocator* pSvcLocator) :
00016 IEFAlgorithm(name, pSvcLocator){
00017 int output = (m_output%100)/10;
00018
00019 MsgStream log(msgSvc(), name);
00020 msgSvc()->setOutputLevel(name,output);
00021 m_ntrk = new CriteriaItemValue;
00022 m_acol = new CriteriaItemValue;
00023 m_mbal = new CriteriaItemValue;
00024 m_pmax1 = new CriteriaItemValue;
00025 m_pmax2 = new CriteriaItemValue;
00026 m_cost1 = new CriteriaItemValue;
00027 m_cost2 = new CriteriaItemValue;
00028 m_vr = new CriteriaItemValue;
00029 m_vz = new CriteriaItemValue;
00030 }
00031
00032 EFChargedTrack::~EFChargedTrack() {
00033 delete m_ntrk;
00034 delete m_acol;
00035 delete m_mbal;
00036 delete m_pmax1;
00037 delete m_pmax2;
00038 delete m_cost1;
00039 delete m_cost2;
00040 delete m_vr;
00041 delete m_vz;
00042 }
00043
00044 StatusCode EFChargedTrack::initialize(){
00045
00046 MsgStream log(msgSvc(), name());
00047 log << MSG::INFO << "in initialize()" << endreq;
00048
00049 IEFAlgorithm::initialize();
00050
00051 StatusCode sc;
00052 sc = m_HltStoreSvc->put("ntrk", m_ntrk);
00053 if ( sc.isFailure() ) {
00054 log << MSG::ERROR << "m_HltStoreSvc->put(ntrk) wrong" << endreq;
00055 return sc;
00056 }
00057 sc = m_HltStoreSvc->put("acol", m_acol);
00058 if ( sc.isFailure() ) {
00059 log << MSG::ERROR << "m_HltStoreSvc->put(acol) wrong" << endreq;
00060 return sc;
00061 }
00062 sc = m_HltStoreSvc->put("mbal", m_mbal);
00063 if ( sc.isFailure() ) {
00064 log << MSG::ERROR << "m_HltStoreSvc->put(mbal) wrong" << endreq;
00065 return sc;
00066 }
00067 sc = m_HltStoreSvc->put("pmax1", m_pmax1);
00068 if ( sc.isFailure() ) {
00069 log << MSG::ERROR << "m_HltStoreSvc->put(pmax1) wrong" << endreq;
00070 return sc;
00071 }
00072 sc = m_HltStoreSvc->put("pmax2", m_pmax2);
00073 if ( sc.isFailure() ) {
00074 log << MSG::ERROR << "m_HltStoreSvc->put(pmax2) wrong" << endreq;
00075 return sc;
00076 }
00077 sc = m_HltStoreSvc->put("cost1", m_cost1);
00078 if ( sc.isFailure() ) {
00079 log << MSG::ERROR << "m_HltStoreSvc->put(cost1) wrong" << endreq;
00080 return sc;
00081 }
00082 sc = m_HltStoreSvc->put("cost2", m_cost2);
00083 if ( sc.isFailure() ) {
00084 log << MSG::ERROR << "m_HltStoreSvc->put(cost2) wrong" << endreq;
00085 return sc;
00086 }
00087 sc = m_HltStoreSvc->put("vr", m_vr);
00088 if ( sc.isFailure() ) {
00089 log << MSG::ERROR << "m_HltStoreSvc->put(vr) wrong" << endreq;
00090 return sc;
00091 }
00092 sc = m_HltStoreSvc->put("vz", m_vz);
00093 if ( sc.isFailure() ) {
00094 log << MSG::ERROR << "m_HltStoreSvc->put(vz) wrong" << endreq;
00095 return sc;
00096 }
00097
00098 log << MSG::DEBUG << "finish initialize()" << endreq;
00099
00100 return StatusCode::SUCCESS;
00101 }
00102
00103 StatusCode EFChargedTrack::execute() {
00104
00105 reset();
00106
00107 MsgStream log(msgSvc(), name());
00108
00109
00110
00111 SmartDataPtr<Event::EventHeader> eventHeader(eventSvc(),"/Event/EventHeader");
00112 if (!eventHeader) {
00113 log << MSG::FATAL << "Could not find Event Header" << endreq;
00114 return( StatusCode::FAILURE);
00115 }
00116
00117
00118 SmartDataPtr<RecMdcTrackCol> mdcTrackCol(eventSvc(),"/Event/Recon/RecMdcTrackCol");
00119 if (!mdcTrackCol) {
00120 log << MSG::FATAL << "Could not find Mdc track collection!!" << endreq;
00121 return( StatusCode::FAILURE);
00122 }
00123
00124 RecMdcTrackCol::iterator iterTrk=mdcTrackCol->begin();
00125 double max1=0.;
00126 double max2=0.;
00127 double cost1=-999,cost2=-999,phi1=-999,phi2=-999;
00128 int hit1=-999,hit2=-999,shit1=-999,shit2=-999,charg1=-999,charg2=-999;
00129 double mdcbalance=-999;
00130
00131 unsigned int ntrk=mdcTrackCol->size();
00132 double kappa=-999,tanl=-999,sint=-999;
00133 double p=-999,theta=-999,phi=-999;
00134 double vz1=-999,vz2=-999,vr1=-999,vr2=-999;
00135
00136
00137 for(;iterTrk!= mdcTrackCol->end();iterTrk++) {
00138
00139
00140
00141 phi=(*iterTrk)->helix(1);
00142 kappa=(*iterTrk)->helix(2);
00143 tanl=(*iterTrk)->helix(4);
00144
00145 theta=0.5*3.1415926-atan(tanl);
00146 sint=sin(theta);
00147
00148
00149
00150 if(abs(kappa)>0.001&&abs(sint)>0.01){
00151 p=abs(1./kappa)/sint;
00152 }
00153 else{
00154 p=1000.;
00155 log << MSG::WARNING << "FastTrk=>"<<" kappa=" <<kappa<<"; sint="<<sint<<endreq;
00156 }
00157 if(p>=max1){
00158 max2=max1;
00159 cost2=cost1;
00160 phi2=phi1;
00161 max1=p;
00162 cost1=cos(theta);
00163 phi1=phi;
00164 vz2=vz1;
00165 vz1=(*iterTrk)->helix(3);
00166 vr2=vr1;
00167 vr1=(*iterTrk)->helix(0);
00168 hit2=hit1;
00169 hit1=(*iterTrk)->getNhits();
00170 shit2=shit1;
00171 shit1=(*iterTrk)->nster();
00172 charg2=charg1;
00173 charg1=(*iterTrk)->charge();
00174 }
00175 else if(p>max2){
00176 max2=p;
00177 cost2=cos(theta);
00178 phi2=phi;
00179 vz2=(*iterTrk)->helix(3);
00180 vr2=(*iterTrk)->helix(0);
00181 hit2=(*iterTrk)->getNhits();
00182 shit2=(*iterTrk)->nster();
00183 charg2=(*iterTrk)->charge();
00184 }
00185 log << MSG::DEBUG << "p=" <<p <<", "<<"theta="<<theta
00186 <<", phi="<<phi<<", vz="<<(*iterTrk)->helix(3)<<", vr="<<(*iterTrk)->helix(0)<<endreq;
00187 if(cos(theta)>0) mdcbalance +=1.;
00188 else if(cos(theta)<0) mdcbalance -=1.;
00189 }
00190 if(ntrk>=2) mdcbalance /= ntrk;
00191 else mdcbalance = 1;
00192
00193 double acol=180.;
00194 if(ntrk>=2){
00195 acol=180.-180./3.1415926*acos(cos(phi1)*sin(acos(cost1))*cos(phi2)*sin(acos(cost2))
00196 +sin(phi1)*sin(acos(cost1))*sin(phi2)*sin(acos(cost2))
00197 +cost1*cost2);
00198 }
00199
00200 log << MSG::INFO << "ntrk=" << ntrk << "; mdc balance=" <<mdcbalance
00201 <<"; pmax1="<< max1 <<"; pmax2="<< max2
00202 <<"; acol="<< acol<<"; cost1="<<cost1<<"; cost2="<<cost2<<endreq;
00203
00204
00205 m_ntrk->setValue(ntrk);
00206 m_acol->setValue(acol);
00207 m_mbal->setValue(mdcbalance);
00208 m_pmax1->setValue(max1);
00209 m_pmax2->setValue(max2);
00210 m_cost1->setValue(cost1);
00211 m_cost2->setValue(cost2);
00212 m_vr->setValue(vr1);
00213 m_vz->setValue(vz1);
00214
00215 m_ef->addToEFVec(ntrk, 6);
00216 m_ef->appToEFVec(max1, 7);
00217 m_ef->appToEFVec(cost1, 8);
00218 m_ef->appToEFVec(vz1, 9);
00219 m_ef->appToEFVec(phi1, 10);
00220 m_ef->appToEFVec(vr1, 11);
00221 m_ef->addToEFVec(hit1, 12);
00222 m_ef->addToEFVec(shit1, 13);
00223 m_ef->appToEFVec(max2, 14);
00224 m_ef->appToEFVec(cost2, 15);
00225 m_ef->appToEFVec(vz2, 16);
00226 m_ef->appToEFVec(phi2, 17);
00227 m_ef->appToEFVec(vr2, 18);
00228 m_ef->addToEFVec(hit2, 19);
00229 m_ef->addToEFVec(shit2, 20);
00230 m_ef->appToEFVec(acol, 21);
00231 m_ef->setVecBit(true, 0, 4);
00232 if(ntrk==0) m_ef->addToEFVec(1, 1);
00233 else if(ntrk==1) m_ef->addToEFVec(8, 1);
00234 else if(ntrk>=2) m_ef->addToEFVec(16,1);
00235 else m_ef->addToEFVec(0, 1);
00236
00237 m_ef->addToEFVec(charg1, 23);
00238 m_ef->addToEFVec(charg2, 24);
00239
00240 m_run=1;
00241
00242 return StatusCode::SUCCESS;
00243 }
00244
00245 StatusCode EFChargedTrack::finalize() {
00246 MsgStream log(msgSvc(), name());
00247 log << MSG::INFO << "in finalize()" << endmsg;
00248 return StatusCode::SUCCESS;
00249 }
00250
00251 void EFChargedTrack::reset() {
00252
00253 if(m_run){
00254 m_ntrk->reset();
00255 m_acol->reset();
00256 m_mbal->reset();
00257 m_pmax1->reset();
00258 m_pmax2->reset();
00259 m_cost1->reset();
00260 m_cost2->reset();
00261 m_vr->reset();
00262 m_vz->reset();
00263 m_run=0;
00264 }
00265 return;
00266 }