Main Page | Namespace List | Class Hierarchy | Alphabetical List | Class List | Directories | File List | Namespace Members | Class Members | File Members | Related Pages

JsiLL Class Reference

#include <JsiLL.h>

List of all members.

Public Member Functions

StatusCode execute ()
StatusCode execute ()
StatusCode finalize ()
StatusCode finalize ()
StatusCode initialize ()
StatusCode initialize ()
 JsiLL (const std::string &name, ISvcLocator *pSvcLocator)
 JsiLL (const std::string &name, ISvcLocator *pSvcLocator)

Private Attributes

NTuple::Item< double > m_chisq
NTuple::Item< double > m_chisq
double m_cthcut
double m_energyThreshold
NTuple::Item< long > m_event
NTuple::Item< long > m_event
double m_gammaAngCut
NTuple::Item< double > m_mLambda
NTuple::Item< double > m_mLambda
NTuple::Item< double > m_mLambdabar
NTuple::Item< double > m_mLambdabar
NTuple::Item< double > m_pLambda
NTuple::Item< double > m_pLambda
NTuple::Item< double > m_pLambdabar
NTuple::Item< double > m_pLambdabar
NTuple::Item< long > m_runNo
NTuple::Item< long > m_runNo
ITHistSvc * m_thsvc
ITHistSvc * m_thsvc
NTuple::Tuple * m_tuple
NTuple::Tuple * m_tuple
double m_vr0cut
double m_vr1cut
double m_vz0cut
double m_vz1cut

Constructor & Destructor Documentation

JsiLL::JsiLL (  const std::string &  name,
ISvcLocator *  pSvcLocator
) 
  

00059                                                             :
00060   Algorithm(name, pSvcLocator) {
00061   
00062   //Declare the properties  
00063   declareProperty("Vr0cut", m_vr0cut=5.0);
00064   declareProperty("Vz0cut", m_vz0cut=20.0);
00065   declareProperty("Vr1cut", m_vr1cut=1.0);
00066   declareProperty("Vz1cut", m_vz1cut=5.0);
00067   declareProperty("Vctcut", m_cthcut=0.93);
00068   declareProperty("EnergyThreshold", m_energyThreshold=0.04);
00069   declareProperty("GammaAngCut", m_gammaAngCut=20.0);
00070 }

JsiLL::JsiLL (  const std::string &  name,
ISvcLocator *  pSvcLocator
) 
 

Member Function Documentation

StatusCode JsiLL::execute (   ) 
 

StatusCode JsiLL::execute (   ) 
  

00124                           {
00125   
00126   MsgStream log(msgSvc(), name());
00127   log << MSG::INFO << "in execute()" << endreq;
00128 
00129   // DQA
00130   // Add the line below at the beginning of execute()
00131   //
00132   setFilterPassed(false);
00133 
00134   SmartDataPtr<Event::EventHeader> eventHeader(eventSvc(),"/Event/EventHeader");
00135   m_runNo=eventHeader->runNumber();
00136   m_event=eventHeader->eventNumber();
00137   log << MSG::DEBUG <<"run, evtnum = "
00138       << m_runNo << " , "
00139       << m_event <<endreq;
00140 
00141 
00142   SmartDataPtr<EvtRecEvent> evtRecEvent(eventSvc(), EventModel::EvtRec::EvtRecEvent);
00143   log << MSG::DEBUG <<"ncharg, nneu, tottks = " 
00144       << evtRecEvent->totalCharged() << " , "
00145       << evtRecEvent->totalNeutral() << " , "
00146       << evtRecEvent->totalTracks() <<endreq;
00147 
00148   SmartDataPtr<EvtRecTrackCol> evtRecTrkCol(eventSvc(), EventModel::EvtRec::EvtRecTrackCol);
00149 
00150   if(evtRecEvent->totalNeutral()>100) {
00151       return StatusCode::SUCCESS;
00152   }
00153 
00154   Vint iGood, iplus, iminus;
00155   iGood.clear();
00156   iplus.clear();
00157   iminus.clear();
00158   Vp4 ppip, ppim;
00159   ppip.clear();
00160   ppim.clear();
00161 
00162   Hep3Vector xorigin(0,0,0);
00163  
00164   IVertexDbSvc*  vtxsvc;
00165   Gaudi::svcLocator()->service("VertexDbSvc", vtxsvc);
00166   if(vtxsvc->isVertexValid()){ 
00167   double* dbv = vtxsvc->PrimaryVertex();
00168   double*  vv = vtxsvc->SigmaPrimaryVertex();
00169     xorigin.setX(dbv[0]);
00170     xorigin.setY(dbv[1]);
00171     xorigin.setZ(dbv[2]);
00172   }
00173 
00174   int nCharge = 0;
00175   for(int i = 0; i < evtRecEvent->totalCharged(); i++){
00176     EvtRecTrackIterator itTrk=evtRecTrkCol->begin() + i;
00177     if(!(*itTrk)->isMdcTrackValid()) continue;
00178     if (!(*itTrk)->isMdcKalTrackValid()) continue;
00179     RecMdcTrack *mdcTrk = (*itTrk)->mdcTrack();
00180 
00181     double pch =mdcTrk->p();
00182     double x0  =mdcTrk->x(); 
00183     double y0  =mdcTrk->y(); 
00184     double z0  =mdcTrk->z(); 
00185     double phi0=mdcTrk->helix(1);
00186     double xv=xorigin.x();
00187     double yv=xorigin.y();
00188     double Rxy=fabs((x0-xv)*cos(phi0)+(y0-yv)*sin(phi0));
00189 
00190     if(fabs(z0) >= m_vz0cut) continue;
00191     if(Rxy >= m_vr0cut) continue;
00192 //    if(fabs(m_Vct)>=m_cthcut) continue;
00193 
00194     // DQA
00195     TH1 *h(0);
00196     if (m_thsvc->getHist("/DQAHist/JsiLL/hrxy", h).isSuccess()) {
00197         h->Fill(Rxy);
00198     } else {
00199         log << MSG::ERROR << "Couldn't retrieve hrxy" << endreq;
00200     }
00201     if (m_thsvc->getHist("/DQAHist/JsiLL/hz", h).isSuccess()) {
00202         h->Fill(z0);
00203     } else {
00204         log << MSG::ERROR << "Couldn't retrieve hz" << endreq;
00205     }
00206 
00207 //     iGood.push_back((*itTrk)->trackId());
00208     iGood.push_back(i);
00209     nCharge += mdcTrk->charge();
00210     if (mdcTrk->charge() > 0) {
00211         iplus.push_back(i);
00212     } else {
00213         iminus.push_back(i);
00214     }
00215   }
00216   
00217   //
00218   // Finish Good Charged Track Selection
00219   //
00220   int nGood = iGood.size();
00221   log << MSG::DEBUG << "ngood, totcharge = " << nGood << " , " << nCharge << endreq;
00222   if((nGood != 4)||(nCharge!=0)){
00223       return StatusCode::SUCCESS;
00224   }
00225 
00226 
00227 //  for(int i = 0; i < nGood; i++) {
00228 //    EvtRecTrackIterator itTrk = evtRecTrkCol->begin() + iGood[i];
00229 //    if(!(*itTrk)->isMdcKalTrackValid())
00230 //      return StatusCode::SUCCESS;
00231 //    }
00232   
00233     int pidp0 = 5, pidp1 = 3, pidm0 = 5, pidm1 = 3;  //PID
00234   
00235     RecMdcKalTrack *itTrkp   = (*(evtRecTrkCol->begin() + iplus[0]))->mdcKalTrack(); 
00236     RecMdcKalTrack *itTrkpip = (*(evtRecTrkCol->begin() + iplus[1]))->mdcKalTrack();
00237     RecMdcKalTrack *itTrkpb  = (*(evtRecTrkCol->begin() + iminus[0]))->mdcKalTrack();
00238     RecMdcKalTrack *itTrkpim = (*(evtRecTrkCol->begin() + iminus[1]))->mdcKalTrack();
00239 
00240 
00241 
00242     if (itTrkp->p() < itTrkpip->p()){
00243         itTrkp   = (*(evtRecTrkCol->begin() + iplus[1]))->mdcKalTrack();
00244         itTrkpip = (*(evtRecTrkCol->begin() + iplus[0]))->mdcKalTrack(); 
00245         pidp0 = 3;
00246         pidp1 = 5;
00247     }     
00248     if (itTrkpb->p() < itTrkpim->p()){
00249         itTrkpb  = (*(evtRecTrkCol->begin() + iminus[1]))->mdcKalTrack();
00250         itTrkpim = (*(evtRecTrkCol->begin() + iminus[0]))->mdcKalTrack(); 
00251         pidm0 = 3;
00252         pidm1 = 5;
00253     }     
00254     if (itTrkp->p()   < 0.7 || itTrkp->p()  >1.1) return StatusCode::SUCCESS; 
00255     if (itTrkpb->p()  < 0.7 || itTrkpb->p() >1.1) return StatusCode::SUCCESS; 
00256     if (itTrkpip->p() > 0.35) return StatusCode::SUCCESS; 
00257     if (itTrkpim->p() > 0.35) return StatusCode::SUCCESS; 
00258     
00259 
00260 
00261     itTrkp->setPidType(RecMdcKalTrack::proton);  
00262     itTrkpb->setPidType(RecMdcKalTrack::proton);  
00263     itTrkpip->setPidType(RecMdcKalTrack::pion);  
00264     itTrkpim->setPidType(RecMdcKalTrack::pion);  
00265 
00266 
00267    m_chisq = 9999.0;
00268   // Vertex Fit
00269   HepPoint3D vx(0., 0., 0.);
00270   HepSymMatrix Evx(3, 0);
00271   double bx = 1E+6; 
00272   double by = 1E+6; 
00273   double bz = 1E+6; 
00274   Evx[0][0] = bx*bx;
00275   Evx[1][1] = by*by;
00276   Evx[2][2] = bz*bz;
00277 
00278   VertexParameter vxpar;
00279   vxpar.setVx(vx);  
00280   vxpar.setEvx(Evx);
00281 
00282   VertexFit* vtxfita0 = VertexFit::instance();
00283   SecondVertexFit *vtxfita = SecondVertexFit::instance();
00284   VertexFit* vtxfitb0 = VertexFit::instance();
00285   SecondVertexFit *vtxfitb = SecondVertexFit::instance();
00286   VertexFit* vtxfit = VertexFit::instance();
00287   
00288     WTrackParameter wpip = WTrackParameter(mpi, itTrkpip->getZHelix(), itTrkpip->getZError());
00289     WTrackParameter wpim = WTrackParameter(mpi, itTrkpim->getZHelix(), itTrkpim->getZError());
00290     WTrackParameter wp   = WTrackParameter(mproton, itTrkp->getZHelixP(), itTrkp->getZErrorP());
00291     WTrackParameter wpb  = WTrackParameter(mproton, itTrkpb->getZHelixP(), itTrkpb->getZErrorP());
00292 
00293   vtxfita0->init();
00294   vtxfita0->AddTrack(0, wp);
00295   vtxfita0->AddTrack(1, wpim);
00296   vtxfita0->AddVertex(0, vxpar, 0, 1);
00297   if(!vtxfita0->Fit(0)) return StatusCode::SUCCESS;
00298  vtxfita0->Swim(0);  
00299  vtxfita0->BuildVirtualParticle(0);
00300   vtxfita->init();
00301   vtxfita->AddTrack(0, vtxfita0->wVirtualTrack(0));
00302   vtxfita->setVpar(vtxfita0->vpar(0));
00303   if(!vtxfita->Fit()) return StatusCode::SUCCESS; 
00304 
00305   WTrackParameter wLambda    = vtxfita->wpar();   
00306 
00307   vtxfitb0->init();
00308   vtxfitb0->AddTrack(0, wpb);
00309   vtxfitb0->AddTrack(1, wpip);
00310   vtxfitb0->AddVertex(0, vxpar, 0, 1);
00311   if(!vtxfitb0->Fit(0)) return StatusCode::SUCCESS;
00312   vtxfitb0->Swim(0); 
00313   vtxfitb0->BuildVirtualParticle(0);
00314 
00315   vtxfitb->init();
00316   vtxfitb->AddTrack(0, vtxfitb0->wVirtualTrack(0));
00317   vtxfitb->setVpar(vtxfitb0->vpar(0));
00318   if(!vtxfitb->Fit()) return StatusCode::SUCCESS; 
00319 
00320   WTrackParameter wLambdabar = vtxfitb->wpar();
00321      
00322   vtxfit->init();
00323   vtxfit->AddTrack(0,   wLambda);
00324   vtxfit->AddTrack(1,   wLambdabar);
00325   vtxfit->AddVertex(0, vxpar,0, 1); 
00326   if(!vtxfit->Fit(0)) return StatusCode::SUCCESS;
00327   vtxfit->Swim(0);
00328   WTrackParameter wwLambda    = vtxfit->wtrk(0);
00329   WTrackParameter wwLambdabar = vtxfit->wtrk(1);
00330 
00331 
00332   // Kinamatic Fit
00333 
00334    KinematicFit* kmfit = KinematicFit::instance();
00335 
00336 
00337     HepLorentzVector ecms(0.034065,0.0,0.0,3.0969);
00338     const Hep3Vector u_cms = -ecms.boostVector();
00339 
00340     kmfit->init();
00341     kmfit->AddTrack(0, wwLambda);
00342     kmfit->AddTrack(1, wwLambdabar);
00343     kmfit->AddFourMomentum(0, ecms);
00344 
00345     if(!kmfit->Fit())  return StatusCode::SUCCESS;
00346       m_chisq = kmfit->chisq();
00347       if(m_chisq > 50)  return StatusCode::SUCCESS;
00348       HepLorentzVector kmf_pLambda    = kmfit->pfit(0);
00349       HepLorentzVector kmf_pLambdabar = kmfit->pfit(1);
00350 
00351   kmf_pLambda.boost(u_cms);
00352   kmf_pLambdabar.boost(u_cms);
00353   m_mLambda    = kmf_pLambda.m();
00354   m_mLambdabar = kmf_pLambdabar.m();
00355   m_pLambda    = kmf_pLambda.rho();
00356   m_pLambdabar = kmf_pLambdabar.rho();
00357 
00358   if(fabs(m_mLambda-1.1157)>0.003||fabs(m_mLambdabar-1.1157)>0.003) return StatusCode::SUCCESS;
00359   // finale selection
00360   
00361   m_tuple->write();
00362 //      return StatusCode::SUCCESS;
00363 
00365   // DQA
00366   // set tag and quality
00367 
00368   // Pid: 1 - electron, 2 - muon, 3 - pion, 4 - kaon, 5 - proton
00369   (*(evtRecTrkCol->begin()+iplus[0]))->setPartId(pidp0);
00370   (*(evtRecTrkCol->begin()+iplus[1]))->setPartId(pidp1);
00371   (*(evtRecTrkCol->begin()+iminus[0]))->setPartId(pidm0);
00372   (*(evtRecTrkCol->begin()+iminus[1]))->setPartId(pidm1);
00373   // Quality: defined by whether dE/dx or TOF is used to identify particle
00374   // 0 - no dE/dx, no TOF (can be used for dE/dx and TOF calibration)
00375   // 1 - only dE/dx (can be used for TOF calibration)
00376   // 2 - only TOF (can be used for dE/dx calibration)
00377   // 3 - Both dE/dx and TOF
00378   (*(evtRecTrkCol->begin()+iplus[0]))->setQuality(0);
00379   (*(evtRecTrkCol->begin()+iplus[1]))->setQuality(0);
00380   (*(evtRecTrkCol->begin()+iminus[0]))->setQuality(0);
00381   (*(evtRecTrkCol->begin()+iminus[1]))->setQuality(0);
00382 
00383   // DQA
00384   // Add the line below at the end of execute(), (before return)
00385   //
00386   setFilterPassed(true);
00388   
00389   return StatusCode::SUCCESS;
00390 
00391 }

StatusCode JsiLL::finalize (   ) 
 

StatusCode JsiLL::finalize (   ) 
  

00394                            {
00395 
00396   MsgStream log(msgSvc(), name());
00397   log << MSG::INFO << "in finalize()" << endmsg;
00398   return StatusCode::SUCCESS;
00399 }

StatusCode JsiLL::initialize (   ) 
 

StatusCode JsiLL::initialize (   ) 
  

00074                             {
00075   MsgStream log(msgSvc(), name());
00076 
00077   log << MSG::INFO << "in initialize()" << endmsg;
00078   StatusCode status;
00079 
00080   // DQA
00081   // The first directory specifier must be "DQAFILE"
00082   // The second is the control sample name, the first letter is in upper format. eg. "Rhopi"
00083   NTuplePtr nt(ntupleSvc(), "DQAFILE/JsiLL");
00084   if ( nt ) m_tuple = nt;
00085    else {
00086       m_tuple = ntupleSvc()->book("DQAFILE/JsiLL", CLID_ColumnWiseTuple, "JsiLL ntuple");
00087       if( m_tuple ) {
00088           status = m_tuple->addItem("runNo", m_runNo);
00089           status = m_tuple->addItem("event", m_event);
00090           status = m_tuple->addItem("chisq", m_chisq);
00091           status = m_tuple->addItem("mLambda", m_mLambda);
00092           status = m_tuple->addItem("mLambdabar", m_mLambdabar);
00093           status = m_tuple->addItem("pLambda", m_pLambda);
00094           status = m_tuple->addItem("pLambdabar", m_pLambdabar);
00095 
00096       } else {
00097           log << MSG::ERROR << "Can not book N-tuple:" << long(m_tuple) << endreq;
00098       }
00099   }
00100 
00101   if(service("THistSvc", m_thsvc).isFailure()) {
00102     log << MSG::ERROR << "Couldn't get THistSvc" << endreq;
00103     return StatusCode::FAILURE;
00104   }
00105 
00106   // "DQAHist" is fixed
00107   // "Rhopi" is control sample name, just as ntuple case.
00108   TH1F* hrxy = new TH1F("Rxy", "Rxy distribution", 110, -1., 10.);
00109   if(m_thsvc->regHist("/DQAHist/JsiLL/hrxy", hrxy).isFailure()) {
00110       log << MSG::ERROR << "Couldn't register Rxy" << endreq;
00111   }
00112   TH1F* hz = new TH1F("z", "z distribution", 200, -100., 100.);
00113   if(m_thsvc->regHist("/DQAHist/JsiLL/hz", hz).isFailure()) {
00114       log << MSG::ERROR << "Couldn't register z" << endreq;
00115   }
00116   
00117   log << MSG::INFO << "successfully return from initialize()" <<endmsg;
00118   return StatusCode::SUCCESS;
00119 
00120 
00121 }

Member Data Documentation

NTuple::Item<double> JsiLL::m_chisq [private]
 

NTuple::Item<double> JsiLL::m_chisq [private]
 

double JsiLL::m_cthcut [private]
 

double JsiLL::m_energyThreshold [private]
 

NTuple::Item<long> JsiLL::m_event [private]
 

NTuple::Item<long> JsiLL::m_event [private]
 

double JsiLL::m_gammaAngCut [private]
 

NTuple::Item<double> JsiLL::m_mLambda [private]
 

NTuple::Item<double> JsiLL::m_mLambda [private]
 

NTuple::Item<double> JsiLL::m_mLambdabar [private]
 

NTuple::Item<double> JsiLL::m_mLambdabar [private]
 

NTuple::Item<double> JsiLL::m_pLambda [private]
 

NTuple::Item<double> JsiLL::m_pLambda [private]
 

NTuple::Item<double> JsiLL::m_pLambdabar [private]
 

NTuple::Item<double> JsiLL::m_pLambdabar [private]
 

NTuple::Item<long> JsiLL::m_runNo [private]
 

NTuple::Item<long> JsiLL::m_runNo [private]
 

ITHistSvc* JsiLL::m_thsvc [private]
 

ITHistSvc* JsiLL::m_thsvc [private]
 

NTuple::Tuple* JsiLL::m_tuple [private]
 

NTuple::Tuple* JsiLL::m_tuple [private]
 

double JsiLL::m_vr0cut [private]
 

double JsiLL::m_vr1cut [private]
 

double JsiLL::m_vz0cut [private]
 

double JsiLL::m_vz1cut [private]
 

The documentation for this class was generated from the following files:
Generated on Wed Feb 2 16:19:05 2011 for BOSS6.5.5 by  doxygen 1.3.9.1