/search.css" rel="stylesheet" type="text/css"/> /search.js">
| Classes | Job Modules | Data Objects | Services | Algorithms | Tools | Packages | Directories | Tracs |

In This Package:

Public Member Functions | Static Public Member Functions | Public Attributes | Private Attributes
MuonCombinedRecTool Class Reference

#include <MuonCombinedRecTool.h>

Inheritance diagram for MuonCombinedRecTool:
Inheritance graph
[legend]
Collaboration diagram for MuonCombinedRecTool:
Collaboration graph
[legend]

List of all members.

Public Member Functions

 MuonCombinedRecTool (const std::string &type, const std::string &name, const IInterface *parent)
virtual ~MuonCombinedRecTool ()
virtual StatusCode reconstruct (const DayaBay::CalibReadout &readout, DayaBay::RecTrigger &recTrigger)
virtual StatusCode initialize ()
virtual StatusCode finalize ()
void Rec (double *, double *, int &, double &, double &, double &, double &, double &, double &, double *, double *, double *)
void CombineRPC (int &, double &, double &, double &, double &, double &, double &, double, double, double)
StatusCode RpcPointToAdLocal (const DayaBay::CalibReadout &readout, Gaudi::XYZPoint &rpcPoint)

Static Public Member Functions

static const InterfaceID & interfaceID ()

Public Attributes

std::map< Site::Site_t,
Gaudi::XYZPoint > 
m_rpcRecOffset

Private Attributes

std::string m_pmtGeomSvcName
IPmtGeomInfoSvcm_pmtGeomSvc
int m_combineRPC
double m_rpcTime
std::vector< const
DayaBay::RecRpcCluster * > 
m_rpcClusters
std::string m_detDataSvcName
IDataProviderSvc * m_detDataSvc
std::string m_rpcRecHeaderLocation

Detailed Description

Definition at line 38 of file MuonCombinedRecTool.h.


Constructor & Destructor Documentation

MuonCombinedRecTool::MuonCombinedRecTool ( const std::string &  type,
const std::string &  name,
const IInterface *  parent 
)

Definition at line 32 of file MuonCombinedRecTool.cc.

  : GaudiTool(type,name,parent)
  , m_pmtGeomSvc(0)
{
    declareInterface< IReconTool >(this) ;  
    declareProperty("CombineRPC", m_combineRPC = 0, "RPC use or not"); 
    declareProperty("PmtGeomSvcName", m_pmtGeomSvcName="PmtGeomInfoSvc",
                    "Name of PMT geometry data service");
    declareProperty("RpcRecHeaderLocation", m_rpcRecHeaderLocation = "/Event/Rec/RpcSimple",
                    "RPC RecHeader location");
    declareProperty("DetDataSvcName", m_detDataSvcName = "DetectorDataSvc",
                    "Name of detector data service");
}
MuonCombinedRecTool::~MuonCombinedRecTool ( ) [virtual]

Definition at line 48 of file MuonCombinedRecTool.cc.

{}

Member Function Documentation

StatusCode MuonCombinedRecTool::reconstruct ( const DayaBay::CalibReadout readout,
DayaBay::RecTrigger recTrigger 
) [virtual]

Implements IReconTool.

Definition at line 66 of file MuonCombinedRecTool.cc.

{
  if( readout.detector().isRPC()){
        m_rpcTime = readout.triggerTime().GetSeconds();
    
        if (!exist<DayaBay::RecRpcHeader>(evtSvc(), m_rpcRecHeaderLocation)) {
            warning() << "Cannot find header at " << m_rpcRecHeaderLocation << endreq;
            m_rpcClusters.clear();
        } else {
            const DayaBay::RecRpcHeader *rpcRecHdr = get<DayaBay::RecRpcHeader>(m_rpcRecHeaderLocation);
            m_rpcClusters = rpcRecHdr->recTrigger().clusters();
        }
  }
  if( !readout.detector().isAD() ){
    debug() << "Not an AD readout; ignoring detector " 
            << readout.detector().detName() << endreq;
    recTrigger.setDirectionStatus( -1 );
    return StatusCode::SUCCESS;
  }
  const CalibReadoutPmtCrate* pmtReadout 
    = dynamic_cast<const CalibReadoutPmtCrate*>(&readout);
  if(!pmtReadout){
    error() << "Incorrect type of readout crate for detector "
            << readout.detector().detName() << endreq;
    recTrigger.setDirectionStatus( -1 );
    return StatusCode::SUCCESS;
  }

  if(recTrigger.rawEvis()<400){
        // Hard to reconstruct since track length is not long enough
        recTrigger.setDirectionStatus( -1 ); 
        return StatusCode::SUCCESS;
  }
  double ADC[192]={0}, ADT[192]={0};
  double PMTx[192]={0}, PMTy[192] = {0}, PMTz[192] ={0};
  // Loop over hits and add up charge
  CalibReadoutPmtCrate::PmtChannelReadouts::const_iterator chanIter, 
    chanEnd = pmtReadout->channelReadout().end();
  for(chanIter=pmtReadout->channelReadout().begin(); chanIter != chanEnd; 
      chanIter++){
    const CalibReadoutPmtChannel& channel = *chanIter;
    AdPmtSensor pmtId(channel.pmtSensorId().fullPackedData());
    if( !pmtId.is8inch() ){ continue; }  // Ignore non 8-inch pmt channels
    int count = pmtId.ring()*24-24+pmtId.column()-1;
    double earliestC = channel.earliestCharge(-1650,-1250); 
    double earliestT = channel.earliestTime(-1650,-1250);
    ADT[count] = earliestT;
    ADC[count] = earliestC;
    if(earliestC<10){
        ADT[count]=-1;
        ADC[count] = 0;
    }
    const CLHEP::Hep3Vector& pmtPos = m_pmtGeomSvc->get(pmtId.fullPackedData())->localPosition();
    PMTx[count] = pmtPos[0]/1000.;
    PMTy[count] = pmtPos[1]/1000.;
    PMTz[count] = pmtPos[2]/1000.; 
  }
  double InX=0, InY=0, InZ=0, OutX=0, OutY=0,OutZ=0;
  int MuonClass=0;
  Rec(ADT,ADC,MuonClass,InX,InY,InZ,OutX,OutY,OutZ,PMTx,PMTy,PMTz);
  double triggerTimeAfterRpc = readout.triggerTime().GetSeconds() - m_rpcTime;
  if(m_combineRPC){
        if(triggerTimeAfterRpc<1e-6&&m_rpcClusters.size()>0){
                const CLHEP::HepLorentzVector &rpcVector = m_rpcClusters[0]->position();
                Gaudi::XYZPoint rpcPoint(rpcVector.x(), rpcVector.y(), rpcVector.z());
                rpcPoint = rpcPoint - m_rpcRecOffset[readout.detector().site()];
                double RPCX=0,RPCY=0,RPCZ=0;
                if (RpcPointToAdLocal(readout, rpcPoint) != StatusCode::SUCCESS) {
                        warning() << "Cannot transform RPC cluster position to AD local coordinates" << endreq;
                        RPCX=-1,RPCY=-1,RPCZ=-1;
                }
                else {
                        RPCX = rpcPoint.X();
                        RPCY = rpcPoint.Y();
                        RPCZ = rpcPoint.Z();
                }
                CombineRPC(MuonClass,InX,InY,InZ,OutX,OutY,OutZ,RPCX/1000,RPCY/1000,RPCZ/1000);
        }
        else{ // no RPC  trigger before
                if(MuonClass==0) MuonClass=5;
                else MuonClass=6;
        }
  }
  CLHEP::HepLorentzVector position(InX*1000, InY*1000, InZ*1000);
  CLHEP::HepLorentzVector direction(OutX*1000,OutY*1000,OutZ*1000);
  recTrigger.setPosition(position);
  recTrigger.setDirection(direction);
  recTrigger.setDirectionStatus(MuonClass);
  return StatusCode::SUCCESS;
}
StatusCode MuonCombinedRecTool::initialize ( ) [virtual]

Definition at line 50 of file MuonCombinedRecTool.cc.

{
  // Get Cable Service
  m_pmtGeomSvc = svc<IPmtGeomInfoSvc>(m_pmtGeomSvcName,true);
  m_rpcRecOffset[Site::kDayaBay] = Gaudi::XYZPoint(2500.0, -500.0, 12500.0);
  m_rpcRecOffset[Site::kLingAo] = Gaudi::XYZPoint(2500.0, -500.0, 12500.0);
  m_rpcRecOffset[Site::kFar] = Gaudi::XYZPoint(5650.0, 500.0, 12500.0);
  m_detDataSvc = svc<IDataProviderSvc>(m_detDataSvcName, true);
  return StatusCode::SUCCESS;
}
StatusCode MuonCombinedRecTool::finalize ( ) [virtual]

Definition at line 61 of file MuonCombinedRecTool.cc.

{
  return StatusCode::SUCCESS;
}
void MuonCombinedRecTool::Rec ( double *  ,
double *  ,
int &  ,
double &  ,
double &  ,
double &  ,
double &  ,
double &  ,
double &  ,
double *  ,
double *  ,
double *   
)

Definition at line 158 of file MuonCombinedRecTool.cc.

                                                                                                                                                                                                                   {
        double max = -1000;
        int earliestPMT=0;
        // Module 1: Remove strange PMT from the event, find the earliest one
        double temp[192]={0};
        for(int i=0;i<192;i++){
                int test1 = i-24, test2=i+24;
                if(test1>=0)
                        if(fabs(ADT[test1]-ADT[i])>10) temp[i]+=1;
                if(test2<192)
                        if(fabs(ADT[test2]-ADT[i])>10) temp[i]+=1;
        }
        for(int i=0;i<192;i++){
                if(temp[i]==1&&(i<24||i>=7*24)) ADT[i] = -1;
                if(temp[i]==2) ADT[i]=-1;
        }
        for(int i=0;i<192;i++){
                if(ADT[i]<max){ max=ADT[i];earliestPMT=i;}
        }
        // Module 1: Finish
        // Module 3: Test the back PMTs on the top ring to the earliest PMT, give a guess to MuonClass
        int earliestColumn = earliestPMT%24;
        int backColumn = earliestColumn + 12;
        if(backColumn>23) backColumn-=24;
        double MeanBack =0, BackCount=0;
        for(int i=0;i<3;i++){
                int currentColumn = backColumn-1+i;
                if(currentColumn<0) currentColumn+=24;
                if(currentColumn>23) currentColumn-=24;
                if(ADT[7*24+currentColumn]!=-1) {
                        MeanBack+=(ADT[7*24+currentColumn]-max);
                        BackCount+=1;
                }
        }
        MeanBack/=BackCount;
        if(MeanBack<20.76-1.5)
                MuonClass=1;
        else
                MuonClass=2;
        // Module 3: Finish
        double maxCharge=0;int  maxChargePMT=0;
        // Module 5: Fine study to MuonClass 1
        if(MuonClass==1){
                // Guess In
                GuessIn(ADT, InX, InY, InZ,max);
                // Guess Out
                GuessOut(ADT, ADC, maxChargePMT,maxCharge,sqrt(InX*InX+InY*InY));
                // Combine charge and timei, if consistent, MuonClass = 0
                int MaxDiffPMT = Combine(InX, InY, InZ, ADT,max,PMTx,PMTy,PMTz);
                if(Near(MaxDiffPMT,maxChargePMT )){
                        MuonClass=0;
                        double OutR = 2;
                        OutZ = PMTz[MaxDiffPMT];
                        OutX = PMTx[MaxDiffPMT]*2/2.3245;
                        OutY = PMTy[MaxDiffPMT]*2/2.3245;
                }
                else{
                        OutX=0;
                        OutY=0;
                        OutZ=0;
                }
        }
        // Module 5: Finish 

        // Module 6: Fine study to MuonClass 2
        if(MuonClass==2){
                double InjectR = 2;
                InZ = PMTz[earliestPMT];
                InX = PMTx[earliestPMT]*2/2.3245;
                InY = PMTy[earliestPMT]*2/2.3245;
                double meanTop=0, meanTopCount=0;
                for(int i=0;i<3;i++){
                        int currentColumn = earliestColumn-1+i;
                        if(currentColumn<0) currentColumn+=24;
                        if(currentColumn>23) currentColumn-=24;
                        if(ADT[7*24+currentColumn]!=-1) {
                                meanTop+=(ADT[7*24+currentColumn]-max);
                                meanTopCount+=1;
                        }
                }
                meanTop/=meanTopCount;
                meanTop/=4.83;
                double guessZ = 1.75-meanTop;
                if(fabs(InZ-guessZ)<0.5) InZ = (InZ+guessZ)/2;

                int MaxDiffPMT = Combine(InX, InY, InZ, ADT,max,PMTx,PMTy,PMTz);
                GuessOut2(ADT, ADC, maxChargePMT,maxCharge,earliestPMT);
                if(Near(MaxDiffPMT,maxChargePMT )){
                        MuonClass=0;
                        OutZ = PMTz[MaxDiffPMT];
                        OutX = PMTx[MaxDiffPMT]*2/2.3245;
                        OutY = PMTy[MaxDiffPMT]*2/2.3245;
                        if(OutX==InX) {
                                MuonClass=1;
                                OutZ=0;
                                OutX=0;
                                OutY=0;
                        }
                }
                else {
                        MuonClass=1;
                        OutZ=0;
                        OutX=0;
                        OutY=0;
                }
        }
        
}
void MuonCombinedRecTool::CombineRPC ( int &  MuonClass,
double &  InX,
double &  InY,
double &  InZ,
double &  OX,
double &  OY,
double &  OZ,
double  RPCX,
double  RPCY,
double  RPCZ 
)

Definition at line 349 of file MuonCombinedRecTool.cc.

                                                                                                                                                                   {
        if(RPCX==-0.001) { // RPC no readout, set as 5/6
                if(MuonClass==0) MuonClass=5;
                if(MuonClass==1) MuonClass=6;
                return;
        }
        double theta=0, phi=0;
        if(OX!=0){ // AD has a track
                double Dis = sqrt((InX-OX)*(InX-OX)+(InY-OY)*(InY-OY)+(InZ-OZ)*(InZ-OZ));
                double RecTheta =  acos(fabs(InZ-OZ)/Dis)*180/3.1415;
                double Dis1 = sqrt((InX-RPCX)*(InX-RPCX)+(InY-RPCY)*(InY-RPCY)+(InZ-RPCZ)*(InZ-RPCZ));
                double RecTheta1 = acos(fabs(InZ-RPCZ)/Dis1)*180/3.1415;
                double vec1 = OX-InX, vec2 = OY-InY;
                double vec3 = InX-RPCX, vec4= InY-RPCY;
                theta = (RecTheta+RecTheta1)/2;
                double Phi1 = acos((vec1)/sqrt(vec1*vec1+vec2*vec2))*180/3.1415;
                double Phi2 = acos((vec3)/sqrt(vec3*vec3+vec4*vec4))*180/3.1415;
                if(vec2<0) Phi1=360-Phi1;
                if(vec4<0) Phi2=360-Phi2;
                double ThetaDiff = fabs(RecTheta-RecTheta1);
                double PhiDiff = fabs(Phi1-Phi2); // theta phi differences at AD and RPC
                phi = (Phi1+Phi2)/2;
                if(ThetaDiff<10&&PhiDiff<10){ //good, mean of RPC and AD are used
                        MuonClass=1;
                        CalculateOut(InX, InY, InZ, OX, OY, OZ, theta, phi);
                }
                else if(ThetaDiff<20&&PhiDiff<20){ // good, mean of RPC and AD are used
                        MuonClass=2;
                        CalculateOut(InX, InY, InZ, OX, OY, OZ, theta, phi);
                }
                else { // too bad, only AD is stored
                        MuonClass=3;
                }
        }
        else{ // AD only give inject point, use RPC to construct a track
                MuonClass=4;
                double Dis1 = sqrt((InX-RPCX)*(InX-RPCX)+(InY-RPCY)*(InY-RPCY)+(InZ-RPCZ)*(InZ-RPCZ));
                theta = acos(fabs(InZ-RPCZ)/Dis1)*180/3.1415;
                double vec3 = InX-RPCX, vec4= InY-RPCY;
                phi = acos((vec3)/sqrt(vec3*vec3+vec4*vec4))*180/3.1415;
                if(vec4<0) phi=360-phi;
                int tag = CalculateOut(InX, InY, InZ, OX, OY, OZ, theta, phi);
                // tag=0: track length is too small, indicating some bad tracks
                if(tag==0) MuonClass=6;
        }
}
StatusCode MuonCombinedRecTool::RpcPointToAdLocal ( const DayaBay::CalibReadout readout,
Gaudi::XYZPoint &  rpcPoint 
)

Definition at line 397 of file MuonCombinedRecTool.cc.

{
    char *siteName = "";
    if (readout.detector().site() == Site::kDayaBay)
        siteName = "db";
    else
    if (readout.detector().site() == Site::kLingAo)
        siteName = "la";
    else
    if (readout.detector().site() == Site::kFar)
        siteName = "far";

    const char *rpcDEName = Form("/dd/Structure/DayaBay/%s-rock/%s-rpc", siteName, siteName);
    SmartDataPtr<IDetectorElement> rpcDE(m_detDataSvc, rpcDEName);
    if (!rpcDE)
        return StatusCode::FAILURE;
    Gaudi::XYZPoint rpcGlobalPoint = rpcDE->geometry()->toGlobal(rpcPoint);

    const char *adDEName = Form("/dd/Structure/DayaBay/%s-rock/%s-ows/%s-curtain/%s-iws/%s-ade%d/%s-sst%d/%s-oil%d",
                                  siteName, siteName, siteName, siteName,
                                  siteName, readout.detector().detectorId(),
                                  siteName, readout.detector().detectorId(),
                                  siteName, readout.detector().detectorId());
    SmartDataPtr<IDetectorElement> adDE(m_detDataSvc, adDEName);
    if (!adDE)
        return StatusCode::FAILURE;

    rpcPoint = adDE->geometry()->toLocal(rpcGlobalPoint);

    return StatusCode::SUCCESS;
}

Member Data Documentation

Definition at line 56 of file MuonCombinedRecTool.h.

Definition at line 59 of file MuonCombinedRecTool.h.

Definition at line 61 of file MuonCombinedRecTool.h.

Definition at line 63 of file MuonCombinedRecTool.h.

Definition at line 65 of file MuonCombinedRecTool.h.

Definition at line 66 of file MuonCombinedRecTool.h.

Definition at line 68 of file MuonCombinedRecTool.h.

IDataProviderSvc* MuonCombinedRecTool::m_detDataSvc [private]

Definition at line 70 of file MuonCombinedRecTool.h.

Definition at line 71 of file MuonCombinedRecTool.h.


The documentation for this class was generated from the following files:
| Classes | Job Modules | Data Objects | Services | Algorithms | Tools | Packages | Directories | Tracs |

Generated on Fri May 16 2014 10:11:52 for MuonCombinedRec by doxygen 1.7.4