XtInteMdcCalib Class Reference

#include <XtInteMdcCalib.h>

Inheritance diagram for XtInteMdcCalib:

MdcCalib List of all members.

Public Member Functions

 XtInteMdcCalib ()
 ~XtInteMdcCalib ()
void initialize (TObjArray *hlist, IMdcGeomSvc *mdcGeomSvc, IMdcCalibFunSvc *mdcFunSvc, IMdcUtilitySvc *mdcUtilitySvc)
void setParam (MdcCalParams &param)
int fillHist (MdcCalEvent *event)
int updateConst (MdcCalibConst *calconst)
void clear ()

Private Attributes

MdcCalParams m_param
TObjArray * m_hlist
IMdcGeomSvcm_mdcGeomSvc
IMdcCalibFunSvcm_mdcFunSvc
IMdcUtilitySvcm_mdcUtilitySvc
bool m_fgIni
int m_nMaxGrPoint
double m_docaMax [MdcCalNLayer]
double m_tbinWid [MdcCalNLayer][3]
double m_tbinLim [MdcCalNLayer][4]
TFolder * m_fdPf
TGraph * m_grXt [MdcCalNLayer][NENTR][2]
TProfile * m_pfNear [MdcCalNLayer][NENTR][2]
TProfile * m_pfMid [MdcCalNLayer][NENTR][2]
TProfile * m_pfFar [MdcCalNLayer][NENTR][2]

Static Private Attributes

static const int NENTR = 18
static const int NTBIN_INN = 72
static const int NTBIN_OUT = 79

Detailed Description

Definition at line 13 of file XtInteMdcCalib.h.


Constructor & Destructor Documentation

XtInteMdcCalib::XtInteMdcCalib (  ) 

Definition at line 19 of file XtInteMdcCalib.cxx.

References m_fgIni, m_nMaxGrPoint, m_tbinLim, m_tbinWid, and MdcCalNLayer.

00019                               {
00020      m_fgIni = false;
00021      m_nMaxGrPoint = 5000;
00022      for(int lay=0; lay<MdcCalNLayer; lay++){
00023           m_tbinWid[lay][0] = 5.0;
00024           m_tbinWid[lay][1] = 10.0;
00025           m_tbinWid[lay][2] = 20.0;
00026 
00027           m_tbinLim[lay][0] = -10.0;
00028           m_tbinLim[lay][1] = 30.0;
00029           if(lay < 8) m_tbinLim[lay][2] = 210.0;
00030           else m_tbinLim[lay][2] = 350.0;
00031           m_tbinLim[lay][3] = 990.0;
00032      }
00033 }

XtInteMdcCalib::~XtInteMdcCalib (  ) 

Definition at line 35 of file XtInteMdcCalib.cxx.

00035                                {
00036 }


Member Function Documentation

void XtInteMdcCalib::clear (  )  [virtual]

Implements MdcCalib.

Definition at line 38 of file XtInteMdcCalib.cxx.

References MdcCalib::clear(), m_fdPf, m_grXt, m_pfFar, m_pfMid, m_pfNear, MdcCalNLayer, and NENTR.

00038                           {
00039      for(int lay=0; lay<MdcCalNLayer; lay++){
00040           for(int iEntr=0; iEntr<NENTR; iEntr++){
00041                for(int lr=0; lr<2; lr++){
00042                     delete m_pfNear[lay][iEntr][lr];
00043                     delete m_pfMid[lay][iEntr][lr];
00044                     delete m_pfFar[lay][iEntr][lr];
00045                     delete m_grXt[lay][iEntr][lr];
00046                }
00047           }
00048      }
00049      delete m_fdPf;
00050      MdcCalib::clear();
00051 }

int XtInteMdcCalib::fillHist ( MdcCalEvent event  )  [virtual]

Implements MdcCalib.

Definition at line 104 of file XtInteMdcCalib.cxx.

References Bes_Common::DEBUG, MdcCalParams::drCut, MdcCalParams::dzCut, MdcCalib::fillHist(), MdcCalRecTrk::getDr(), MdcCalRecTrk::getDz(), MdcCalEvent::getEsCutFlag(), MdcCalEvent::getRecTrk(), genRecEmupikp::i, m_docaMax, m_fgIni, m_grXt, m_mdcFunSvc, m_nMaxGrPoint, m_param, m_pfFar, m_pfMid, m_pfNear, m_tbinLim, MdcCalParams::maxDocaInner, MdcCalParams::maxDocaOuter, MdcCalNLayer, msgSvc(), MdcCalParams::nHitLayCut, and MdcCalParams::resiCut.

00104                                               {
00105      IMessageSvc* msgSvc;
00106      Gaudi::svcLocator() -> service("MessageSvc", msgSvc);
00107      MsgStream log(msgSvc, "XtInteMdcCalib");
00108      log << MSG::DEBUG << "XtInteMdcCalib::fillHist()" << endreq;
00109 
00110      MdcCalib::fillHist(event);
00111 
00112      // get EsTimeCol
00113      bool esCutFg = event->getEsCutFlag();
00114      if( ! esCutFg ) return -1;
00115 
00116      int i;
00117      int k;
00118      int lay;
00119      int iLR;
00120      int iEntr;
00121 
00122      double dr;
00123      double dz;
00124      double doca;
00125      double tdr;
00126      double resi;
00127      double entrance;
00128 
00129      int nhitlay;
00130      bool fgHitLay[MdcCalNLayer];
00131      bool fgTrk;
00132      
00133      if(! m_fgIni){
00134           for(lay=0; lay<MdcCalNLayer; lay++){
00135                if(lay < 8) m_docaMax[lay] = m_param.maxDocaInner;
00136                else m_docaMax[lay] = m_param.maxDocaOuter;
00137           }
00138           m_fgIni = true;
00139      }
00140 
00141      MdcCalRecTrk* rectrk;
00142      MdcCalRecHit* rechit;
00143 
00144      int nhit;
00145      int ntrk = event -> getNTrk();
00146      for(i=0; i<ntrk; i++){
00147           fgTrk = true;
00148           rectrk = event->getRecTrk(i);
00149           nhit = rectrk -> getNHits();
00150 
00151           // dr cut
00152           dr = rectrk->getDr();
00153           if(fabs(dr) > m_param.drCut) continue;
00154         
00155           // dz cut
00156           dz = rectrk->getDz();
00157           if(fabs(dz) > m_param.dzCut) continue;
00158 
00159           for(lay=0; lay<MdcCalNLayer; lay++){
00160                fgHitLay[lay] = false;
00161           }
00162 
00163           for(k=0; k<nhit; k++){
00164                rechit = rectrk -> getRecHit(k);
00165                lay = rechit -> getLayid();
00166                doca = rechit -> getDocaInc();
00167                resi = rechit -> getResiInc();
00168                fgHitLay[lay] = true;
00169 
00170 //             if( (fabs(doca) > m_docaMax[lay]) || 
00171 //                 (fabs(resi) > m_param.resiCut[lay]) ){
00172 //                  fgTrk = false;
00173 //             }
00174           }
00175           if(! fgTrk) continue;
00176 
00177           nhitlay = 0;
00178           for(lay=0; lay<MdcCalNLayer; lay++){
00179                if(fgHitLay[lay]) nhitlay++;
00180           }
00181           if(nhitlay < m_param.nHitLayCut) continue;
00182 
00183           for(k=0; k<nhit; k++){
00184                rechit = rectrk -> getRecHit(k);
00185                lay = rechit -> getLayid();
00186                doca = rechit -> getDocaInc();
00187                resi = rechit -> getResiInc();
00188                iLR = rechit -> getLR();
00189                entrance = rechit -> getEntra();
00190                tdr = rechit -> getTdrift();
00191 
00192                if( (fabs(doca) > m_docaMax[lay]) || 
00193                    (fabs(resi) > m_param.resiCut[lay]) ){
00194                     continue;
00195                }
00196 
00197                if(0 == lay){
00198                     if( ! fgHitLay[1] ) continue;
00199                } else if(42 == lay){
00200                     if( ! fgHitLay[41] ) continue;
00201                } else{
00202                     if( (!fgHitLay[lay-1]) && (!fgHitLay[lay+1]) ) continue;
00203                }
00204 
00205                iEntr = m_mdcFunSvc -> getXtEntrIndex(entrance);
00206                int npoint = m_grXt[lay][iEntr][iLR]->GetN();
00207                if(npoint < m_nMaxGrPoint) m_grXt[lay][iEntr][iLR]->SetPoint(npoint, tdr, fabs(doca));
00208 
00209                if(tdr<m_tbinLim[lay][0]) continue;
00210                else if(tdr<m_tbinLim[lay][1])  m_pfNear[lay][iEntr][iLR]->Fill(tdr, fabs(doca));
00211                else if(tdr<m_tbinLim[lay][2])  m_pfMid[lay][iEntr][iLR]->Fill(tdr, fabs(doca));
00212                else if(tdr<m_tbinLim[lay][3]) m_pfFar[lay][iEntr][iLR]->Fill(tdr, fabs(doca));
00213           } // hit loop
00214      } // track loop
00215      return 1;
00216 }

void XtInteMdcCalib::initialize ( TObjArray *  hlist,
IMdcGeomSvc mdcGeomSvc,
IMdcCalibFunSvc mdcFunSvc,
IMdcUtilitySvc mdcUtilitySvc 
) [virtual]

Implements MdcCalib.

Definition at line 53 of file XtInteMdcCalib.cxx.

References Bes_Common::INFO, MdcCalib::initialize(), m_fdPf, m_grXt, m_hlist, m_mdcFunSvc, m_mdcGeomSvc, m_mdcUtilitySvc, m_pfFar, m_pfMid, m_pfNear, m_tbinLim, m_tbinWid, MdcCalNLayer, msgSvc(), and NENTR.

00054                                                                                            {
00055      IMessageSvc* msgSvc;
00056      Gaudi::svcLocator() -> service("MessageSvc", msgSvc);
00057      MsgStream log(msgSvc, "XtInteMdcCalib");
00058      log << MSG::INFO << "XtInteMdcCalib::initialize()" << endreq;
00059 
00060      m_hlist = hlist;
00061      m_mdcGeomSvc = mdcGeomSvc;
00062      m_mdcFunSvc = mdcFunSvc;
00063      m_mdcUtilitySvc = mdcUtilitySvc;
00064 
00065      MdcCalib::initialize(m_hlist, m_mdcGeomSvc, m_mdcFunSvc, m_mdcUtilitySvc);
00066 
00067      m_fdPf = new TFolder("fdProfile", "fdProfile");
00068      m_hlist -> Add(m_fdPf);
00069 
00070      char hname[200];
00071      for(int lay=0; lay<MdcCalNLayer; lay++){
00072           for(int iEntr=0; iEntr<NENTR; iEntr++){
00073                for(int lr=0; lr<2; lr++){
00074                     sprintf(hname, "xt%02d_%02d_%d_gr", lay, iEntr, lr);
00075                     m_grXt[lay][iEntr][lr] = new TGraph();
00076                     m_grXt[lay][iEntr][lr]->SetName(hname);
00077                     m_fdPf->Add(m_grXt[lay][iEntr][lr]);
00078 
00079                     int xbinN = (int)((m_tbinLim[lay][1] - m_tbinLim[lay][0])/m_tbinWid[lay][0] + 0.5);
00080                     sprintf(hname, "xt%02d_%02d_%d_near", lay, iEntr, lr);
00081                     m_pfNear[lay][iEntr][lr] = new TProfile(hname, hname, xbinN, m_tbinLim[lay][0], m_tbinLim[lay][1]);
00082                     m_fdPf->Add(m_pfNear[lay][iEntr][lr]);
00083 
00084                     int xbinM = (int)((m_tbinLim[lay][2] - m_tbinLim[lay][1])/m_tbinWid[lay][1] + 0.5);
00085                     sprintf(hname, "xt%02d_%02d_%d_mid", lay, iEntr, lr);
00086                     m_pfMid[lay][iEntr][lr] = new TProfile(hname, hname, xbinM, m_tbinLim[lay][1], m_tbinLim[lay][2]);
00087                     m_fdPf->Add(m_pfMid[lay][iEntr][lr]);
00088 
00089                     int xbinF = (int)((m_tbinLim[lay][3] - m_tbinLim[lay][2])/m_tbinWid[lay][2] + 0.5);
00090                     sprintf(hname, "xt%02d_%02d_%d_far", lay, iEntr, lr);
00091                     m_pfFar[lay][iEntr][lr] = new TProfile(hname, hname, xbinF, m_tbinLim[lay][2], m_tbinLim[lay][3]);
00092                     m_fdPf->Add(m_pfFar[lay][iEntr][lr]);
00093 
00094 //                  if((0==iEntr)&&(0==lr))
00095 //                  cout << setw(5) << lay
00096 //                       << setw(5) << xbinN << setw(10) << m_tbinLim[lay][0] << setw(10) << m_tbinLim[lay][1]
00097 //                       << setw(5) << xbinM << setw(10) << m_tbinLim[lay][1] << setw(10) << m_tbinLim[lay][2]
00098 //                       << setw(5) << xbinF << setw(10) << m_tbinLim[lay][2] << setw(10) << m_tbinLim[lay][3] << endl;
00099                }
00100           }
00101      }
00102 }

void XtInteMdcCalib::setParam ( MdcCalParams param  )  [inline, virtual]

Implements MdcCalib.

Definition at line 50 of file XtInteMdcCalib.h.

References m_param, and MdcCalib::setParam().

00050                                                        {
00051      MdcCalib::setParam(param);
00052      m_param = param;
00053 }

int XtInteMdcCalib::updateConst ( MdcCalibConst calconst  )  [virtual]

Implements MdcCalib.

Definition at line 218 of file XtInteMdcCalib.cxx.

References Bes_Common::INFO, msgSvc(), and MdcCalib::updateConst().

00218                                                       {
00219      IMessageSvc* msgSvc;
00220      Gaudi::svcLocator() -> service("MessageSvc", msgSvc);
00221      MsgStream log(msgSvc, "XtInteMdcCalib");
00222      log << MSG::INFO << "XtInteMdcCalib::updateConst()" << endreq;
00223 
00224      MdcCalib::updateConst(calconst);
00225 }


Member Data Documentation

double XtInteMdcCalib::m_docaMax[MdcCalNLayer] [private]

Reimplemented from MdcCalib.

Definition at line 39 of file XtInteMdcCalib.h.

Referenced by fillHist().

TFolder* XtInteMdcCalib::m_fdPf [private]

Definition at line 43 of file XtInteMdcCalib.h.

Referenced by clear(), and initialize().

bool XtInteMdcCalib::m_fgIni [private]

Reimplemented from MdcCalib.

Definition at line 37 of file XtInteMdcCalib.h.

Referenced by fillHist(), and XtInteMdcCalib().

TGraph* XtInteMdcCalib::m_grXt[MdcCalNLayer][NENTR][2] [private]

Definition at line 44 of file XtInteMdcCalib.h.

Referenced by clear(), fillHist(), and initialize().

TObjArray* XtInteMdcCalib::m_hlist [private]

Reimplemented from MdcCalib.

Definition at line 32 of file XtInteMdcCalib.h.

Referenced by initialize().

IMdcCalibFunSvc* XtInteMdcCalib::m_mdcFunSvc [private]

Reimplemented from MdcCalib.

Definition at line 34 of file XtInteMdcCalib.h.

Referenced by fillHist(), and initialize().

IMdcGeomSvc* XtInteMdcCalib::m_mdcGeomSvc [private]

Reimplemented from MdcCalib.

Definition at line 33 of file XtInteMdcCalib.h.

Referenced by initialize().

IMdcUtilitySvc* XtInteMdcCalib::m_mdcUtilitySvc [private]

Reimplemented from MdcCalib.

Definition at line 35 of file XtInteMdcCalib.h.

Referenced by initialize().

int XtInteMdcCalib::m_nMaxGrPoint [private]

Definition at line 38 of file XtInteMdcCalib.h.

Referenced by fillHist(), and XtInteMdcCalib().

MdcCalParams XtInteMdcCalib::m_param [private]

Reimplemented from MdcCalib.

Definition at line 30 of file XtInteMdcCalib.h.

Referenced by fillHist(), and setParam().

TProfile* XtInteMdcCalib::m_pfFar[MdcCalNLayer][NENTR][2] [private]

Definition at line 47 of file XtInteMdcCalib.h.

Referenced by clear(), fillHist(), and initialize().

TProfile* XtInteMdcCalib::m_pfMid[MdcCalNLayer][NENTR][2] [private]

Definition at line 46 of file XtInteMdcCalib.h.

Referenced by clear(), fillHist(), and initialize().

TProfile* XtInteMdcCalib::m_pfNear[MdcCalNLayer][NENTR][2] [private]

Definition at line 45 of file XtInteMdcCalib.h.

Referenced by clear(), fillHist(), and initialize().

double XtInteMdcCalib::m_tbinLim[MdcCalNLayer][4] [private]

Definition at line 41 of file XtInteMdcCalib.h.

Referenced by fillHist(), initialize(), and XtInteMdcCalib().

double XtInteMdcCalib::m_tbinWid[MdcCalNLayer][3] [private]

Definition at line 40 of file XtInteMdcCalib.h.

Referenced by initialize(), and XtInteMdcCalib().

const int XtInteMdcCalib::NENTR = 18 [static, private]

Definition at line 26 of file XtInteMdcCalib.h.

Referenced by clear(), and initialize().

const int XtInteMdcCalib::NTBIN_INN = 72 [static, private]

Definition at line 27 of file XtInteMdcCalib.h.

const int XtInteMdcCalib::NTBIN_OUT = 79 [static, private]

Definition at line 28 of file XtInteMdcCalib.h.


Generated on Tue Nov 29 23:36:22 2016 for BOSS_7.0.2 by  doxygen 1.4.7