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

In This Package:

Classes | Public Member Functions | Protected Types | Protected Attributes | Private Member Functions | Private Attributes
UnObserverStepAction Class Reference

#include <UnObserverStepAction.h>

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

List of all members.

Classes

struct  NamedStat
struct  StatGroup

Public Member Functions

 UnObserverStepAction (const std::string &type, const std::string &name, const IInterface *parent)
virtual ~UnObserverStepAction ()
virtual StatusCode initialize ()
virtual StatusCode finalize ()
virtual void UserSteppingAction (const G4Step *)
virtual StatusCode GetTrackParameterList (RuleParser::ParameterList &)
 To get the initial list of Queriables.
virtual StatusCode GetVertexParameterList (RuleParser::ParameterList &)
void Reset (const G4Step *, const G4Track *, const G4StepPoint *)
 To prepare a Step for querying:
virtual void queryParam (int id, double &output) const
virtual void queryParam (int id, std::string &output) const
const IDetectorElement * getDetectorElement () const
int getDetectorElementMatch () const
const ILVolume * getLogicalVolume () const
const IPVolume * getPhysicalVolume () const
const DayaBay::SimProcessgetProcess () const
double getQuenchedEnergy () const
const G4HistoryUserTrackInfogetUserTrackInfo () const
unsigned int getDetectorId () const
unsigned int getDetectorId (const IDetectorElement *) const

Protected Types

enum  EParamIds {
  kPar_Unknown = 0, kPar_START_VERTEX, kPar_Step_dE, kPar_Step_dE_Ion,
  kPar_Step_qdE, kPar_Step_dx, kPar_Step_dt, kPar_Step_dAngle,
  kPar_capTargetZ, kPar_capTargetA, kPar_Step_E_weighted_x, kPar_Step_E_weighted_y,
  kPar_Step_E_weighted_z, kPar_Step_E_weighted_t, kPar_Step_qE_weighted_x, kPar_Step_qE_weighted_y,
  kPar_Step_qE_weighted_z, kPar_Step_qE_weighted_t, kPar_IsStopping, kPar_IsStarting,
  kPar_StepNumber, kPar_VolumeChanged, kPar_MaterialChanged, kPar_START_TRACK,
  kPar_t, kPar_x, kPar_y, kPar_z,
  kPar_r, kPar_local_x, kPar_local_y, kPar_local_z,
  kPar_local_r, kPar_LogicalVolumeName, kPar_MaterialName, kPar_DetectorElement,
  kPar_DetectorElementName, kPar_DetectorElementMatch, kPar_NicheId, kPar_DetectorId,
  kPar_SiteId, kPar_Site, kPar_AD, kPar_momentum,
  kPar_TotEnergy, kPar_KineticEnergy, kPar_vx, kPar_vy,
  kPar_vz, kPar_local_vx, kPar_local_vy, kPar_local_vz,
  kPar_ProcessType, kPar_ProcessName, kPar_Pdg, kPar_Charge,
  kPar_TrackId, kPar_CreatorPdg, kPar_AncestorPdg, kPar_mass,
  kPar_ParticleName, kPar_Prescale, kPar_CreatorProcessName, kPar_END_VERTEX,
  kPar_END_QUERIABLE
}
typedef std::map< const
IDetectorElement *, unsigned
int > 
DetectorIdCache_t

Protected Attributes

std::vector< std::string > m_DetectorElementSearchPath
std::string m_TouchableToDetelem_name
std::vector< std::string > m_IdParameterNames
const G4Step * mCurrentStep
const G4Track * mCurrentTrack
const G4StepPoint * mCurrentStepPoint
ITouchableToDetectorElementmTouchToDetElem
const IDetectorElement * mDetElement
int mDetElementMatch
const ILVolume * mLogicVolume
const IPVolume * mPhysVolume
DayaBay::SimProcess mProcess
double mQuenchedEnergy
G4HistoryUserTrackInfomUserTrackInfo
unsigned int mDetectorId
DetectorIdCache_t mDetectorIdCache

Private Member Functions

 UnObserverStepAction ()
 no default constructor
 UnObserverStepAction (const UnObserverStepAction &)
 no copy
UnObserverStepActionoperator= (const UnObserverStepAction &)
 no =

Private Attributes

std::vector< std::vector
< std::string > > 
m_Stats
std::vector< StatGroupmGroups
IHistoryKeepermKeeper
int mBail

Detailed Description

Definition at line 11 of file UnObserverStepAction.h.


Member Typedef Documentation

typedef std::map<const IDetectorElement*,unsigned int> QueriableStepAction::DetectorIdCache_t [protected, inherited]

Definition at line 88 of file QueriableStepAction.h.


Member Enumeration Documentation

enum QueriableStepAction::EParamIds [protected, inherited]
Enumerator:
kPar_Unknown 
kPar_START_VERTEX 
kPar_Step_dE 
kPar_Step_dE_Ion 
kPar_Step_qdE 
kPar_Step_dx 
kPar_Step_dt 
kPar_Step_dAngle 
kPar_capTargetZ 
kPar_capTargetA 
kPar_Step_E_weighted_x 
kPar_Step_E_weighted_y 
kPar_Step_E_weighted_z 
kPar_Step_E_weighted_t 
kPar_Step_qE_weighted_x 
kPar_Step_qE_weighted_y 
kPar_Step_qE_weighted_z 
kPar_Step_qE_weighted_t 
kPar_IsStopping 
kPar_IsStarting 
kPar_StepNumber 
kPar_VolumeChanged 
kPar_MaterialChanged 
kPar_START_TRACK 
kPar_t 
kPar_x 
kPar_y 
kPar_z 
kPar_r 
kPar_local_x 
kPar_local_y 
kPar_local_z 
kPar_local_r 
kPar_LogicalVolumeName 
kPar_MaterialName 
kPar_DetectorElement 
kPar_DetectorElementName 
kPar_DetectorElementMatch 
kPar_NicheId 
kPar_DetectorId 
kPar_SiteId 
kPar_Site 
kPar_AD 
kPar_momentum 
kPar_TotEnergy 
kPar_KineticEnergy 
kPar_vx 
kPar_vy 
kPar_vz 
kPar_local_vx 
kPar_local_vy 
kPar_local_vz 
kPar_ProcessType 
kPar_ProcessName 
kPar_Pdg 
kPar_Charge 
kPar_TrackId 
kPar_CreatorPdg 
kPar_AncestorPdg 
kPar_mass 
kPar_ParticleName 
kPar_Prescale 
kPar_CreatorProcessName 
kPar_END_VERTEX 
kPar_END_QUERIABLE 

Definition at line 93 of file QueriableStepAction.h.

                  {
     kPar_Unknown = 0
     , kPar_START_VERTEX
     // Applies to end-of-step (Vertex) only
     , kPar_Step_dE                    // double
     , kPar_Step_dE_Ion                // double
     , kPar_Step_qdE                   // double
     , kPar_Step_dx                    // double
     , kPar_Step_dt                    // double
     , kPar_Step_dAngle                // double
     , kPar_capTargetZ                 // int
     , kPar_capTargetA                 // int | No correct A info. 8/2008 Wei
     
     , kPar_Step_E_weighted_x          // double
     , kPar_Step_E_weighted_y          // double
     , kPar_Step_E_weighted_z          // double
     , kPar_Step_E_weighted_t          // double
     , kPar_Step_qE_weighted_x         // double
     , kPar_Step_qE_weighted_y         // double
     , kPar_Step_qE_weighted_z         // double
     , kPar_Step_qE_weighted_t         // double

     , kPar_IsStopping                 // int
     , kPar_IsStarting                 // int
     , kPar_StepNumber                 // int

     , kPar_VolumeChanged              // int as bool
     , kPar_MaterialChanged            // int as bool

     // Both Vertex and Track:
     , kPar_START_TRACK
     , kPar_t                          // double   |
     , kPar_x                          // double   |
     , kPar_y                          // double   |
     , kPar_z                          // double   |
     , kPar_r                          // double   |
     , kPar_local_x                    // double   |
     , kPar_local_y                    // double   |
     , kPar_local_z                    // double   | all of these guys refer to prepoint for tracks
     , kPar_local_r                    // double   | or postpoint for vertices.
     , kPar_LogicalVolumeName          // string
     , kPar_MaterialName               // string
     , kPar_DetectorElement            // custom
     , kPar_DetectorElementName        // string           
     , kPar_DetectorElementMatch       // int
     , kPar_NicheId                    // int
     , kPar_DetectorId                 // int
     , kPar_SiteId                     // int
     , kPar_Site                       // int
     , kPar_AD                         // int
     , kPar_momentum                   // double   |
     , kPar_TotEnergy                  // double   |
     , kPar_KineticEnergy              // double   |
     , kPar_vx                         // double   |
     , kPar_vy                         // double   |
     , kPar_vz                         // double   |
     , kPar_local_vx                   // double   |
     , kPar_local_vy                   // double   |
     , kPar_local_vz                   // double   |
     , kPar_ProcessType                // int      |
     , kPar_ProcessName                // string   |

     , kPar_Pdg                        // int 
     , kPar_Charge                    // int
     , kPar_TrackId                    // int
     , kPar_CreatorPdg                 // int
     , kPar_AncestorPdg                // int
     , kPar_mass                       // double
     , kPar_ParticleName               // string
     , kPar_Prescale                   // custom
     , kPar_CreatorProcessName         // string
     , kPar_END_VERTEX
     , kPar_END_QUERIABLE
   }; 

Constructor & Destructor Documentation

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

Definition at line 26 of file UnObserverStepAction.cc.

  : QueriableStepAction( type , name , parent )
  , m_Stats()
{ 
   // Declare some good starting stats.
  m_Stats.push_back(AddableVector().add("photon_created_energy" ).add( "E" ).add( "StepNumber==1 and pdg==20022" ) );
  m_Stats.push_back(AddableVector().add("photon_backscatter_r" ).add( "local_r" ).add( "pdg==20022 and dAngle >= 90" ) );
  m_Stats.push_back(AddableVector().add("photon_forescatter_r" ).add( "local_r" ).add( "pdg==20022 and (dAngle <= 90 and dAngle > 0)" ) );
  m_Stats.push_back(AddableVector().add("photon_stop_r" ).add(        "local_r" ).add( "(pdg==20022) and ((IsStopping == 1) and (LogicalVolumeName == '/dd/Geometry/PMT/lvPmtHemi'))" ) );
  m_Stats.push_back(AddableVector().add("pmt_hit" ).add( "t" ).add( "pdg==20022 and (IsStopping == 1 and (LogicalVolumeName=='/dd/Geometry/PMT/lvPmtHemi'))" ) );
  m_Stats.push_back(AddableVector().add("edep-water").add(  "dE").add( "pdg!=20022 and (MaterialName == '/dd/Materials/Water')") );
  m_Stats.push_back(AddableVector().add("edep-oil").add(    "dE").add( "pdg!=20022 and (MaterialName == '/dd/Materials/Oil')") );
  m_Stats.push_back(AddableVector().add("edep-ls").add(     "dE").add( "pdg!=20022 and (MaterialName == '/dd/Materials/LiquidScintillator')") );
  m_Stats.push_back(AddableVector().add("edep-gdls").add(   "dE").add( "pdg!=20022 and (MaterialName == '/dd/Materials/GdDopedLS')") );
  m_Stats.push_back(AddableVector().add("edep-acryilc").add("dE").add( "pdg!=20022 and (MaterialName == '/dd/Materials/Acrylic')") );
  m_Stats.push_back(AddableVector().add("edep-ad1"  ).add("dE"  ).add("pdg!=20022 and ((MaterialName == '/dd/Materials/LiquidScintillator' or MaterialName == '/dd/Materials/GdDopedLS') and AD==1)"  ) );
  m_Stats.push_back(AddableVector().add("edep-ad2"  ).add("dE"  ).add("pdg!=20022 and ((MaterialName == '/dd/Materials/LiquidScintillator' or MaterialName == '/dd/Materials/GdDopedLS') and AD==2)"  ) );
  m_Stats.push_back(AddableVector().add("edep-ad3"  ).add("dE"  ).add("pdg!=20022 and ((MaterialName == '/dd/Materials/LiquidScintillator' or MaterialName == '/dd/Materials/GdDopedLS') and AD==3)"  ) );
  m_Stats.push_back(AddableVector().add("edep-ad4"  ).add("dE"  ).add("pdg!=20022 and ((MaterialName == '/dd/Materials/LiquidScintillator' or MaterialName == '/dd/Materials/GdDopedLS') and AD==4)"  ) );
  m_Stats.push_back(AddableVector().add("qedep-ad1" ).add("qdE" ).add("pdg!=20022 and ((MaterialName == '/dd/Materials/LiquidScintillator' or MaterialName == '/dd/Materials/GdDopedLS') and AD==1)"   ) );
  m_Stats.push_back(AddableVector().add("qedep-ad2" ).add("qdE" ).add("pdg!=20022 and ((MaterialName == '/dd/Materials/LiquidScintillator' or MaterialName == '/dd/Materials/GdDopedLS') and AD==2)"   ) );
  m_Stats.push_back(AddableVector().add("qedep-ad3" ).add("qdE" ).add("pdg!=20022 and ((MaterialName == '/dd/Materials/LiquidScintillator' or MaterialName == '/dd/Materials/GdDopedLS') and AD==3)"   ) );
  m_Stats.push_back(AddableVector().add("qedep-ad4" ).add("qdE" ).add("pdg!=20022 and ((MaterialName == '/dd/Materials/LiquidScintillator' or MaterialName == '/dd/Materials/GdDopedLS') and AD==4)"   ) );
  m_Stats.push_back(AddableVector().add("mu-path-water").add( "dx").add( "(pdg==13 || pdg==-13) and MaterialName == '/dd/Materials/Water'") );
  m_Stats.push_back(AddableVector().add("mu-path-oil").add(   "dx").add( "(pdg==13 || pdg==-13) and MaterialName == '/dd/Materials/Oil'") );
  m_Stats.push_back(AddableVector().add("mu-path-ls").add(    "dx").add( "(pdg==13 || pdg==-13) and MaterialName == '/dd/Materials/LiquidScintillator'") );
  m_Stats.push_back(AddableVector().add("mu-path-gdls").add(  "dx").add( "(pdg==13 || pdg==-13) and MaterialName == '/dd/Materials/GdDopedLS'") );
  m_Stats.push_back(AddableVector()
    .add("ad"    ).add( "AD" )
    .add("expos" ).add("Ex"  )
    .add("eypos" ).add("Ey"  )
    .add("ezpos" ).add("Ez"  )
    .add("et"    ).add("Et"  )
    .add("qexpos").add("qEx" )
    .add("qeypos").add("qEy" )
    .add("qezpos").add("qEz" )
    .add("qet"   ).add("qEt" )
    .add( "(pdg!=20022) and (MaterialName == '/dd/Materials/LiquidScintillator' or MaterialName == '/dd/Materials/GdDopedLS')"  ) );
   
   declareProperty("Stats",m_Stats);
   
   
   // Also, change the default search path:
   m_DetectorElementSearchPath =  AddableVector()
     .add("/dd/Structure/DayaBay")
     .add("/dd/Structure/Sites")
     .add("/dd/Structure/Pool")
     .add("/dd/Structure/AD")
     .add("/dd/Structure/RPC");
   
   
}
UnObserverStepAction::~UnObserverStepAction ( ) [virtual]

Definition at line 82 of file UnObserverStepAction.cc.

{
  for(size_t i=0;i<mGroups.size();i++) {
    if(mGroups[i].mRule) delete mGroups[i].mRule;
  }   
}
UnObserverStepAction::UnObserverStepAction ( ) [private]

no default constructor

UnObserverStepAction::UnObserverStepAction ( const UnObserverStepAction ) [private]

no copy


Member Function Documentation

StatusCode UnObserverStepAction::initialize ( ) [virtual]

Reimplemented from QueriableStepAction.

Definition at line 90 of file UnObserverStepAction.cc.

{
  info() << "UnObserverStepAction::initialize()" << endreq;

  StatusCode sc = QueriableStepAction::initialize();
  if (sc.isFailure()) return sc;
  
  sc = service("HistoryKeeper",mKeeper,true);
  if (sc.isFailure()) return sc;

  // Build the selection rules.
  
  using namespace RuleParser;
  using std::string;

  ParameterList plist;
  QueriableStepAction::GetVertexParameterList(plist);

  // Now loop through the various statistics that have been offered.
  for(size_t istat = 0; istat < m_Stats.size(); istat++ ){
    std::vector<std::string>& stat = m_Stats[istat];

    // Check basic layout. Odd number of strings.
    if(stat.size()%2 != 1) {
      err() << "Bad configuration given to the UnObserverStepAction:" << endreq;
      err() << " Expected strings of the form <name> : <parameter> [: <name> : <parameter> ..]: <cut>" << endreq;
      err() << " but got " << stat.size() << " strings"
      //<< " \"" << line << << "\""
       <<  endreq;
      return StatusCode::FAILURE;
    }

    StatGroup statgroup;
    
    // Compile the cut
    statgroup.mRule = 0;
    bool compiled = RuleParser::CreateRules(stat[stat.size()-1],plist,statgroup.mRule);
    if(!compiled || !statgroup.mRule) {
      err() << "Bad configuration given to the UnObserverStepAction:" << endreq;
      err() << " Expected 3 strings of the form <name>,<parameter>,<cut>" << endreq;
      err() << " but I failed to interpret selection string: " << stat[stat.size()-1] << endreq;
      return StatusCode::FAILURE;      
    }
    info() << "Cut group : " << stat[stat.size()-1] << endreq;
    
    for(size_t isub=0; isub<(stat.size()-1)/2 ; isub+=2) {
      // Check that the paramter name matches something, and record the match.
      std::string parname = stat[isub+1];
      for(std::string::iterator c = parname.begin(); c!= parname.end(); ++c)  *c = tolower(*c);
    
      bool foundit = false;
      int parid = 0;
      for(size_t ipar = 0; ipar < plist.size(); ipar++) {
        std::string temp = plist[ipar].name();
        for(std::string::iterator c = temp.begin(); c!= temp.end(); ++c)  *c = tolower(*c);
        if(parname==temp) {
          foundit = true;
          parid = plist[ipar].id();
        }
      }
      if(!foundit) {   
        err() << "Bad configuration given to the UnObserverStepAction:" << endreq; 
        err() << "Didn't recognize parameter name " << parname << endreq;
        return StatusCode::FAILURE;
      }
      // Report success and record the info we'll need at runtime.
      info() << "   Statistic " << stat[isub] 
             << " param=" << parname << " (id=" << parid <<")"
             << endreq;
       
      NamedStat ns;
      ns.mName = stat[isub];
      ns.mParam = parid;

      statgroup.push_back(ns);
    }
    if(statgroup.size())
     mGroups.push_back(statgroup);
  }

  // keep track of bail-out messages
  mBail = 0; 

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

Reimplemented from QueriableStepAction.

Definition at line 178 of file UnObserverStepAction.cc.

{
  info() << "UnObserverStepAction::finalize()" << endreq;
  info() << "There were "<<mBail<<" bail-out warnings encountered in UserSteppingAction."<<endreq;
  return QueriableStepAction::finalize(); 
}
void UnObserverStepAction::UserSteppingAction ( const G4Step *  g4step) [virtual]

This is the meat of the history-constructing code. Although the TrackAction is responsible for setting up primary tracks, every daughter track is set up through this interface.

Implements QueriableStepAction.

Definition at line 187 of file UnObserverStepAction.cc.

{

  // Use the PreStepPoint to correctly evaluate geometry-, material-, volume-related information
  // Use the PostStepPoint to evaluate 'track'-related infomation (ie, energy loss, step length)

  // Reference material:
  // The algorithm to handle a single step in G4:
  // http://geant4.web.cern.ch/geant4/UserDocumentation/UsersGuides/ForApplicationDeveloper/html/ch05.html#sect.Track.Basic
  // G4 stepping info from http://geant4.cern.ch/support/faq.shtml#TRACK-1 :
  // "Hereafter we call current volume the volume where the step has just 
  //   gone through. Geometrical informations are available from preStepPoint."
  // 
  // Reset temporary things:
  QueriableStepAction::Reset(g4step,g4step->GetTrack(),g4step->GetPostStepPoint());
  if ( g4step->GetTrack()->GetCurrentStepNumber() == 1 )
    // This is the first step.  
    // The interesting point is the first one.
    QueriableStepAction::Reset(g4step,g4step->GetTrack(),g4step->GetPreStepPoint());

  // The following is useful for intensive debugging of each step on every track.
  // Make it so it does not get executed without re-compilation.
  if ( 0 && g4step->GetTrack()->GetDefinition()->GetParticleName() != "opticalphoton" )
    {
      G4StepPoint* pre = g4step->GetPreStepPoint();
      G4String mpre = "NONE", mpost = "NONE", status1 = "", status2 = "";
      G4String ppre = "NONE",ppost = "NONE",lpre = "NONE",lpost = "NONE";
      G4ThreeVector localPosition = 0, worldPosition = 0;
      if (pre) {
        if ( pre->GetStepStatus() == fGeomBoundary ) status1 = "Entering";
        G4Material* mm = pre->GetMaterial();
        if (mm) mpre = mm->GetName();
        G4TouchableHandle touch = pre->GetTouchableHandle();
        G4VPhysicalVolume* volume = touch->GetVolume();
        ppre = volume->GetName();
        G4LogicalVolume* lvolume = volume->GetLogicalVolume();
        lpre = lvolume->GetName();
      }
      else {
        verbose() << " THERE IS NO PRE-STEP-POINT FOR THIS STEP " << endreq;
      }
      pre = g4step->GetPostStepPoint();
      if (pre) {
        if (pre->GetStepStatus() == fGeomBoundary ) status2 = "Leaving";
        G4Material* mm = pre->GetMaterial();
        if (mm) mpost = mm->GetName();
        G4TouchableHandle touch = pre->GetTouchableHandle();
        if (touch) {
          worldPosition = pre->GetPosition();
          localPosition = touch->GetHistory()->GetTopTransform().TransformPoint( worldPosition );
          G4VPhysicalVolume* volume = touch->GetVolume();
          if (volume) {
            ppost = volume->GetName();
            G4LogicalVolume* lvolume = volume->GetLogicalVolume();
            if (lvolume)          lpost = lvolume->GetName();
          }
        }
      }
      else {
        verbose() << " THERE IS NO POST-STEP-POINT FOR THIS STEP " << endreq;
      }

      verbose() << "Track in " << g4step->GetTrack()->GetVolume()->GetName() 
              << ", particle is " << g4step->GetTrack()->GetDefinition()->GetParticleName()
              << " Etot(MeV) " << g4step->GetTrack()->GetTotalEnergy()/MeV
              << " Ekin(MeV) " << g4step->GetTrack()->GetKineticEnergy()/MeV
              << " worldX,Y,Z " << worldPosition.x() << ", " << worldPosition.y() << ", " << worldPosition.z()
              << " localX,Y,Z " << localPosition.x() << ", " << localPosition.y() << ", " << localPosition.z()
              << " DirX,Y,Z " << g4step->GetTrack()->GetMomentumDirection().x()
              << ", " << g4step->GetTrack()->GetMomentumDirection().y()
              << ", " << g4step->GetTrack()->GetMomentumDirection().z()
              << " G4TrackID "<< g4step->GetTrack()->GetTrackID()
              << " current step# " << g4step->GetTrack()->GetCurrentStepNumber()
              << " step length(mm) " << g4step->GetTrack()->GetStepLength()/mm
              << " status " << g4step->GetTrack()->GetTrackStatus() << endreq;
      verbose() << " ...." << status1 << " " << status2 
              << " prestep (current) material " << mpre << " pVol " << ppre << " lVol " << lpre 
              << " poststep material " << mpost << " pVol " << ppost << " lVol " << lpost
              << endreq;
      }
        
  // djaffe Bail with warning if we are in uppermost 'universe' volume
  if ( g4step->GetTrack()->GetVolume()->GetMotherLogical() == 0x0 ) {
    mBail++;
    if (mBail <= 10) { 
      warning() << "Bailing out because track is in " << g4step->GetTrack()->GetVolume()->GetName() 
                << ", particle is " << g4step->GetTrack()->GetDefinition()->GetParticleName()
                << " with total energy(MeV) " << g4step->GetTrack()->GetTotalEnergy()/MeV
                << " & kinetic energy(MeV) " << g4step->GetTrack()->GetKineticEnergy()/MeV
                << " G4TrackID "<< g4step->GetTrack()->GetTrackID()
                << " current step# " << g4step->GetTrack()->GetCurrentStepNumber()
                << " step length(mm) " << g4step->GetTrack()->GetStepLength()/mm
                << " status " << g4step->GetTrack()->GetTrackStatus()
                << endreq;
    }
    if (mBail > 10) {
      debug()   << "Bailing out because track is in " << g4step->GetTrack()->GetVolume()->GetName() 
                << ", particle is " << g4step->GetTrack()->GetDefinition()->GetParticleName()
                << " with total energy(MeV) " << g4step->GetTrack()->GetTotalEnergy()/MeV
                << " & kinetic energy(MeV) " << g4step->GetTrack()->GetKineticEnergy()/MeV
                << " G4TrackID "<< g4step->GetTrack()->GetTrackID()
                << " current step# " << g4step->GetTrack()->GetCurrentStepNumber()
                << " step length(mm) " << g4step->GetTrack()->GetStepLength()/mm
                << " status " << g4step->GetTrack()->GetTrackStatus()
                << endreq;
    }
    if (mBail == 10) { 
      warning() << "Additional bail-out messages will be suppressed" << endreq;
    }
    return;
  }
  
  // Get the current set of stats.
  SimUnobservableStatisticsHeader* ush = 0;
  StatusCode sc = mKeeper->GetCurrentUnobservable(ush);
  if(sc.isFailure() || (!ush)) {
    err() << "Problem with Keeper." << endreq;
    return;
  }
  assert(ush);
  
  SimUnobservableStatisticsHeader::stat_map& statmap = ush->stats();
  if(ush->stats().size() == 0) {

    // This is a first call.  Make sure that all the stats get intialized, by simply looking them up.
    for(size_t igrp = 0; igrp < mGroups.size(); igrp++) {
      StatGroup& statgroup = mGroups[igrp];
      for(size_t istat = 0; istat < statgroup.size(); istat++ ){
        statmap[statgroup[istat].mName];
        verbose() << "Initializing group " << igrp << " stat " << istat 
               << " name=" <<statgroup[istat].mName << " paramid=" << statgroup[istat].mParam << endreq;
      }
    }  
  }

  for(size_t igrp = 0; igrp < mGroups.size(); igrp++) {
    StatGroup& statgroup = mGroups[igrp];
    
    // The magic line.. does this step pass the cut?
    verbose() << " Group " << igrp << " statgroup.mRule->name() " << statgroup.mRule->name() << endreq;
    if(statgroup.mRule->select(this)) {
      // yes, so loop over the parameters.
      for(size_t istat = 0; istat < statgroup.size(); istat++ ){
        double x = 0;
        queryParam(statgroup[istat].mParam,x);
        SimStatistic& stat = statmap[statgroup[istat].mName];
        verbose() << "Querying group " << igrp << " stat " << istat 
                << " name=" <<statgroup[istat].mName << " paramid=" << statgroup[istat].mParam 
                << " and incrementing by " << x << endreq;
        stat.increment(x);
      }
    }
  }

}
UnObserverStepAction& UnObserverStepAction::operator= ( const UnObserverStepAction ) [private]

no =

StatusCode QueriableStepAction::GetTrackParameterList ( RuleParser::ParameterList tplist) [virtual, inherited]

To get the initial list of Queriables.

Usually called at initialize()

Definition at line 107 of file QueriableStepAction.cc.

{
  using namespace RuleParser;
  using std::string;
  
  tplist.add<double>(kPar_t,"time","t");
  tplist.add<double>(kPar_x,"x","global_x");
  tplist.add<double>(kPar_y,"y","global_y");
  tplist.add<double>(kPar_z,"z","global_z");
  tplist.add<double>(kPar_r,"r","radius","pos_r");
  tplist.add<double>(kPar_local_x,"lx","local_x","det_x");             
  tplist.add<double>(kPar_local_y,"ly","local_y","det_y");
  tplist.add<double>(kPar_local_z,"lz","local_z","det_z");
  tplist.add<double>(kPar_local_r,"lr","local_r","det_r");
  tplist.add<string>(kPar_LogicalVolumeName,"Volume","LogicalVolumeName","LogicalVolume","VolumeName");
  tplist.add<string>(kPar_MaterialName,"Material","MaterialName");
  tplist.add<string>(kPar_DetectorElementName,"DetectorElementName");
  tplist.add<double>(kPar_DetectorElementMatch,"Match","DetectorElementMatch");
  tplist.add<double>(kPar_NicheId,"NicheId","Niche");
  tplist.add<double>(kPar_DetectorId,"DetectorId");
  tplist.add<double>(kPar_SiteId,"SiteId");
  tplist.add<double>(kPar_Site,"Site");
  tplist.add<double>(kPar_AD,"AD","AdNumber");

  tplist.add<double>(kPar_momentum,"momentum","p");
  tplist.add<double>(kPar_TotEnergy,"E","totEnergy","TotalEnergy");
  tplist.add<double>(kPar_KineticEnergy,"KE","kineticEnergy");
  tplist.add<double>(kPar_vx,"vx","global_dir_x","global_u");
  tplist.add<double>(kPar_vy,"vy","global_dir_y","global_v");
  tplist.add<double>(kPar_vz,"vz","global_dir_z","global_w");
  tplist.add<double>(kPar_local_vx,"lvx","local_dir_x","local_u");
  tplist.add<double>(kPar_local_vy,"lvy","local_dir_y","local_v");
  tplist.add<double>(kPar_local_vz,"lvz","local_dir_z","local_w");
  tplist.add<double>(kPar_ProcessType,"ProcessType");
  tplist.add<string>(kPar_ProcessName,"Process","ProcessName");
  tplist.add<double>(kPar_Pdg,"pdg","pdgcode","particle");
  tplist.add<double>(kPar_Charge,"charge","ParticleCharge","q");
  tplist.add<double>(kPar_TrackId,"id","trackid");
  tplist.add<double>(kPar_CreatorPdg,"creatorPdg","creator");
  tplist.add<double>(kPar_AncestorPdg,"ParentPdg","AncestorPdg","Ancestor");
  tplist.add<double>(kPar_mass,"mass","m");
  tplist.add<string>(kPar_ParticleName,"ParticleName");
  tplist.add<string>(kPar_CreatorProcessName,"CreatorProcessName","CreatorProcess");



  boost::shared_ptr<RuleFactory> prescaleFactory(new PrescaleRuleFactory);
  boost::shared_ptr<RuleFactory> customFactory(new DetElemContainsRuleFactory<QueriableStepAction>(this));
  tplist.add<int>(kPar_Prescale,"Prescale","by",prescaleFactory);
  tplist.add<const IDetectorElement*>(kPar_DetectorElement,"DetElem","in",customFactory);
  tplist.add<const IDetectorElement*>(kPar_DetectorElement,"DetectorElement","in",customFactory);


  return StatusCode::SUCCESS;
}
StatusCode QueriableStepAction::GetVertexParameterList ( RuleParser::ParameterList vplist) [virtual, inherited]

Definition at line 164 of file QueriableStepAction.cc.

{
  GetTrackParameterList(vplist); // Everything valid for a track is valid for a vertex.

  using namespace RuleParser;
  using std::string;


  // And add vertex-specific stuff:
  vplist.add<double>(kPar_Step_dE,"Step_dE","dE");
  vplist.add<double>(kPar_Step_dE_Ion,"Step_dE_Ion","de_ion","ionization");
  vplist.add<double>(kPar_Step_qdE,"Step_qDE","quenched_dE","qdE");
  vplist.add<double>(kPar_Step_dx,"Step_dx","StepLength","dx");
  vplist.add<double>(kPar_Step_dt,"Step_dt","StepDuration","dt");
  vplist.add<double>(kPar_Step_dAngle,"Step_dAngle","dAngle");
  vplist.add<double>(kPar_capTargetZ,"capTargetZ");
  vplist.add<double>(kPar_capTargetA,"capTargetA");
  vplist.add<double>(kPar_Step_E_weighted_x,"Ex","E_weighted_x");
  vplist.add<double>(kPar_Step_E_weighted_y,"Ey","E_weighted_y"); 
  vplist.add<double>(kPar_Step_E_weighted_z,"Ez","E_weighted_z");    
  vplist.add<double>(kPar_Step_E_weighted_t,"Et","E_weighted_t");    
  vplist.add<double>(kPar_Step_qE_weighted_x,"qEx","qE_weighted_x","quenched_weighted_x");
  vplist.add<double>(kPar_Step_qE_weighted_y,"qEy","qE_weighted_y","quenched_weighted_y");    
  vplist.add<double>(kPar_Step_qE_weighted_z,"qEz","qE_weighted_z","quenched_weighted_z");    
  vplist.add<double>(kPar_Step_qE_weighted_t,"qEt","qE_weighted_t","quenched_weighted_t");                                     


  vplist.add<double>(kPar_IsStopping,"IsStopping","stop","End");
  vplist.add<double>(kPar_IsStarting,"IsStarting","start","begin");
  vplist.add<double>(kPar_StepNumber,"StepNumber");                    
  vplist.add<double>(kPar_VolumeChanged,"VolumeChanged","NewVolume");
  vplist.add<double>(kPar_MaterialChanged,"MaterialChanged","NewMaterial");


  verbose() << " kPar_MaterialName " << kPar_MaterialName << " kPar_MaterialChanged " << kPar_MaterialChanged
          << " kPar_VolumeChanged " << kPar_VolumeChanged << " kPar_Pdg " << kPar_Pdg 
          << " kPar_LogicalVolumeName " << kPar_LogicalVolumeName 
          << endreq;
  
  return StatusCode::SUCCESS;
}                                                   
void QueriableStepAction::Reset ( const G4Step *  step,
const G4Track *  track,
const G4StepPoint *  point 
) [inherited]

To prepare a Step for querying:

Definition at line 85 of file QueriableStepAction.cc.

{
  mCurrentStep = step;
  mCurrentTrack = track;
  mCurrentStepPoint = point;
  
  mDetElement=0;
  mDetElementMatch=-1;
  mLogicVolume=0;
  mPhysVolume=0;
  mProcess=SimProcess();
  mQuenchedEnergy = -1e9;
  mUserTrackInfo = 0;
  mDetectorId = 0xFFFFFFF;

  verbose() << " The current track has been reset " << endreq;
}
void QueriableStepAction::queryParam ( int  id,
double &  output 
) const [virtual, inherited]

Reimplemented from RuleParser::Queriable.

Reimplemented in HistorianStepAction.

Definition at line 211 of file QueriableStepAction.cc.

{

  switch(id) {
    case kPar_IsStopping: 
      output = (  ( mCurrentTrack->GetTrackStatus() == fStopAndKill)
               || ( mCurrentTrack->GetTrackStatus() == fKillTrackAndSecondaries ) );
      break;

    case kPar_IsStarting: 
      output = (mCurrentTrack->GetCurrentStepNumber() == 1);
      break;

    case kPar_StepNumber:
      output = mCurrentTrack->GetCurrentStepNumber();
      break;
 
    case kPar_MaterialChanged:
      output = (mCurrentStep->GetPostStepPoint()->GetMaterial() != mCurrentStep->GetPreStepPoint()->GetMaterial());
      break;

    case kPar_VolumeChanged:
      output = (mCurrentStep->GetPostStepPoint()->GetStepStatus() == fGeomBoundary);
      break;
    
    case kPar_DetectorElementMatch:
      output = getDetectorElementMatch();
      break;

    case kPar_ProcessType:
      output = getProcess().type();
      break;

    case kPar_Pdg:
      output = mCurrentTrack->GetDefinition()->GetPDGEncoding();
      break;

    case kPar_Charge:
      output = mCurrentTrack->GetDefinition()->GetPDGCharge();
      break;
  
    case kPar_TrackId:
      output = mCurrentStep->GetTrack()->GetTrackID();
      break;
      
    // wangzhe
    // Help me track things.
    case kPar_CreatorPdg: {
      G4CompositeTrackInfo* composite=dynamic_cast<G4CompositeTrackInfo*>(mCurrentStep->GetTrack()->GetUserInformation());
      G4HistoryUserTrackInfo* userinfo = composite?dynamic_cast<G4HistoryUserTrackInfo*>
                                                     ( composite->GetHistoryTrackInfo() ):0;
      if(userinfo) 
        output = userinfo->parentPdg();
      break;
      }
      
    case kPar_AncestorPdg: {
      G4CompositeTrackInfo* composite=dynamic_cast<G4CompositeTrackInfo*>(mCurrentStep->GetTrack()->GetUserInformation());
      G4HistoryUserTrackInfo* userinfo = composite?dynamic_cast<G4HistoryUserTrackInfo*>
                                                     ( composite->GetHistoryTrackInfo() ):0;
      if(userinfo)
            output = userinfo->track().track()->particle(); // This will bring me to ancestor track
      break;
    }
    // wz
      
    case kPar_Step_dE:
      output = mCurrentStep->GetTotalEnergyDeposit();
      break;
    case kPar_Step_dE_Ion:
      output = mCurrentStep->GetTotalEnergyDeposit() - mCurrentStep->GetNonIonizingEnergyDeposit();
      break;
    case kPar_Step_qdE:
      output = getQuenchedEnergy();
      break;
    case kPar_Step_dx:
      output = mCurrentStep->GetStepLength();
      break;
    case kPar_Step_dt:
      output = mCurrentStep->GetDeltaTime();
      break;
    case kPar_Step_dAngle:
      output = acos(mCurrentStep->GetPostStepPoint()->GetMomentumDirection().dot(
                     mCurrentStep->GetPreStepPoint()->GetMomentumDirection())
                     ) * CLHEP::radian;
      break;
    case kPar_capTargetZ:
      // Return the neutron capture target Z
      output = m_capinfo->getCapture().GetCapTargetZ();
      //info()<<"stepAction: Z "<<m_capinfo->getCapture().GetCapTargetZ()
       //        <<" A: "<<m_capinfo->getCapture().GetCapTargetA()
        //       <<" time: "<<m_capinfo->getCapture().GetCapTime()
         //      <<" gammaNum: "<<m_capinfo->getCapture().GetCapGammaN()
          //     <<endreq;
      break;
    case kPar_capTargetA:
      // Return the neutron capture target A
      output = m_capinfo->getCapture().GetCapTargetA();
      //info()<<"stepAction: A"<<m_capinfo->getCapture().GetCapTargetA()
       //        <<endreq;               
      break;
      // For energy-weighted quantities, use the average position or time of the step rather than
      // the final position or time of the step.
    case kPar_Step_E_weighted_x:
      output = mCurrentStep->GetPreStepPoint()->GetPosition().x();
      if ( mCurrentStep->GetPostStepPoint() ) output = 0.5*(output + mCurrentStep->GetPostStepPoint()->GetPosition().x());
      output = mCurrentStep->GetTotalEnergyDeposit() * output;
      break;
    case kPar_Step_E_weighted_y:
      output = mCurrentStep->GetPreStepPoint()->GetPosition().y();
      if ( mCurrentStep->GetPostStepPoint() ) output = 0.5*(output + mCurrentStep->GetPostStepPoint()->GetPosition().y());
      output = mCurrentStep->GetTotalEnergyDeposit() * output;
      break;
    case kPar_Step_E_weighted_z:
      output = mCurrentStep->GetPreStepPoint()->GetPosition().z();
      if ( mCurrentStep->GetPostStepPoint() ) output = 0.5*(output + mCurrentStep->GetPostStepPoint()->GetPosition().z());
      output = mCurrentStep->GetTotalEnergyDeposit() * output;
      break;
    case kPar_Step_E_weighted_t:
      output = mCurrentStep->GetTotalEnergyDeposit() * mCurrentStepPoint->GetGlobalTime();
      output = mCurrentStep->GetPreStepPoint()->GetGlobalTime();
      if ( mCurrentStep->GetPostStepPoint() ) output = 0.5*(output + mCurrentStep->GetPostStepPoint()->GetGlobalTime());
      output = mCurrentStep->GetTotalEnergyDeposit() * output;
      break;

    case kPar_Step_qE_weighted_x:
      output = mCurrentStep->GetPreStepPoint()->GetPosition().x();
      if ( mCurrentStep->GetPostStepPoint() ) output = 0.5*(output + mCurrentStep->GetPostStepPoint()->GetPosition().x());
      output = getQuenchedEnergy() * output ; 
      break;
    case kPar_Step_qE_weighted_y:
      output = mCurrentStep->GetPreStepPoint()->GetPosition().y();
      if ( mCurrentStep->GetPostStepPoint() ) output = 0.5*(output + mCurrentStep->GetPostStepPoint()->GetPosition().y());
      output = getQuenchedEnergy() * output ; 
      break;    
    case kPar_Step_qE_weighted_z:
      output = mCurrentStep->GetPreStepPoint()->GetPosition().z();
      if ( mCurrentStep->GetPostStepPoint() ) output = 0.5*(output + mCurrentStep->GetPostStepPoint()->GetPosition().z());
      output = getQuenchedEnergy() * output ; 
      break;
    case kPar_Step_qE_weighted_t:
      output = mCurrentStep->GetPreStepPoint()->GetGlobalTime();
      if ( mCurrentStep->GetPostStepPoint() ) output = 0.5*(output + mCurrentStep->GetPostStepPoint()->GetGlobalTime());
      output = getQuenchedEnergy() * output ; 
      break;

    case kPar_t:
      output = mCurrentStepPoint->GetGlobalTime();
      break;
    case kPar_x:
      output = mCurrentStepPoint->GetPosition().x();
      break;    
    case kPar_y:
      output = mCurrentStepPoint->GetPosition().y();
      break;
    case kPar_z:
      output = mCurrentStepPoint->GetPosition().z();
      break;
    case kPar_r:
      output = sqrt(mCurrentStepPoint->GetPosition().x()*mCurrentStepPoint->GetPosition().x() +
                    mCurrentStepPoint->GetPosition().y()*mCurrentStepPoint->GetPosition().y());
      break;

    case kPar_local_x:
    case kPar_local_y:
    case kPar_local_z:
    case kPar_local_r:
    {
      const G4AffineTransform transformation =
        mCurrentStepPoint->GetTouchable()->GetHistory()->GetTopTransform();
      G4ThreeVector local = transformation.TransformPoint(mCurrentStepPoint->GetPosition());
      
      switch(id){
        case kPar_local_x: output = local.x(); break;
        case kPar_local_y: output = local.y(); break;
        case kPar_local_z: output = local.z(); break;
        case kPar_local_r: output = sqrt(local.x()*local.x() + local.y()*local.y()); break;        
      }
    }
    break;
    
    case kPar_NicheId:
      output = getDetectorId();
      break;

    case kPar_DetectorId:
      output = getDetectorId() & 0xFFFF0000;
      break;
    case kPar_SiteId:
      output = getDetectorId() & 0xFF000000;
      break;
    case kPar_AD:
      output = (getDetectorId() & 0x00FF0000) >> 16;
      break;
    case kPar_Site:
      output = (getDetectorId() & 0xFF000000) >> 24;
      break;

    case kPar_momentum:
       output = mCurrentStepPoint->GetMomentum().mag();
       break;
    case kPar_TotEnergy:
      output = mCurrentStepPoint->GetTotalEnergy();
      break;
    case kPar_KineticEnergy:
      output = mCurrentStepPoint->GetKineticEnergy();
      break;
    case kPar_mass:
      output = mCurrentStep->GetTrack()->GetDynamicParticle()->GetMass();
      break;
    case kPar_vx:
      output = mCurrentStepPoint->GetMomentumDirection().x();
      break;
    case kPar_vy:
      output = mCurrentStepPoint->GetMomentumDirection().y();
      break;
    case kPar_vz:
      output = mCurrentStepPoint->GetMomentumDirection().z();
      break;

    case kPar_local_vx:
    case kPar_local_vy:
    case kPar_local_vz:
    {
      const G4AffineTransform transformationV =
        mCurrentStepPoint->GetTouchable()->GetHistory()->GetTopTransform();
      G4ThreeVector localV = transformationV.TransformAxis(mCurrentStepPoint->GetMomentumDirection());
      
      switch(id){
        case kPar_local_vx: output = localV.x(); break;
        case kPar_local_vy: output = localV.y(); break;
        case kPar_local_vz: output = localV.z(); break;
      }
    }
    break;
        
    default:
      err() << "Unknown real parameter " << id << endreq;
      output = 0;
  }
  verbose() << " id " << id << " output " << output << endreq;

}
void QueriableStepAction::queryParam ( int  id,
std::string &  output 
) const [virtual, inherited]

Reimplemented from RuleParser::Queriable.

Reimplemented in HistorianStepAction.

Definition at line 455 of file QueriableStepAction.cc.

{
  const ILVolume* volume = 0;
  const IDetectorElement* ide = 0;

  switch(id) {
  case kPar_CreatorProcessName:
    if(mCurrentTrack->GetParentID()==0) // This is a primary.
        mProcess = SimProcess(SimProcess::kPrimaryVertex,"");
    else if(mCurrentTrack->GetCreatorProcess())
        mProcess = SimProcess(SimProcess::kParticleStart,mCurrentTrack->GetCreatorProcess()->GetProcessName());
    else
        mProcess = SimProcess(SimProcess::kParticleStart,""); 
    output = mProcess.name();
    verbose()<<" id " << id << " creator process name "<<output<<endreq;
    return;
    
  case kPar_ProcessName:
    switch(mCurrentStep->GetPostStepPoint()->GetStepStatus()) {
        case fWorldBoundary:   mProcess = SimProcess(SimProcess::kWorldBoundary,""); break;
        case fGeomBoundary:    mProcess = SimProcess(SimProcess::kGeomBoundary,""); break;
        default:
        SimProcessFromG4Process(mCurrentStep->GetPostStepPoint()->GetProcessDefinedStep(),mProcess);
    }
    output = mProcess.name();
    verbose()<<" id " << id << " process name "<<output<<endreq;
    return;

  case kPar_DetectorElementName:
    ide = getDetectorElement();
    if(ide)  output = ide->name();
    else output = "";
    verbose() << " id " << id << " detector element output " << output << endreq;
    return;

  case kPar_LogicalVolumeName:
    volume = getLogicalVolume();
    if(volume)  output = volume->name();
    else output = "";
    verbose() << " id " << id << " logical volume name output " << output << endreq;
    return;
    
  case kPar_MaterialName:
    volume = getLogicalVolume();
    if(volume)output = volume->materialName();
    else output ="";
    verbose() << " id " << id << " material name output " << output << endreq;
    return;
    
  case kPar_ParticleName:
    output = mCurrentStep->GetTrack()->GetDefinition()->GetParticleName();
    verbose() << " id " << id << " particle name output " << output << endreq;
    return;
    
  default:
    err() << "Unknown string parameter " << id << endreq;
    output = "";
  }
}
const IDetectorElement * QueriableStepAction::getDetectorElement ( ) const [inherited]

Definition at line 520 of file QueriableStepAction.cc.

{
  if(!mDetElement) {
    // wangzhe:
    // (GetPhysicalVolume() == 0) at the boundary of the world.
    // (MotherLogical pointer == 0) means that current step point is
    // in the outmost volume, "Universe". 
    // Currently no touchable history info is available here, so neglect
    // any query in this volume.
    if(mCurrentStep->GetPreStepPoint()->GetPhysicalVolume() == 0) {
      return 0;
    }
    if(mCurrentStep->GetPreStepPoint()->GetPhysicalVolume()->GetMotherLogical() == 0) {
      return 0;
    }
    // wz
    G4VTouchable* touch = mCurrentStep->GetPreStepPoint()->GetTouchableHandle()();
    G4TouchableHistory* hist = dynamic_cast<G4TouchableHistory*>(touch);
    assert(hist);

    if (!hist->GetHistoryDepth()) {
        err() << "getDetectorElement(): given empty touchable history" << endreq;
        return 0;
    }

    StatusCode sc = mTouchToDetElem->GetBestDetectorElement(hist,
                                                            m_DetectorElementSearchPath,
                                                            mDetElement,
                                                            mDetElementMatch);
    if (sc.isFailure()) {
        err() << "Failed to get best detector element" << endreq;
        return 0;
    }
  }
  return mDetElement;  
}
int QueriableStepAction::getDetectorElementMatch ( ) const [inherited]

Definition at line 557 of file QueriableStepAction.cc.

const ILVolume * QueriableStepAction::getLogicalVolume ( ) const [inherited]

Definition at line 563 of file QueriableStepAction.cc.

{
  if(!mLogicVolume) {
    const IPVolume* pvol = getPhysicalVolume();
    if(pvol) mLogicVolume = pvol->lvolume();
  }
  return mLogicVolume; 
}
const IPVolume * QueriableStepAction::getPhysicalVolume ( ) const [inherited]

Definition at line 572 of file QueriableStepAction.cc.

{
  if(!mPhysVolume) {
     mTouchToDetElem->G4VolumeToDetDesc(mCurrentStep->GetPreStepPoint()->GetPhysicalVolume(),
                                        mPhysVolume
                                        );
   }
   return mPhysVolume;  
}
const SimProcess & QueriableStepAction::getProcess ( ) const [inherited]

Definition at line 582 of file QueriableStepAction.cc.

{
  if(!mProcess.isValid()){
    if(mCurrentStepPoint == mCurrentStep->GetPreStepPoint()) {
      if(mCurrentTrack->GetParentID()==0) // This is a primary.
         mProcess = SimProcess(SimProcess::kPrimaryVertex,"");
      else if(mCurrentTrack->GetCreatorProcess())
        mProcess = SimProcess(SimProcess::kParticleStart,mCurrentTrack->GetCreatorProcess()->GetProcessName());
      else
        mProcess = SimProcess(SimProcess::kParticleStart,"");
    } else {
      switch(mCurrentStep->GetPostStepPoint()->GetStepStatus()) {
        case fWorldBoundary:   mProcess = SimProcess(SimProcess::kWorldBoundary,""); break;
        case fGeomBoundary:    mProcess = SimProcess(SimProcess::kGeomBoundary,""); break;
        default:
        SimProcessFromG4Process(mCurrentStepPoint->GetProcessDefinedStep(),mProcess);
      }
    }
  }
  return mProcess;
}
double QueriableStepAction::getQuenchedEnergy ( ) const [inherited]

Definition at line 604 of file QueriableStepAction.cc.

{
  if(mQuenchedEnergy>=0) return mQuenchedEnergy;
  G4Track *aTrack = mCurrentStep->GetTrack();
  G4ParticleDefinition* aParticle = aTrack->GetDefinition();
  G4String particleName = aParticle->GetParticleName();
  //info()<<"particle name "<<particleName<<endreq;

  double dE = mCurrentStep->GetTotalEnergyDeposit();
  double dx = mCurrentStep->GetStepLength();

  // Find quenched energy deposit.
  mQuenchedEnergy = 0;
  if(dE > 0) {
    if(aParticle == G4Gamma::Gamma()) // It is a gamma
    {
      G4LossTableManager* manager = G4LossTableManager::Instance();
      dx = manager->GetRange(G4Electron::Electron(), dE, aTrack->GetMaterialCutsCouple());
      //info()<<"dE_dx = "<<dE/dx/(MeV/mm)<<"MeV/mm"<<endreq;
    } 
    G4Material* aMaterial = mCurrentStep->GetPreStepPoint()->GetMaterial();
    G4MaterialPropertiesTable* aMaterialPropertiesTable =
        aMaterial->GetMaterialPropertiesTable();
    if (aMaterialPropertiesTable) {

      // There are some properties. Is there a scintillator property?
      const G4MaterialPropertyVector* Fast_Intensity = 
          aMaterialPropertiesTable->GetProperty("FASTCOMPONENT"); 
      const G4MaterialPropertyVector* Slow_Intensity =
          aMaterialPropertiesTable->GetProperty("SLOWCOMPONENT");
  
      if (Fast_Intensity || Slow_Intensity ) {
        // It's a scintillator.
        double delta = dE/dx/aMaterial->GetDensity(); 
        //double birk1 = 0.0125*g/cm2/MeV;
        double birk1 = m_BirksConstant1; 
        if(mCurrentTrack->GetDefinition()->GetPDGCharge()>1.1)//for particle charge greater than 1.
           birk1 = 0.57*birk1;
        //double birk2 = (0.0031*g/MeV/cm2)*(0.0031*g/MeV/cm2);
        double birk2 = m_BirksConstant2;
        mQuenchedEnergy= dE /(1+birk1*delta+birk2*delta*delta);
      }
    }
  }
  return mQuenchedEnergy;  
}
const G4HistoryUserTrackInfo * QueriableStepAction::getUserTrackInfo ( ) const [inherited]

Definition at line 651 of file QueriableStepAction.cc.

{
  if(!mUserTrackInfo) {
    G4CompositeTrackInfo* composite=dynamic_cast<G4CompositeTrackInfo*>(mCurrentTrack->GetUserInformation());
    mUserTrackInfo = composite?dynamic_cast<G4HistoryUserTrackInfo*>
                                                 ( composite->GetHistoryTrackInfo() ):0;
    // mUserTrackInfo = dynamic_cast<G4HistoryUserTrackInfo*>(mCurrentTrack->GetUserInformation());
  }
  return mUserTrackInfo;  
}
unsigned int QueriableStepAction::getDetectorId ( ) const [inherited]

Definition at line 662 of file QueriableStepAction.cc.

{
  if(mDetectorId==0xFFFFFFF) {
    const IDetectorElement* de = getDetectorElement();
    mDetectorId = getDetectorId(de);
  }
  return mDetectorId;
}
unsigned int QueriableStepAction::getDetectorId ( const IDetectorElement *  de) const [inherited]

Definition at line 671 of file QueriableStepAction.cc.

{
  DetectorIdCache_t::iterator it = mDetectorIdCache.find(de);
  if(it!=mDetectorIdCache.end()) {
    return it->second;
  }

  const ParamValidDataObject* params = de->params();
  
  // Check each type of user param that the user could have provided:
  for(unsigned int ip = 0; ip<m_IdParameterNames.size(); ip++) {
    if (params->exists(m_IdParameterNames[ip])) {
      if (unsigned int id = (unsigned int)(params->param<int>(m_IdParameterNames[ip]))) {
        // Got it!
        mDetectorIdCache[de] = id;  // cache it for next time
        //info() << "Found id for " << getDetectorElement()->name() << " --> id " << std::hex << mDetectorId << std::dec << endreq;
        return id;
      }
    }
  }
  
  // Ok, we didnt' find it.
  // So, move up the support tree recursively.

  // --this does't work because geo->detElem() is private
  // const IGeometryInfo* igeo = de->geometry();
  //   ILVolume::ReplicaPath dummy;
  //   IGeometryInfo* nextGeo = 0;
  //   igeo->location(nextGeo, dummy); 
  //   if(nextGeo==0) return 0x0; // We failed to find any kind of ID.
  // const GeometryInfoPlus* geoplus = dynamic_cast<GeometryInfoPlus*>(nextGeo);
  //   assert(geoplus); // What's up with this?
  //   return getDetectorId(geoplus->detElem());

  // -- this works if the child elements are all set correctly in the XML.
  const IDetectorElement* nextElem = de->parentIDetectorElement();
  if(nextElem==0) {
    const IDetectorElement* ide = getDetectorElement();
    std::string DeName;
    if(ide) DeName = ide->name();
    else DeName = "";
    warning() << "Found no next-element for " << DeName << " --> id " << std::hex << mDetectorId << std::dec << endreq;
    mDetectorIdCache[de] = 0;
    return 0;
  }
  unsigned int id = getDetectorId(nextElem);
  mDetectorIdCache[de] = id;
  return id;
}

Member Data Documentation

std::vector< std::vector<std::string> > UnObserverStepAction::m_Stats [private]

Definition at line 25 of file UnObserverStepAction.h.

std::vector<StatGroup> UnObserverStepAction::mGroups [private]

Definition at line 36 of file UnObserverStepAction.h.

Definition at line 39 of file UnObserverStepAction.h.

Definition at line 46 of file UnObserverStepAction.h.

std::vector<std::string> QueriableStepAction::m_DetectorElementSearchPath [protected, inherited]

Definition at line 68 of file QueriableStepAction.h.

std::string QueriableStepAction::m_TouchableToDetelem_name [protected, inherited]

Definition at line 69 of file QueriableStepAction.h.

std::vector<std::string> QueriableStepAction::m_IdParameterNames [protected, inherited]

Definition at line 70 of file QueriableStepAction.h.

const G4Step* QueriableStepAction::mCurrentStep [protected, inherited]

Definition at line 73 of file QueriableStepAction.h.

const G4Track* QueriableStepAction::mCurrentTrack [protected, inherited]

Definition at line 74 of file QueriableStepAction.h.

const G4StepPoint* QueriableStepAction::mCurrentStepPoint [protected, inherited]

Definition at line 75 of file QueriableStepAction.h.

Definition at line 76 of file QueriableStepAction.h.

const IDetectorElement* QueriableStepAction::mDetElement [mutable, protected, inherited]

Definition at line 80 of file QueriableStepAction.h.

int QueriableStepAction::mDetElementMatch [mutable, protected, inherited]

Definition at line 81 of file QueriableStepAction.h.

const ILVolume* QueriableStepAction::mLogicVolume [mutable, protected, inherited]

Definition at line 82 of file QueriableStepAction.h.

const IPVolume* QueriableStepAction::mPhysVolume [mutable, protected, inherited]

Definition at line 83 of file QueriableStepAction.h.

DayaBay::SimProcess QueriableStepAction::mProcess [mutable, protected, inherited]

Definition at line 84 of file QueriableStepAction.h.

double QueriableStepAction::mQuenchedEnergy [mutable, protected, inherited]

Definition at line 85 of file QueriableStepAction.h.

Definition at line 86 of file QueriableStepAction.h.

unsigned int QueriableStepAction::mDetectorId [mutable, protected, inherited]

Definition at line 87 of file QueriableStepAction.h.

Definition at line 89 of file QueriableStepAction.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:33 for Historian by doxygen 1.7.4