/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
ExpQCalcTool Class Reference

#include <ExpQCalcTool.h>

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

List of all members.

Public Member Functions

 ExpQCalcTool (const std::string &type, const std::string &name, const IInterface *parent)
virtual ~ExpQCalcTool ()
virtual vector< float > expqcalc (Site::Site_t site, DetectorId::DetectorId_t detid, CLHEP::Hep3Vector vertex)
virtual vector< float > geomcalc (Site::Site_t site, DetectorId::DetectorId_t detid, CLHEP::Hep3Vector vertex, string name)
virtual StatusCode initialize ()
virtual StatusCode finalize ()
StatusCode getOpPara (std::string &data_src)
StatusCode getGeomPara (std::string &data_src)
void printPara ()

Static Public Member Functions

static const InterfaceID & interfaceID ()

Public Attributes

std::string m_cableSvcName
ICableSvcm_cableSvc
std::string m_pmtGeomSvcName
IPmtGeomInfoSvcm_pmtGeomSvc
OpParam_opPara
std::string m_opLocation
GeomParam_geomPara
std::string m_geomLocation
double m_ATTL
double m_TopRef
double m_BotRef
double m_LightYield
double m_PhotoCathodeArea
double m_TopRefZpos
double m_BotRefZpos

Detailed Description

Definition at line 34 of file ExpQCalcTool.h.


Constructor & Destructor Documentation

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

Definition at line 24 of file ExpQCalcTool.cc.

  :  GaudiTool(type, name, parent)
  , m_cableSvc(0)
  , m_pmtGeomSvc(0)
  , m_opPara(0)
  , m_geomPara(0)
{
  declareInterface<IReconHelperTool>(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 Information Service");

  // declare property for those optical parameters
  declareProperty("opLocation", m_opLocation="DetDesc", "Optical parameters location");
  declareProperty("geomLocation", m_geomLocation="DetDesc", "Optical parameters location");
  declareProperty("AbsLength", m_ATTL=0., "Effective absorption length");
  declareProperty("TopRefZ", m_TopRefZpos=0., "Z position of top reflector");
  declareProperty("BotRefZ", m_BotRefZpos=0., "Z position of bottom reflector");
  declareProperty("TopReflectivity", m_TopRef=0., "Reflectivity of top reflector");
  declareProperty("BotReflectivity", m_BotRef=0., "Reflectivity of bottom reflector");
}
ExpQCalcTool::~ExpQCalcTool ( ) [virtual]

Definition at line 50 of file ExpQCalcTool.cc.

{}

Member Function Documentation

vector< float > ExpQCalcTool::expqcalc ( Site::Site_t  site,
DetectorId::DetectorId_t  detid,
CLHEP::Hep3Vector  vertex 
) [virtual]

Implements IReconHelperTool.

Definition at line 76 of file ExpQCalcTool.cc.

                             {

  // Always calculate the value of the function, FVAL,
  // which is usually a chisquare or log likelihood.
    vector<float> qvec;
    qvec.clear();
  
    double S_angle, dist;
    double expected_q, cosa_tmp;
    double AA, BB;
    double ZF, dx, dy, dz;
    double sumQ = 0.0;

    // Loop over calib readout channels
    // Remember: only non-zero readout channels are stored the map
    for(int local_id=0; local_id<192; local_id++) {
      int ring = local_id/24 + 1;
      int col  = local_id%24 + 1;

      AdPmtSensor adPmtId(ring, col, site, detid);
      unsigned int pmtid = adPmtId.fullPackedData();
      CLHEP::Hep3Vector pmtpos =  m_pmtGeomSvc->get(pmtid)->localPosition();
      CLHEP::Hep3Vector pmtnorm = m_pmtGeomSvc->get(pmtid)->localDirection();

      AA = pmtnorm.x();
      BB = pmtnorm.y();

      dx = vertex.x() - pmtpos.x();
      dy = vertex.y() - pmtpos.y();
      dz = vertex.z() - pmtpos.z();

      dist = dx*dx + dy*dy + dz*dz;
      if(dist<0.0001) dist = 0.0001;
      dist = sqrt(dist);

      cosa_tmp = ( AA*dx + BB*dy ) / dist;
      if(cosa_tmp>1.0) cosa_tmp = 1.0;

      S_angle = 0.378 + 0.5099*cosa_tmp + 0.1131*pow(cosa_tmp,2);
      if(S_angle<1.0e-5) S_angle = 1.0e-5;

      expected_q = S_angle*exp(-dist/m_opPara->m_attl)/TMath::Power(dist, 2);


      int ref_deg;
      for(ref_deg=0; ref_deg<10; ref_deg++){
     
        if(ref_deg==0)
          ZF =  2*m_geomPara->m_botRefZ - vertex.z();   // 1st-order reflection
        else if(ref_deg==1)
          ZF =  2*m_geomPara->m_topRefZ - vertex.z();   // 1st-order reflection
        else if(ref_deg==2)
          ZF = 2*m_geomPara->m_botRefZ - 2*m_geomPara->m_topRefZ + vertex.z(); // 2nd-order reflection
        else if(ref_deg==3)
          ZF = 2*m_geomPara->m_topRefZ - 2*m_geomPara->m_botRefZ + vertex.z(); // 2nd-order reflection
        else if(ref_deg==4)
          ZF = 4*m_geomPara->m_botRefZ - 2*m_geomPara->m_topRefZ - vertex.z(); // 3rd-order reflection
        else if(ref_deg==5)
          ZF = 4*m_geomPara->m_topRefZ - 2*m_geomPara->m_botRefZ - vertex.z(); // 3rd-order reflection
        else if(ref_deg==6)
          ZF = 4*m_geomPara->m_botRefZ - 4*m_geomPara->m_topRefZ + vertex.z(); // 4th-order reflection
        else if(ref_deg==7)
          ZF = 4*m_geomPara->m_topRefZ - 4*m_geomPara->m_botRefZ + vertex.z(); // 4th-order reflection
        else if(ref_deg==8)
          ZF = 6*m_geomPara->m_botRefZ - 4*m_geomPara->m_topRefZ - vertex.z(); // 5th-order reflection
        else
          ZF = 6*m_geomPara->m_topRefZ - 4*m_geomPara->m_botRefZ - vertex.z(); // 5th-order reflection

        dz = ZF - pmtpos.z();

        dist = dx*dx + dy*dy + dz*dz;
        if(dist<0.0001) dist = 0.0001;
        dist = sqrt(dist);

        cosa_tmp = ( AA*dx + BB*dy ) / dist;
        if(cosa_tmp>1.0) cosa_tmp = 1.0;

        S_angle = 0.378 + 0.5099*cosa_tmp + 0.1131*pow(cosa_tmp,2);
        if(S_angle<1.0e-5) S_angle = 1.0e-5;

        int index = (int)(ref_deg/2.0) + 1;
        int index_top, index_bot;
        if(ref_deg == 0) 
          index_top = 0;
        else 
          index_top = (int)((ref_deg-1)/4.0) + 1;

        index_bot = index - index_top;

        expected_q += TMath::Power(m_opPara->m_topRef, index_top)*TMath::Power(m_opPara->m_botRef, index_bot)*S_angle*TMath::Exp(-dist/m_opPara->m_attl)/TMath::Power(dist, 2);
      }

      sumQ += expected_q;
      qvec.push_back(expected_q);
    }

    for(int i=0; i<192; i++) {
      qvec[i] = qvec[i]/sumQ;
    }

    return qvec;
}
vector< float > ExpQCalcTool::geomcalc ( Site::Site_t  site,
DetectorId::DetectorId_t  detid,
CLHEP::Hep3Vector  vertex,
string  name 
) [virtual]

Implements IReconHelperTool.

Definition at line 183 of file ExpQCalcTool.cc.

                                          {

  // Always calculate the value of the function, FVAL,
  // which is usually a chisquare or log likelihood.
    vector<float> distvec;
    vector<float> anglevec;
    distvec.clear();
    anglevec.clear();
  
    double dist, cosa_tmp;
    double AA, BB, dx, dy, dz;

    // Loop over calib readout channels
    // Remember: only non-zero readout channels are stored the map
    for(int local_id=0; local_id<192; local_id++) {
      int ring = local_id/24 + 1;
      int col  = local_id%24 + 1;

      AdPmtSensor adPmtId(ring, col, site, detid);
      unsigned int pmtid = adPmtId.fullPackedData();
      CLHEP::Hep3Vector pmtpos =  m_pmtGeomSvc->get(pmtid)->localPosition();
      CLHEP::Hep3Vector pmtnorm = m_pmtGeomSvc->get(pmtid)->localDirection();

      AA = pmtnorm.x();
      BB = pmtnorm.y();

      dx = vertex.x() - pmtpos.x();
      dy = vertex.y() - pmtpos.y();
      dz = vertex.z() - pmtpos.z();

      dist = dx*dx + dy*dy + dz*dz;
      if(dist<0.0001) dist = 0.0001;
      dist = sqrt(dist);

      cosa_tmp = ( AA*dx + BB*dy ) / dist;
      if(cosa_tmp>1.0) cosa_tmp = 1.0;

      distvec.push_back(dist);
      anglevec.push_back(cosa_tmp);
    }

    if(name == "dist") 
      return distvec;
    else if(name == "angle") 
      return anglevec;
    else {
      error() << "Unknown geomcalc Type!" << endl;
      exit(0);
    }
}
StatusCode ExpQCalcTool::initialize ( ) [virtual]

Definition at line 52 of file ExpQCalcTool.cc.

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

  // Get Optical Parameters
  getOpPara(m_opLocation);

  // Get Geometry Parameters
  getGeomPara(m_geomLocation);

  // Show Optical and Geom Parameters
  printPara();

  return StatusCode::SUCCESS;
}
StatusCode ExpQCalcTool::finalize ( ) [virtual]

Definition at line 71 of file ExpQCalcTool.cc.

{
    return StatusCode::SUCCESS;
}
StatusCode ExpQCalcTool::getOpPara ( std::string &  data_src)

Definition at line 237 of file ExpQCalcTool.cc.

                                                      {

  // construct OpPara, default value
  m_opPara = new OpPara(1.25E+04, 0.98, 0.98, 9000, 3.14*103*103);

  if(data_src=="DetDesc") {

    // Get Optical Parameters from Det/DetDesc
    // A tabular was defined in DetDesc, but only one value was get here
    // FIXME
    // Get the Detector Data Service to access TDS:
    // This should be done in your tools initialize()
    // and save m_DDS as a data member
    IDataProviderSvc* dds = 0;
    StatusCode sc = service( "DetectorDataSvc", dds, true );
    if (sc.isFailure()) {
        // ... handle error ...
    }
    // Later, look up some material.  If this material is always used you
    // should look it up once in initialize() and save it for later use

    // absorption length
    std::string location = "/dd/Materials/GdDopedLS";
    Material* gdls = GaudiCommon<AlgTool>::get<Material>(dds,location);
    
    Material::Tables& gdls_tab = gdls->tabulatedProperties();
    Material::Tables::const_iterator mtIter = gdls_tab.begin(), 
      done = gdls_tab.end();
    TabulatedProperty::Table attl_vec;
    for(; mtIter!=done; mtIter++) {
      std::string type = (*mtIter)->type();
      if(type=="ABSLENGTH") {
        attl_vec = (*mtIter)->table(); 
      }
    }
    TabulatedProperty::Table::const_iterator tb_iter;
    int count = 0;
    double attl_tmp = 0.0, ref_tmp = 0.0;
    for(tb_iter=attl_vec.begin(); tb_iter!=attl_vec.end(); tb_iter++){
      if(tb_iter->first<3.18e-06 && tb_iter->first>2.07e-06 ){ // 390nm ~ 600nm
        count++; attl_tmp += tb_iter->second;
      }
    }
    attl_tmp = attl_tmp/count; // average value between 390nm~600nm
    m_opPara->m_attl = attl_tmp;
      
    // ESR surface location
    std::string topESR_location = "/dd/Geometry/AdDetails/AdSurfacesAll/ESRAirSurfaceTop";
    std::string botESR_location = "/dd/Geometry/AdDetails/AdSurfacesAll/ESRAirSurfaceBot";

    Surface* esrtop = GaudiCommon<AlgTool>::get<Surface>(dds, topESR_location);
    Surface* esrbot = GaudiCommon<AlgTool>::get<Surface>(dds, botESR_location);
    
    Surface::Tables& esrtop_tab = esrtop->tabulatedProperties();
    Surface::Tables& esrbot_tab = esrbot->tabulatedProperties();

    Surface::Tables::const_iterator sfIter; 
    TabulatedProperty::Table ref_vec;

    // look up Tables, get an average top ESR reflectivity
    for(sfIter=esrtop_tab.begin(); sfIter!=esrtop_tab.end(); sfIter++) {
      std::string type = (*sfIter)->type();
      if(type=="REFLECTIVITY") {
        ref_vec = (*sfIter)->table(); 
      }
    }
    count = 0;
    for(tb_iter=ref_vec.begin(); tb_iter!=ref_vec.end(); tb_iter++){
      if(tb_iter->first<3.18e-06 && tb_iter->first>2.07e-06 ){ // 390nm ~ 600nm
        count++; ref_tmp += tb_iter->second;
      }
    }
    ref_tmp = ref_tmp/count; // average value between 390nm~600nm
    m_opPara->m_topRef = ref_tmp;
 
    // look up Tables, get an average bot ESR reflectivity
    for(sfIter=esrbot_tab.begin(); sfIter!=esrbot_tab.end(); sfIter++) {
      std::string type = (*sfIter)->type();
      if(type=="REFLECTIVITY") {
        ref_vec = (*sfIter)->table(); 
      }
    }
    count = 0; ref_tmp = 0.0;
    for(tb_iter=ref_vec.begin(); tb_iter!=ref_vec.end(); tb_iter++){
      if(tb_iter->first<3.18e-06 && tb_iter->first>2.07e-06 ){ // 390nm ~ 600nm
        count++; ref_tmp += tb_iter->second;
      }
    }
    ref_tmp = ref_tmp/count; // average value between 390nm~600nm
    m_opPara->m_botRef = ref_tmp;

  } else if(data_src=="OpCalibSvc"){

  } else {

    warning() <<"No external optical data source, use default value."<< endreq;
  }

  // if set new values from declareProperty
  if(m_ATTL!=0.) m_opPara->m_attl = m_ATTL;
  if(m_TopRef!=0.) m_opPara->m_topRef = m_TopRef;
  if(m_BotRef!=0.) m_opPara->m_botRef = m_BotRef;

  return StatusCode::SUCCESS;
}
StatusCode ExpQCalcTool::getGeomPara ( std::string &  data_src)

Definition at line 343 of file ExpQCalcTool.cc.

                                                        {

  // construct GeomPara, default value
  m_geomPara = new GeomPara(2010., -2010.);

  if(data_src=="DetDesc") {

    IDataProviderSvc* dds = 0;
    StatusCode sc = service( "DetectorDataSvc", dds, true );
    if (sc.isFailure()) {
       error() << "No DetectorDataSvc available!" << endreq;
    }
    
    // ESR structure location
    std::string topEsrPath = "/dd/Structure/AdReflectorStructure/db-ad/db-ad1-topesr";
    std::string botEsrPath = "/dd/Structure/AdReflectorStructure/db-ad/db-ad1-botesr";
    SmartDataPtr<IDetectorElement> topEsr(dds,topEsrPath);
    SmartDataPtr<IDetectorElement> botEsr(dds,botEsrPath);

    // define a pointer to ICoordSysSvc to
    // use to cache the service during the job.
    ICoordSysSvc* css;
    sc = service( "CoordSysSvc", css, true );

    Gaudi::XYZPoint zero(0,0,0);
    // Set top ESR position
    Gaudi::XYZPoint globalPoint = topEsr->geometry()->toGlobal(zero);
    IDetectorElement* theAdDe = css->coordSysDE(globalPoint);
    Gaudi::XYZPoint adLocalPoint = theAdDe->geometry()->toLocal(globalPoint);
    m_geomPara->m_topRefZ = adLocalPoint.z();

    // Set bot ESR position
    globalPoint = botEsr->geometry()->toGlobal(zero);
    theAdDe = css->coordSysDE(globalPoint);
    adLocalPoint = theAdDe->geometry()->toLocal(globalPoint);
    m_geomPara->m_botRefZ = adLocalPoint.z();
        
  } else if(data_src=="OnsiteSuveyData"){
    // Get geometry info from onsite suvey data

  } else {

    warning() <<"Unknown geometry data source, use default value."<< endreq;
  }

  // if set new values from declareProperty
  if(m_TopRefZpos!=0.) m_geomPara->m_topRefZ = m_TopRefZpos;
  if(m_BotRefZpos!=0.) m_geomPara->m_botRefZ = m_BotRefZpos;

  return StatusCode::SUCCESS;
}
void ExpQCalcTool::printPara ( )

Definition at line 395 of file ExpQCalcTool.cc.

                            {
  info() << "Op and Geom parameters used in reconstruction." << endreq;
  info() << "Attenuation length : " << m_opPara->m_attl <<"mm"<< endreq;
  info() << "Top reflectivity : " << m_opPara->m_topRef << endreq;
  info() << "Bottom reflectivity : " << m_opPara->m_botRef << endreq;
  info() << "Top reflector Z pos :" << m_geomPara->m_topRefZ <<"mm"<< endreq;
  info() << "Bot reflector Z pos :" << m_geomPara->m_botRefZ <<"mm"<< endreq;
}

Member Data Documentation

Definition at line 60 of file ExpQCalcTool.h.

Definition at line 63 of file ExpQCalcTool.h.

Definition at line 66 of file ExpQCalcTool.h.

Definition at line 69 of file ExpQCalcTool.h.

Definition at line 71 of file ExpQCalcTool.h.

Definition at line 72 of file ExpQCalcTool.h.

Definition at line 75 of file ExpQCalcTool.h.

Definition at line 76 of file ExpQCalcTool.h.

Definition at line 78 of file ExpQCalcTool.h.

Definition at line 79 of file ExpQCalcTool.h.

Definition at line 80 of file ExpQCalcTool.h.

Definition at line 81 of file ExpQCalcTool.h.

Definition at line 82 of file ExpQCalcTool.h.

Definition at line 85 of file ExpQCalcTool.h.

Definition at line 86 of file ExpQCalcTool.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:26 for AdRec by doxygen 1.7.4