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

#include <GtDiffuserBallTool.h>

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

List of all members.

Public Member Functions

 GtDiffuserBallTool (const std::string &type, const std::string &name, const IInterface *parent)
virtual ~GtDiffuserBallTool ()
virtual StatusCode initialize ()
virtual StatusCode finalize ()
virtual StatusCode mutate (HepMC::GenEvent &event)
 Modify the event.

Static Public Member Functions

static const InterfaceID & interfaceID ()
 Retrieve interface ID.

Private Member Functions

StatusCode oneVertex (HepMC::GenEvent &event)
double GetWave ()
double GetTime ()
double ConvertCdfRand (const double &rand, const std::vector< double > &cdf, const std::vector< double > &edges)

Private Attributes

int m_particlesPerEvent
std::string m_particlesPerEventMode
int m_particlesPerEventSpread
std::string m_wavelengthMode
std::string m_timingMode
double m_wavelength
double m_wavelengthSpread
double m_pulseFrequency
double m_radius
bool m_anisotropy
double m_bright_theta
double m_bright_phi
double m_photonScaleWeight
std::string m_geomPosToolName
std::string m_particleName
int m_pid
Rndm::Numbers m_uni
Rndm::Numbers m_randWave
Rndm::Numbers m_randTime
Rndm::Numbers m_randPartNo
std::vector< double > m_pdfWave
std::vector< double > m_edgesWave
std::vector< double > m_cdfWave
std::vector< double > m_pdfTime
std::vector< double > m_edgesTime
std::vector< double > m_cdfTime

Detailed Description

Definition at line 49 of file GtDiffuserBallTool.h.


Constructor & Destructor Documentation

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

Definition at line 20 of file GtDiffuserBallTool.cc.

  :GaudiTool(type,name,parent)
{
  // Initialization
  m_pid = 0;
  m_pdfWave.clear();   m_pdfWave.push_back(1);
  m_edgesWave.clear(); m_edgesWave.push_back(429*nm); m_edgesWave.push_back(431*nm);
  m_pdfTime.clear();   m_pdfTime.push_back(1);
  m_edgesTime.clear(); m_edgesTime.push_back(0*ns);   m_edgesTime.push_back(1*ns);

  declareInterface<IHepMCEventMutator>(this);
  declareProperty("PhotonsPerEvent",m_particlesPerEvent=3500,"Number of photons per event");
  declareProperty("PhotonsPerEventMode",m_particlesPerEventMode="Fixed","Type of distribution of number of photons per event");
  declareProperty("PhotonsPerEventSpread",m_particlesPerEventSpread=0,"Sigma/FWHM of number of photons per event");

  declareProperty("WavelengthMode",m_wavelengthMode="Fixed","Wavelength Mode");
  declareProperty("Wavelength", m_wavelength=430.0*nm,"Fixed (Mean) Photon Wavelength");
  declareProperty("WavelengthSpread",m_wavelengthSpread=5.0*nm,"Wavelength Spectrum Spread");
  declareProperty("PdfWave",m_pdfWave,"User-Defined PDF for Wavelength Spectrum");
  declareProperty("PdfEdgesWave",m_edgesWave,"Bin Edges of PdfWave");

  declareProperty("TimingMode",m_timingMode="Instant","Timing Mode");
  declareProperty("PdfTime",m_pdfTime,"User-Defined PDF for Timing Distribution");
  declareProperty("PdfEdgesTime",m_edgesTime,"Bin Edges of PdfTime"); 

  declareProperty("Radius",m_radius=0.9525*cm,"Radius of Diffuser Ball");
  declareProperty("Anisotropy",m_anisotropy=true,"Turn on/off anisotropy");
  declareProperty("BrightSpotTheta",m_bright_theta=0.0,"Polar angle of bright spot");
  declareProperty("BrightSpotPhi",m_bright_phi=0.0,"Azimuthal angle of bright spot");
  declareProperty("PhotonScaleWeight",m_photonScaleWeight = 3.74,
                  "Reweight photons based on maximum quantum efficiency");
  declareProperty("GeomPosTool",m_geomPosToolName="diffuserBallGeomPos",
          "Name of tool used to dynamically place the diffuser ball geometry");

}
GtDiffuserBallTool::~GtDiffuserBallTool ( ) [virtual]

Definition at line 58 of file GtDiffuserBallTool.cc.

{
}

Member Function Documentation

StatusCode GtDiffuserBallTool::initialize ( ) [virtual]

Definition at line 62 of file GtDiffuserBallTool.cc.

{
  // Particle identification
  IParticlePropertySvc * ppSvc = 
    svc< IParticlePropertySvc >( "ParticlePropertySvc" , true ) ;
  ParticleProperty* particle = 0;

  m_particleName = "opticalphoton";
  particle = ppSvc->find(m_particleName);
  if(particle) {
    m_pid = particle->pdgID();
  }

  // if failed to find particle, return failure
  if (!particle) {
    fatal() << "Failed to find particle named \"" 
            << m_particleName << "\" and with ID " 
            << m_pid << endreq;
    return StatusCode::FAILURE;
  }

  // check modes
  if ((m_timingMode != "Instant") && (m_timingMode != "Defined")) {
    fatal() << "Not a recognized vertex timing distribution mode" << endreq;
    return StatusCode::FAILURE;
  }
  if ((m_wavelengthMode != "Fixed")   && (m_wavelengthMode != "Uniform") &&
      (m_wavelengthMode != "Smeared") && (m_wavelengthMode != "Defined")) {
    fatal() << "Not a recognized wavelength spectrum mode" << endreq;
    return StatusCode::FAILURE;
  }

  // Random number service
  IRndmGenSvc *rgs = 0;
  if (service("RndmGenSvc",rgs,true).isFailure()) {
    fatal() << "Failed to get random service" << endreq;
    return StatusCode::FAILURE;        
  }
  
  StatusCode sc;
  sc = m_uni.initialize(rgs, Rndm::Flat(0,1));
  if (sc.isFailure()) {
    fatal() << "Failed to initialize uniform random numbers" << endreq;
    return StatusCode::FAILURE;
  }

  // random numbers for particle number
  if (m_particlesPerEventMode == "Gaus") {
    sc = m_randPartNo.initialize(rgs, Rndm::Gauss(0,1));
    if (sc.isFailure()) {
      fatal() << "Failed to initialize random numbers for particle number distribution" << endreq;
      return StatusCode::FAILURE;
    }
  }
  else if (m_particlesPerEventMode == "Lorentz"){
    sc = m_randPartNo.initialize(rgs, Rndm::BreitWigner(0,1));
    if (sc.isFailure()) {
      fatal() << "Failed to initialize random numbers for particle number distribution" << endreq;
      return StatusCode::FAILURE;
    }
  }
  
  // random numbers for timing
  if (m_timingMode == "Defined") {
    sc = m_randTime.initialize(rgs, Rndm::DefinedPdf(m_pdfTime,0));
    if (sc.isFailure()) {
      fatal() << "Failed to initialize random numbers for timing distribution" << endreq;
      return StatusCode::FAILURE;
    }
  }
  
  // random numbers for wavelength spectrum
  if (m_wavelengthMode == "Uniform") {
    sc = m_randWave.initialize(rgs, Rndm::Flat(0,1));
    if (sc.isFailure()) {
      fatal() << "Failed to initialize random numbers for wavelength spectrum" << endreq;
      return StatusCode::FAILURE;
    }
  }
  else if (m_wavelengthMode == "Smeared") {
    sc = m_randWave.initialize(rgs, Rndm::Gauss(0,1));
    if (sc.isFailure()) {
      fatal() << "Failed to initialize random numbers for wavelength spectrum" << endreq;
      return StatusCode::FAILURE;
    }
  }
  else if (m_wavelengthMode == "Defined") {
    sc = m_randWave.initialize(rgs, Rndm::DefinedPdf(m_pdfWave,0));
    if (sc.isFailure()) {
      fatal() << "Failed to initialize random numbers for wavelength spectrum" << endreq;
      return StatusCode::FAILURE;
    }
  }

  // check user input for user-defined PDFs
  if (m_pdfWave.size()+1 != m_edgesWave.size() || 
      m_pdfTime.size()+1 != m_edgesTime.size() ) {
    fatal() << "There should be 1 more number of bin edges than bins when defining pdf" << endreq;
    return StatusCode::FAILURE; 
  }

  // Fill pre-calculated Cumulative Distribution Functions
  if(m_pdfTime.size()>0){
    m_cdfTime.push_back( 0 );
    for (unsigned int bin=0; bin<m_pdfTime.size(); bin++)
      m_cdfTime.push_back( m_cdfTime[bin] + m_pdfTime[bin] );
    // Normalize
    for (unsigned int bin=0; bin<m_cdfTime.size(); bin++)
      m_cdfTime[bin] /= m_cdfTime.back();
  }
  if(m_pdfWave.size()>0){
    m_cdfWave.push_back( 0 );
    for (unsigned int bin=0; bin<m_pdfWave.size(); bin++)
      m_cdfWave.push_back( m_cdfWave[bin] + m_pdfWave[bin] );
    // Normalize
    for (unsigned int bin=0; bin<m_cdfWave.size(); bin++)
      m_cdfWave[bin] /= m_cdfWave.back();
  }
  
  if( m_photonScaleWeight > 1.0 ){
    // Reduce number of photons, and scale each photon weight to compensate
    info() << "Scaling each photon weight by " << m_photonScaleWeight 
           << " and reducing the number of photons from "
           << m_particlesPerEvent << " to " 
           << int( m_particlesPerEvent / m_photonScaleWeight )
           << endreq;
    m_particlesPerEvent = int( m_particlesPerEvent / m_photonScaleWeight );
  }

  if( m_geomPosToolName == "" ){
    info() << "Placement of Diffuser Ball geometry not requested.  "
           << "Optical photons will still be produced at the requested "
           << "location." << endreq;
  }else{
    IPositionerTool* geomPosTool = 0;
    try{
      geomPosTool = tool<IPositionerTool>(m_geomPosToolName);
    }
    catch(const GaudiException& exg) {
      fatal() << "Failed to get geometry tool: \"" << m_geomPosToolName << "\"" 
              << endreq;
      return StatusCode::FAILURE;
    }
    info () << "Placing diffuser ball geometry using " << m_geomPosToolName 
            << endreq;
    sc = geomPosTool->placeVolume();
    if(sc.isFailure()) return sc;
    geomPosTool->release();
    geomPosTool = 0;
  }

  return this->GaudiTool::initialize();
}
StatusCode GtDiffuserBallTool::finalize ( ) [virtual]

Definition at line 216 of file GtDiffuserBallTool.cc.

{
  return this->GaudiTool::finalize();
}
StatusCode GtDiffuserBallTool::mutate ( HepMC::GenEvent event) [virtual]

Modify the event.

Implements IHepMCEventMutator.

Definition at line 221 of file GtDiffuserBallTool.cc.

{ 
  verbose() << "Mode: " << m_particlesPerEventMode << endreq;
  int particleNumber = m_particlesPerEvent;
  debug() << " particleNumber " << particleNumber << " mode " << m_particlesPerEventMode << endreq;
  if(m_particlesPerEventMode != "Fixed"){
    // Smear particle number
    particleNumber = m_particlesPerEvent + int(m_particlesPerEventSpread/m_photonScaleWeight * m_randPartNo()); 
    debug() << " smeared # of particles " << particleNumber << endreq;
    // Make sure that at least one dummy photon is generated
    if ( particleNumber < 1 ) { 
      particleNumber = 1; 
      debug() << " after requiring >0 photons. particleNumber is " << particleNumber << endreq;
    }
    else {
      // Cap photon number
      if ( particleNumber > int(m_particlesPerEvent * 4)) { particleNumber = int(m_particlesPerEvent * 4); }   
      debug() << " after cap. particleNumber is " << particleNumber << endreq;
    }
    debug() << "Generate " << particleNumber << " photons" << endreq;
  }
  for (int ind=0; ind<particleNumber; ++ind) {
    if (this->oneVertex(event).isFailure()) 
      return StatusCode::FAILURE;
  }
  return StatusCode::SUCCESS;
}
StatusCode GtDiffuserBallTool::oneVertex ( HepMC::GenEvent event) [private]

Definition at line 249 of file GtDiffuserBallTool.cc.

{
  // generate vertex time (within single event)
  double vertex_t = GetTime();

  // generate vertex position (random point on ball from which photon is emitted)
  double vertex_costh;
  double vertex_phi = m_uni()*CLHEP::twopi;  // range from 0 to 2*pi radians uniformly
  
  if (m_anisotropy){
    // create bright spot at top of diffuser ball, theta=0
    // theta=0 is 50% brighter than theta=pi (brightness taken to be linear)
    // see DocDB-1212

    bool accept = false;
    double vertex_theta;
    while (!accept){
      vertex_theta = m_uni()*CLHEP::pi; 
      double test  = m_uni()*1.25; // function below has max of 1.25 
      if (test <= (-vertex_theta/twopi + 1.5)*sin(vertex_theta)) 
        accept = true;
    }

    // rotate vertex to random bright spot location
    CLHEP::Hep3Vector vertex_bright;
    vertex_bright.setRThetaPhi(1.0,vertex_theta,vertex_phi);
    vertex_bright.rotateY(m_bright_theta);
    vertex_bright.rotateZ(m_bright_phi);
    vertex_bright.unit();

    vertex_costh = vertex_bright.cosTheta();
    vertex_phi   = vertex_bright.phi();
  }else{
    vertex_costh = 2*m_uni()-1;  // cos(theta) ranges from -1 to 1 uniformly
  }

  double vertex_sinth = sqrt(1-vertex_costh*vertex_costh);
  double vertex_x = m_radius*vertex_sinth*cos(vertex_phi);
  double vertex_y = m_radius*vertex_sinth*sin(vertex_phi);
  double vertex_z = m_radius*vertex_costh;
  HepMC::GenVertex* vertex = new HepMC::GenVertex(HepMC::FourVector(vertex_x,
                                                                    vertex_y,
                                                                    vertex_z,
                                                                    vertex_t));
  // Catch scaling of photon weight
  if( m_photonScaleWeight > 1.0 ){
     vertex->weights().push_back( m_photonScaleWeight );
  }
  if( !event.add_vertex(vertex) ){
    error() << "Failed to add vertex to event" << endreq;
    return StatusCode::FAILURE;
  }

  // generate momentum direction
  CLHEP::Hep3Vector dir;
  double dir_costh = m_uni(); // cos(theta) range 0 to 1
                              // do not want photon traveling into diffuser ball
  double dir_phi = m_uni()*twopi; // range from 0 to 2*pi radians
  
  dir.setRThetaPhi(1.0,acos(dir_costh),dir_phi);

  // rotate along Y-axis by vertex direction's theta
  dir.rotateY(acos(vertex_costh));
  // rotate along Z-axis by vertex direction's phi
  dir.rotateZ(vertex_phi);
  
  dir = dir.unit();
  
  // generate wavelength/momentum/energy
  double momentum = CLHEP::twopi*CLHEP::hbarc / GetWave();
  double energy = momentum;
  HepMC::FourVector fourmom(momentum*dir.x(),
                            momentum*dir.y(),
                            momentum*dir.z(),
                            energy);
  
  HepMC::GenParticle* particle = new HepMC::GenParticle(fourmom,m_pid,1/*=status*/);
  double pol_phi = m_uni()*CLHEP::twopi;
  particle->set_polarization(HepMC::Polarization(0.5*pi,pol_phi));
    
  vertex->add_particle_out(particle);

  return StatusCode::SUCCESS;
}
double GtDiffuserBallTool::GetWave ( ) [private]

Definition at line 345 of file GtDiffuserBallTool.cc.

{
  double wave = 0.0;
  double current_rand = m_randWave();
  double mean = m_wavelength;
  double spread = m_wavelengthSpread;

  if (m_wavelengthMode == "Fixed")        wave = m_wavelength;
  else if (m_wavelengthMode == "Uniform") wave = (2*current_rand-1)*spread + mean;
  else if (m_wavelengthMode == "Smeared") wave = current_rand*spread + mean;
  else if (m_wavelengthMode == "Defined") wave = ConvertCdfRand(current_rand,m_cdfWave,m_edgesWave); 

  return wave;
}
double GtDiffuserBallTool::GetTime ( ) [private]

Definition at line 334 of file GtDiffuserBallTool.cc.

{
  double time = -1.0;
  double current_rand = m_randTime();

  if (m_timingMode == "Instant")      time = 0.0;
  else if (m_timingMode == "Defined") time = ConvertCdfRand(current_rand,m_cdfTime,m_edgesTime);

  return time;
}
double GtDiffuserBallTool::ConvertCdfRand ( const double &  rand,
const std::vector< double > &  cdf,
const std::vector< double > &  edges 
) [private]

Definition at line 360 of file GtDiffuserBallTool.cc.

                                                                          {
  // Convert a uniform random number in [0,1] to a random sample from
  // the given cumulative distribution function.
  int left = 0;
  int right = cdf.size();
  const double* cdfStart = &cdf[0];
  const double* edgesStart = &edges[0];

  while( right-left > 1 ){
    int mid = (right+left)/2;
    if( rand < *(cdfStart+mid) ) right = mid;
    else left = mid;
  }
  double value = *(edgesStart + left);
  double slope = (*(edgesStart+right) - *(edgesStart+left))
                /(*(cdfStart+right) - *(cdfStart+left));
  value += slope * (rand - *(cdfStart+left));
  return value;
}
const InterfaceID & IHepMCEventMutator::interfaceID ( ) [static, inherited]

Retrieve interface ID.

Definition at line 8 of file IHepMCEventMutator.cc.

{ 
    return IID_IHepMCEventMutator; 
}

Member Data Documentation

Definition at line 68 of file GtDiffuserBallTool.h.

Definition at line 74 of file GtDiffuserBallTool.h.

Definition at line 76 of file GtDiffuserBallTool.h.

std::string GtDiffuserBallTool::m_wavelengthMode [private]

Definition at line 83 of file GtDiffuserBallTool.h.

std::string GtDiffuserBallTool::m_timingMode [private]

Definition at line 88 of file GtDiffuserBallTool.h.

Definition at line 91 of file GtDiffuserBallTool.h.

Definition at line 94 of file GtDiffuserBallTool.h.

Definition at line 97 of file GtDiffuserBallTool.h.

double GtDiffuserBallTool::m_radius [private]

Definition at line 100 of file GtDiffuserBallTool.h.

Definition at line 103 of file GtDiffuserBallTool.h.

Definition at line 105 of file GtDiffuserBallTool.h.

Definition at line 106 of file GtDiffuserBallTool.h.

Definition at line 111 of file GtDiffuserBallTool.h.

Definition at line 114 of file GtDiffuserBallTool.h.

std::string GtDiffuserBallTool::m_particleName [private]

Definition at line 116 of file GtDiffuserBallTool.h.

Definition at line 117 of file GtDiffuserBallTool.h.

Rndm::Numbers GtDiffuserBallTool::m_uni [private]

Definition at line 119 of file GtDiffuserBallTool.h.

Rndm::Numbers GtDiffuserBallTool::m_randWave [private]

Definition at line 119 of file GtDiffuserBallTool.h.

Rndm::Numbers GtDiffuserBallTool::m_randTime [private]

Definition at line 119 of file GtDiffuserBallTool.h.

Rndm::Numbers GtDiffuserBallTool::m_randPartNo [private]

Definition at line 119 of file GtDiffuserBallTool.h.

std::vector<double> GtDiffuserBallTool::m_pdfWave [private]

Definition at line 123 of file GtDiffuserBallTool.h.

std::vector<double> GtDiffuserBallTool::m_edgesWave [private]

Definition at line 124 of file GtDiffuserBallTool.h.

std::vector<double> GtDiffuserBallTool::m_cdfWave [private]

Definition at line 125 of file GtDiffuserBallTool.h.

std::vector<double> GtDiffuserBallTool::m_pdfTime [private]

Definition at line 129 of file GtDiffuserBallTool.h.

std::vector<double> GtDiffuserBallTool::m_edgesTime [private]

Definition at line 130 of file GtDiffuserBallTool.h.

std::vector<double> GtDiffuserBallTool::m_cdfTime [private]

Definition at line 131 of file GtDiffuserBallTool.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:18:51 for GenTools by doxygen 1.7.4