/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 | Private Member Functions | Private Attributes
CenterOfChargePosTool Class Reference

#include <CenterOfChargePosTool.h>

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

List of all members.

Public Member Functions

 CenterOfChargePosTool (const std::string &type, const std::string &name, const IInterface *parent)
virtual ~CenterOfChargePosTool ()
virtual StatusCode reconstruct (const DayaBay::CalibReadout &readout, DayaBay::RecTrigger &recTrigger)
virtual StatusCode initialize ()
virtual StatusCode finalize ()

Static Public Member Functions

static const InterfaceID & interfaceID ()

Private Member Functions

void GetVertexShift (double r_rec, double z_rec, double &dr_int, double &dz_int)

Private Attributes

std::string m_correctionFileName
std::string m_cableSvcName
ICableSvcm_cableSvc
std::string m_pmtGeomSvcName
IPmtGeomInfoSvcm_pmtGeomSvc
bool m_applyScaling
bool m_applyNewScaling
bool m_output_correction_table
double r [rbins+1]
double z [zbins]
double dz [rbins+1][zbins]
double dr [rbins+1][zbins]

Detailed Description

Definition at line 30 of file CenterOfChargePosTool.h.


Constructor & Destructor Documentation

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

Definition at line 25 of file CenterOfChargePosTool.cc.

  : GaudiTool(type,name,parent)
  , m_cableSvc(0)
  , m_pmtGeomSvc(0)
{
    declareInterface< IReconTool >(this) ;   
    declareProperty("CableSvcName",m_cableSvcName="CableSvc",
                    "Name of service to map between detector, hardware, and electronic IDs");
    declareProperty("PmtGeomSvcName", m_pmtGeomSvcName="PmtGeomInfoSvc",
                    "Name of PMT geometry data service");
    declareProperty("ApplyScaling", m_applyScaling=false,
           "Apply empirical scaling factor to get a 'reconstructed' position");
    declareProperty("ApplyNewScaling", m_applyNewScaling=true,
           "Apply new empirical scaling factor to get a 'reconstructed' position");

    std::string Path = getenv ("DATASVCROOT");
    declareProperty("CorrectionFileName", m_correctionFileName = Path+"/share/VertexCorrection.txt",
                    "File name of the vertex correction table");

    declareProperty("OutputCorrectionTable",m_output_correction_table=false,
                    "Output the correction table for debugging.");
    
}
CenterOfChargePosTool::~CenterOfChargePosTool ( ) [virtual]

Definition at line 51 of file CenterOfChargePosTool.cc.

{}

Member Function Documentation

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

Implements IReconTool.

Definition at line 110 of file CenterOfChargePosTool.cc.

{
  if( !readout.detector().isAD() ){
    debug() << "Not an AD readout; ignoring detector " 
            << readout.detector().detName() << endreq;
    recTrigger.setPositionStatus( ReconStatus::kNotProcessed );
    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.setPositionStatus( ReconStatus::kBadReadout );
    return StatusCode::FAILURE;
  }

  // Context for this data
  int task = 0;
  ServiceMode svcMode(readout.header()->context(), task);

  // Loop over hits and determine charge-weighted mean position
  double pos[3];
  double posSquare[3];
  for(int idx = 0; idx<3; idx++){ pos[idx]=0; posSquare[idx]=0; }
  double time = 0;
  double qSum = 0;
  int nHits = 0;
  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
  
    const CLHEP::Hep3Vector& pmtPos 
      = m_pmtGeomSvc->get(pmtId.fullPackedData())->localPosition();
    double peakAdc=channel.earliestCharge(-1650,-1250);
    double firstTdc = channel.earliestTime();
    for(int idx = 0; idx<3; idx++){
      pos[idx] += peakAdc*pmtPos[idx];
      posSquare[idx] += peakAdc*pmtPos[idx]*pmtPos[idx];
    }
    if (peakAdc > 0){
      time += firstTdc * peakAdc;
      qSum += peakAdc;
      nHits++;
    }
  }

  if(nHits<1){
    CLHEP::HepLorentzVector position(1e10, 1e10, 1e10, 1e10); // set dummy number in case someone try to use it
    recTrigger.setPosition( position );
    recTrigger.setPositionQuality( 1e10 );
    recTrigger.setPositionStatus( ReconStatus::kNoHits );
    return StatusCode::SUCCESS;
  }
  double posStdDev = 0;
  for(int idx = 0; idx<3; idx++){
    posStdDev += posSquare[idx] - pos[idx]*pos[idx]/qSum;
    pos[idx] /= qSum;
  }
  posStdDev /= qSum;
  posStdDev = sqrt(posStdDev);
  time /= qSum;

  if(m_applyScaling){
    // Apply empirical rho, z scaling to vertex position
    // This is not very good, but is useful as a prefit
    double rho = sqrt(pos[0]*pos[0] + pos[1]*pos[1]);
    double rhoNorm = rho/2340.;
    double zNorm = pos[2]/1750.;
    if(zNorm<0) zNorm *= -1;
    double rhoScaled = 1-(1-rhoNorm)/TMath::Power((1+rhoNorm), 3.13);
    double zScaled = 
      1-(1-zNorm)/TMath::Power((1+zNorm), 
                               10.0*(1-sqrt(rhoNorm)));
    //warning() << "Scaling by: " << rhoScaled << " , " << zScaled << endreq;
    pos[0] *= rhoScaled/rhoNorm;
    pos[1] *= rhoScaled/rhoNorm;
    pos[2] *= zScaled/zNorm;
  }

  if(m_applyNewScaling){
    double delta_r = 0;
    double delta_z = 0;
    double rho = sqrt(pos[0]*pos[0] + pos[1]*pos[1]);
    GetVertexShift(rho,pos[2],delta_r,delta_z);
    debug() << "(r,z,dr,dz) =  " << rho << " " << pos[2]
              << " " << delta_r << " " << delta_z << endreq;

    pos[0] += delta_r * pos[0]/rho;
    pos[1] += delta_r * pos[1]/rho;
    pos[2] += delta_z;
    
  }

  CLHEP::HepLorentzVector position(pos[0], pos[1], pos[2], time);
  recTrigger.setPosition( position );
  // Set quality to charge-weighted standard deviation of hit positions
  recTrigger.setPositionQuality( posStdDev );
  recTrigger.setPositionStatus( ReconStatus::kGood );

  return StatusCode::SUCCESS;
}
StatusCode CenterOfChargePosTool::initialize ( ) [virtual]

Definition at line 53 of file CenterOfChargePosTool.cc.

{
  // Get Cable Service
  m_cableSvc = svc<ICableSvc>(m_cableSvcName,true);
  m_pmtGeomSvc = svc<IPmtGeomInfoSvc>(m_pmtGeomSvcName,true);


  // July 1, 2011: Changed to set dr = 0 at r = 0
  if (m_applyNewScaling){
    info() << "Reading vertex correction teble from " << m_correctionFileName.c_str() << endreq;
    std::ifstream drdz_table(m_correctionFileName.c_str());
    //        std::ifstream drdz_table("/u/ynakajim/mywork/vertex_reconstruction/coc_correction_wgt.txt");
    for (int ir = 0; ir < rbins; ir++){
      for (int iz = 0; iz < zbins; iz++){
        drdz_table >> r[ir+1] >> z[iz]
                   >> dr[ir+1][iz]
                   >> dz[ir+1][iz];
        debug() << r[ir+1] << " "
                  << z[iz] << " "
                  << dr[ir+1][iz] << " "
                  << dz[ir+1][iz] << endreq;
      }
    }
    r[0] = 0;
    for (int iz = 0; iz < zbins; iz++){
      dr[0][iz] = 0;
      dz[0][iz] = dz[1][iz];
    }
    if (m_output_correction_table){
      TFile * f_out = new TFile ("vertex_corrections.root","recreate");
      TH2D * h_dr =new TH2D("h_dr","h_dr",1000,0,2000,1000,-2000,2000);
      TH2D * h_dz =new TH2D("h_dz","h_dz",1000,0,2000,1000,-2000,2000);
      for (Int_t i = 0; i < 1000; i++){
        Double_t r = h_dr->GetXaxis()->GetBinCenter(i+1);
        for (Int_t j = 0; j < 1000; j++){
          Double_t z = h_dr->GetYaxis()->GetBinCenter(j+1);
          Double_t dr_int;
          Double_t dz_int;
          GetVertexShift(r,z,dr_int,dz_int);
          h_dr->SetBinContent(i+1,j+1,dr_int);
          h_dz->SetBinContent(i+1,j+1,dz_int);
        }
      }
      f_out->Write();
      f_out->Close();
    }
      
  }
  
  return StatusCode::SUCCESS;
}
StatusCode CenterOfChargePosTool::finalize ( ) [virtual]

Definition at line 105 of file CenterOfChargePosTool.cc.

{
  return StatusCode::SUCCESS;
}
void CenterOfChargePosTool::GetVertexShift ( double  r_rec,
double  z_rec,
double &  dr_int,
double &  dz_int 
) [private]

Definition at line 220 of file CenterOfChargePosTool.cc.

                                                                                                    {

  double r_binwidth = 2000./(double)rbins;
  double z_binwidth = 4000./(double)zbins;
  int rid = (int)(r_rec/r_binwidth+0.5); // it is OK, since rid = 0 is used for r = 0
  int zid = (int)((z_rec+2000)/z_binwidth-0.5);
  if (rid < 0) {
    rid = 0;
    r_rec = r[0];
  }
  if (zid < 0) {
    zid = 0;
    z_rec = z[0];
  }
  if (rid > rbins -1){
    rid = rbins-1;
    r_rec = r[rbins];
  }
  if (zid > zbins -2){
    zid = zbins-2;
    z_rec = z[zbins-1];
  }
  

  dr_int
    = dr[rid][zid]*(r[rid+1]-r_rec)*(z[zid+1]-z_rec)
    + dr[rid+1][zid]*(r_rec-r[rid])*(z[zid+1]-z_rec)
    + dr[rid][zid+1]*(r[rid+1]-r_rec)*(z_rec-z[zid])
    + dr[rid+1][zid+1]*(r_rec-r[rid])*(z_rec-z[zid]);
  dr_int = dr_int/(r[rid+1]-r[rid])/(z[zid+1]-z[zid]);

  dz_int
    = dz[rid][zid]*(r[rid+1]-r_rec)*(z[zid+1]-z_rec)
    + dz[rid+1][zid]*(r_rec-r[rid])*(z[zid+1]-z_rec)
    + dz[rid][zid+1]*(r[rid+1]-r_rec)*(z_rec-z[zid])
    + dz[rid+1][zid+1]*(r_rec-r[rid])*(z_rec-z[zid]);
  dz_int = dz_int/(r[rid+1]-r[rid])/(z[zid+1]-z[zid]);

}

Member Data Documentation

Definition at line 49 of file CenterOfChargePosTool.h.

Definition at line 52 of file CenterOfChargePosTool.h.

Definition at line 54 of file CenterOfChargePosTool.h.

Definition at line 57 of file CenterOfChargePosTool.h.

Definition at line 59 of file CenterOfChargePosTool.h.

Definition at line 62 of file CenterOfChargePosTool.h.

Definition at line 64 of file CenterOfChargePosTool.h.

Definition at line 66 of file CenterOfChargePosTool.h.

double CenterOfChargePosTool::r[rbins+1] [private]

Definition at line 70 of file CenterOfChargePosTool.h.

double CenterOfChargePosTool::z[zbins] [private]

Definition at line 71 of file CenterOfChargePosTool.h.

double CenterOfChargePosTool::dz[rbins+1][zbins] [private]

Definition at line 72 of file CenterOfChargePosTool.h.

double CenterOfChargePosTool::dr[rbins+1][zbins] [private]

Definition at line 73 of file CenterOfChargePosTool.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:13:03 for CenterOfChargePos by doxygen 1.7.4