/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 | Protected Types | Protected Attributes | Private Types | Private Member Functions | Private Attributes
HistorianStepAction Class Reference

#include <HistorianStepAction.h>

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

List of all members.

Public Member Functions

 HistorianStepAction (const std::string &type, const std::string &name, const IInterface *parent)
virtual ~HistorianStepAction ()
virtual StatusCode initialize ()
virtual StatusCode finalize ()
virtual void UserSteppingAction (const G4Step *)
virtual void queryParam (int id, double &output) const
virtual void queryParam (int id, std::string &output) const
DayaBay::SimTrackReference getAncestorTrack () const
DayaBay::SimVertexReference getAncestorVertex () const
virtual bool IsInterestingTrack (const G4Track *)
virtual bool IsInterestingVertex (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:
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 Types

enum  EHistorianStepParams {
  kPar_HISTORIAN_START = kPar_END_QUERIABLE, kPar_ParentPdg, kPar_ParentIndirection, kPar_GrandParentPdg,
  kPar_GrandParentIndirection, kPar_distanceFromLastVertex, kPar_timeSinceLastVertex, kPar_energyLossSinceLastVertex,
  kPar_angleFromLastVertex
}

Private Member Functions

StatusCode CreateTrack (const G4Track *, DayaBay::SimTrack *&track)
StatusCode CreateVertex (DayaBay::SimVertex *&v, const DayaBay::SimTrackReference &parent)
 HistorianStepAction ()
 no default constructor
 HistorianStepAction (const HistorianStepAction &)
 no copy
HistorianStepActionoperator= (const HistorianStepAction &)
 no =

Private Attributes

G4bool idflag
std::string m_TrackSelection
std::string m_VertexSelection
G4bool m_UseFastMuEnergyCut
IHistoryKeepermHistoryKeeper
RuleParser::RulemTrackRule
RuleParser::RulemVertexRule
DayaBay::SimParticleHistorymCurrentHistory
DayaBay::SimTrackReference mCurrentTrackRef

Detailed Description

Definition at line 30 of file HistorianStepAction.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

Enumerator:
kPar_HISTORIAN_START 
kPar_ParentPdg 
kPar_ParentIndirection 
kPar_GrandParentPdg 
kPar_GrandParentIndirection 
kPar_distanceFromLastVertex 
kPar_timeSinceLastVertex 
kPar_energyLossSinceLastVertex 
kPar_angleFromLastVertex 

Definition at line 81 of file HistorianStepAction.h.

                            {
    kPar_HISTORIAN_START = kPar_END_QUERIABLE
    // Tracky things specific to Histories:
    , kPar_ParentPdg                  // int
    , kPar_ParentIndirection          // int
    , kPar_GrandParentPdg             // int
    , kPar_GrandParentIndirection     // int
  
     // Vertex-y things only applicable if we have current history information.
    , kPar_distanceFromLastVertex     // double
    , kPar_timeSinceLastVertex        // double
    , kPar_energyLossSinceLastVertex  // double
    , kPar_angleFromLastVertex        // double
    
  };
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

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

Definition at line 36 of file HistorianStepAction.cc.

  : QueriableStepAction( type, name, parent )
  , mHistoryKeeper(0)
  , mTrackRule(0)
  , mVertexRule(0)
{ 
  declareProperty("TrackSelection",m_TrackSelection = "none");
  declareProperty("VertexSelection",m_VertexSelection = "none");
  declareProperty("UseFastMuEnergyCut",m_UseFastMuEnergyCut = false);
}
virtual HistorianStepAction::~HistorianStepAction ( ) [inline, virtual]

Definition at line 37 of file HistorianStepAction.h.

{};
HistorianStepAction::HistorianStepAction ( ) [private]

no default constructor

HistorianStepAction::HistorianStepAction ( const HistorianStepAction ) [private]

no copy


Member Function Documentation

StatusCode HistorianStepAction::initialize ( ) [virtual]

Reimplemented from QueriableStepAction.

Definition at line 50 of file HistorianStepAction.cc.

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

  StatusCode sc = QueriableStepAction::initialize();
  if (sc.isFailure()) return sc;
  
  // Set up default search path for detector elements. Note that this throws away PMT-specific information.
  // Also, change the default search path:
  m_DetectorElementSearchPath.clear();
  m_DetectorElementSearchPath.push_back("/dd/Structure/DayaBay");
  m_DetectorElementSearchPath.push_back("/dd/Structure/Sites");
  m_DetectorElementSearchPath.push_back("/dd/Structure/Pool");
  m_DetectorElementSearchPath.push_back("/dd/Structure/AD");
  m_DetectorElementSearchPath.push_back("/dd/Structure/RPC");
  
  sc = service("HistoryKeeper",mHistoryKeeper,true);
  if (sc.isFailure()) return sc;
  

  // Build the selection rules.
  
  using namespace RuleParser;
  using std::string;
    
  ParameterList tplist;
  QueriableStepAction::GetTrackParameterList(tplist);
  tplist.add<double>(kPar_ParentPdg,"ParentPdg","AncestorPdg","Ancestor");
  tplist.add<double>(kPar_ParentIndirection,"ParentIndirection","AncestorIndirection");
  tplist.add<double>(kPar_GrandParentPdg,"GrandParentPdg","GrandParent");
  tplist.add<double>(kPar_GrandParentIndirection,"GrandParentIndirection");


  ParameterList vplist;
  QueriableStepAction::GetVertexParameterList(vplist);
  vplist.add<double>(kPar_ParentPdg,"ParentPdg","AncestorPdg","Ancestor");
  vplist.add<double>(kPar_ParentIndirection,"ParentIndirection","AncestorIndirection");
  vplist.add<double>(kPar_GrandParentPdg,"GrandParentPdg","GrandParent");
  vplist.add<double>(kPar_GrandParentIndirection,"GrandParentIndirection");
  
  vplist.add<double>(kPar_distanceFromLastVertex,"distanceFromLastVertex");
  vplist.add<double>(kPar_timeSinceLastVertex,"TimeSinceLastVertex");
  vplist.add<double>(kPar_energyLossSinceLastVertex,"EnergyLostSinceLastVertex");
  vplist.add<double>(kPar_angleFromLastVertex,"AngleFromLastVertex");
                                                     
                                                     
  bool result;
  result = RuleParser::CreateRules(m_TrackSelection,tplist,mTrackRule);
  
  if(mTrackRule==0 || result==false) {
    err() << "Failed to interpret selection string: " << m_TrackSelection << endreq;
    return StatusCode::FAILURE;
  } else {
    info() << "Configured with track selection:  " << mTrackRule->name() << endreq;
  }
  

  result = RuleParser::CreateRules(m_VertexSelection,vplist,mVertexRule);
  
  if(mVertexRule==0 || result==false) {
    err() << "Failed to interpret selection string: " << m_VertexSelection << endreq;
    return StatusCode::FAILURE;
  } else {
    info() << "Configured with vertex selection: " << mVertexRule->name() << endreq;
  }
  

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

Reimplemented from QueriableStepAction.

Definition at line 123 of file HistorianStepAction.cc.

void HistorianStepAction::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.

wangzhe The following code is for trac #748 waiting for a long lifetime particle, like > 1 year. They will temporarily stay here until we identify them.

wz

Implements QueriableStepAction.

Definition at line 376 of file HistorianStepAction.cc.

{
  
  // Reset temporary things:
  mHistoryKeeper->GetCurrentHistory(mCurrentHistory);
  assert(mCurrentHistory);
 
  //Add the fast muon simulation energy cut for different particle, intend to increase the simulation speed. By Haoqi Jan,27, 2011.
  
  G4Track* g4track = g4step->GetTrack();
  if( g4track->GetGlobalTime()> 3e16 )  {  //3e16 /* a year */
      G4CompositeTrackInfo* composite=dynamic_cast<G4CompositeTrackInfo*>(g4track->GetUserInformation());
      G4HistoryUserTrackInfo* userinfo = composite?dynamic_cast<G4HistoryUserTrackInfo*>
                                                   ( composite->GetHistoryTrackInfo() ):0;
      int parentPdg = 0;
      if(userinfo) {
        parentPdg = userinfo->parentPdg();
      }
      warning()<<"A 1 year old track orignated by particle "<<parentPdg<<endreq;
  }

  G4int trackId = g4track ->GetTrackID();
  G4StepPoint* thePrePoint;
  G4VPhysicalVolume* physV;
  G4Material* material;
  G4ParticleDefinition* particle = g4track->GetDefinition();
  G4String pname = particle->GetParticleName();
  G4ProcessVector* processVector; 
  thePrePoint = g4step->GetPreStepPoint();
  physV = thePrePoint->GetPhysicalVolume();
  material = physV->GetLogicalVolume()->GetMaterial(); 
  G4double EnergyCut;
 
  idflag = false ; // variables must be initialized before use
  if(g4track->GetCurrentStepNumber()==1&&m_UseFastMuEnergyCut==true)
    { idflag=true;}
 //G4cout<<"FastMuCut ="<<m_UseFastMuEnergyCut <<G4endl;
  if(idflag==true)
    {
       
      if(material->GetName()=="/dd/Materials/Rock"||material->GetName()=="/dd/Materials/Air")
      {  
       //  G4cout<<"Rock or air aterial Name ="<<material->GetName()<<G4endl;
         if(particle->GetPDGCharge()!=0)
         {
            EnergyCut=500.0;
         }
         else
         {
            EnergyCut=30.0; 
         }
      }//rock
      else if(material->GetName()=="/dd/Materials/OwsWater"||material->GetName()=="/dd/Materials/IwsWater")
      {
         // G4cout<<"Ows or IWs aterial Name ="<<material->GetName()<<G4endl;
         if(particle->GetPDGCharge()!=0)
         {
            EnergyCut=50.0;
         }
         else
         {
            EnergyCut=3.0;   
          }
       }
       else 
       {
        EnergyCut=0.;
     //    if(material->GetName()=="/dd/Materials/LiquidScintillator"||material->GetName()=="/dd/Materials/GdDopedLS") G4cout<<"Material ="<<material->GetName()<<G4endl;

        }//material
      if(g4track->GetKineticEnergy()/MeV<EnergyCut)
      {
          
          //aTrack->SetKineticEnergy(0.00001/MeV);
          g4track->SetTrackStatus(fStopAndKill);
          idflag=false;
      }
   }//idflag
  
  // First, do things related to SimTracks.

  if(g4track->GetCurrentStepNumber() == 1) {
    // This is the first time we've seen this track. 
    // The interesting point is the first one, the track initial vertex.
    QueriableStepAction::Reset(g4step,g4track,g4step->GetPreStepPoint());

    // Does a SimTrack already exist for this G4track?
    int id = g4track->GetTrackID();
    SimTrackReference ref = mCurrentHistory->track(id);
    SimTrack* track = ref.track();

    if(track==0) { // There is no track or reference.
      // Maybe we need to make one. Is it a primary track?
      if(g4track->GetParentID()==0) {
        // This is a primary track, so this action is responsible for starting it.
        // Create the SimTrack
        CreateTrack(g4track,track);
        ref = SimTrackReference(track,0);
      } else {
        if(IsInterestingTrack(g4track)) {
          CreateTrack(g4track,track);        
          ref = SimTrackReference(track,0);
        }
      }
    }

    if(ref.track()==0) {
      // OK, there is still no valid reference. Let's build an indirect one.
      SimTrackReference  parentTrack = getAncestorTrack();
      SimVertexReference parentVertex = getAncestorVertex();

       // Increment indirection.
       ref = SimTrackReference(parentTrack.track(),parentTrack.indirection()+1);
       mCurrentHistory->addTrackRef(g4track->GetTrackID(),ref);
     }


    // currentTrack should now be set correctly, one way or another
    if(! ref.isDirect()) {
      // The reference is indirect, meaning we aren't recording the current one.
      // So, instead, we want to store some 'missing' information.
      ref.track()->incrementUnrecordedDescendants(pdgFromG4Track(g4track));
    }

    // The current track reference, for convenience.
    //mCurrentHistory->setCurrentTrack(ref);
  }
  
  // Second, do things related to SimVertex.
  
  // The point of interest is the post-step if we're looking at the step/vertex
  QueriableStepAction::Reset(g4step,g4track,g4step->GetPostStepPoint());
  
  mCurrentTrackRef = mCurrentHistory->track(g4track->GetTrackID());
  if(!mCurrentTrackRef.track()) {
    // This is an error and shouldn't happen.
    err() << "Somehow we got a bad track reference on track id " 
          << g4track->GetTrackID() << endreq;
    err() << " This shouldn't happen - everything should have a parent." << endreq;
    assert(0);
  }
  

  // Is the current track getting recorded?
  if(mCurrentTrackRef.isDirect()){
    // Yes. Is this interaction interesting on it's own?
    if(IsInterestingVertex(g4step)) {
      verbose() << " InterestingVertex. Track in " << g4step->GetTrack()->GetVolume()->GetName() 
              << " particle " << g4step->GetTrack()->GetDefinition()->GetParticleName()
              << " G4TrackID "<< g4step->GetTrack()->GetTrackID()
              << " current step# " << g4step->GetTrack()->GetCurrentStepNumber()
              << " step length(mm) " << g4step->GetTrack()->GetStepLength()/mm
              << endreq;
      SimVertex* vertex = 0;
      CreateVertex(vertex,mCurrentTrackRef);
      // Give ownership to History.
      mCurrentHistory->addVertex(vertex); 
      // We have a recorded track and recorded vertex. So, add the vertex.
      mCurrentTrackRef.track()->addVertex(vertex);
    }
  }

  // Get current PDG code.
  int pdg = pdgFromG4Track(g4track);
  
  // Get the current best reference to a vertex
  SimVertexReference curVertexRef(mCurrentTrackRef.track()->vertices().back(),mCurrentTrackRef.indirection());

  assert(fpSteppingManager);
  assert(fpSteppingManager->GetStep() == g4step);
  
  // Loop over secondaries, push information forward.
  int nSecAtRest = fpSteppingManager->GetfN2ndariesAtRestDoIt();
  int nSecAlong  = fpSteppingManager->GetfN2ndariesAlongStepDoIt();
  int nSecPost   = fpSteppingManager->GetfN2ndariesPostStepDoIt();
  int nSecTotal  = nSecAtRest+nSecAlong+nSecPost;
  G4TrackVector* secVec = fpSteppingManager->GetfSecondary();

  size_t end = (*secVec).size();
  size_t start = end-nSecTotal;
  for(size_t itrack=start; itrack < end; itrack++)
  {
    G4Track* g4track = (*secVec)[itrack];
    G4HistoryUserTrackInfo* user
      = new G4HistoryUserTrackInfo(curVertexRef,mCurrentTrackRef,pdg);
    G4VUserTrackInformation* composite=g4track->GetUserInformation();
    if( composite ) {
        G4CompositeTrackInfo* compositei=dynamic_cast<G4CompositeTrackInfo*>(composite);
        if( !compositei || compositei->GetHistoryTrackInfo() ) {
          warning() << "Already found user info!  itrack=" << itrack
            << " out of " << (*secVec).size() << endreq;
        }
        if (compositei) compositei->SetHistoryTrackInfo(user);
    } else {
        G4CompositeTrackInfo* comp=new G4CompositeTrackInfo();
        comp->SetHistoryTrackInfo(user);
        g4track->SetUserInformation(comp);
    }
  }

  // G4CompositeTrackInfo* composite=dynamic_cast<G4CompositeTrackInfo*>(mCurrentStep->GetTrack()->GetUserInformation());
  // G4HistoryUserTrackInfo* userinfo = composite?dynamic_cast<G4HistoryUserTrackInfo*>
  //                                                ( composite->GetHistoryTrackInfo() ):0;
  // printf("comp %p user %p parent %i\n", composite, userinfo, userinfo?userinfo->parentPdg():0);
}
void HistorianStepAction::queryParam ( int  id,
double &  output 
) const [virtual]

Reimplemented from QueriableStepAction.

Definition at line 132 of file HistorianStepAction.cc.

{
  SimTrackReference ref;
  G4StepPoint* postPoint = mCurrentStep->GetPostStepPoint();
  const SimVertex* v = 0;
  Gaudi::XYZPoint pos;
  Gaudi::XYZVector dPos;
  
  switch(id) {
    case kPar_ParentPdg:
      ref = getAncestorTrack();
      if(ref.track()) output = ref.track()->particle();
      else output = 0;
      break;
      
    case kPar_ParentIndirection:
      output = getAncestorTrack().indirection();
      break;

    case kPar_GrandParentPdg:
      output = 0;
      ref = getAncestorTrack();
      if(ref.track()) ref = ref.track()->ancestorTrack();
      if(ref.track()) output = ref.track()->particle();
      break;
    
    case kPar_GrandParentIndirection:
      output = -999;
      ref = getAncestorTrack();
      if(ref.track()) output = ref.track()->ancestorTrack().indirection();
      break;

    // Special history stuff.
    case kPar_distanceFromLastVertex:      
    case kPar_timeSinceLastVertex:
    case kPar_energyLossSinceLastVertex:
    case kPar_angleFromLastVertex:
      assert(mCurrentTrackRef.track());
      assert(mCurrentTrackRef.track()->vertices().size());
      v = mCurrentTrackRef.track()->vertices().back();

      switch(id) {
        case kPar_distanceFromLastVertex:
          pos = Gaudi::XYZPoint(postPoint->GetPosition().x()
                                ,postPoint->GetPosition().y()
                                ,postPoint->GetPosition().z());
          dPos = pos - v->position();
          output = dPos.R();
          return;
        case kPar_timeSinceLastVertex:
          output = postPoint->GetGlobalTime() - v->time();
          return;
        case kPar_energyLossSinceLastVertex:
          output = postPoint->GetTotalEnergy() - v->totalEnergy();
          return;
        case kPar_angleFromLastVertex:
          output = ( postPoint->GetMomentumDirection().x() * v->momentum().x()
                   + postPoint->GetMomentumDirection().y() * v->momentum().y()
                   + postPoint->GetMomentumDirection().z() * v->momentum().z() )
                   / (v->momentum().R());
          output = acos(output)*CLHEP::radian;
          return;
      }
      break;

    default:
      QueriableStepAction::queryParam(id,output);
  }
  verbose() << " id " << id << " output " << output << endreq;

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

Reimplemented from QueriableStepAction.

Definition at line 204 of file HistorianStepAction.cc.

{
  QueriableStepAction::queryParam(id,output);
  verbose() << " id " << id << " output " << output << endreq;

}
SimTrackReference HistorianStepAction::getAncestorTrack ( ) const

Definition at line 212 of file HistorianStepAction.cc.

{
  const G4HistoryUserTrackInfo* info = getUserTrackInfo();
  if(!info) {
    err() << "Found track id " << mCurrentTrack->GetTrackID() << " with no user info. Shouldn't happen." << endreq;
    err() << "(Are you running with both the HistorianTrackAction AND the HistorianStepAction?)" << endreq;
    assert(info);    
  }
  return info->track();
}
SimVertexReference HistorianStepAction::getAncestorVertex ( ) const

Definition at line 223 of file HistorianStepAction.cc.

{
  const G4HistoryUserTrackInfo* info = getUserTrackInfo();
  if(!info) {
    err() << "Found track id " << mCurrentTrack->GetTrackID() << " with no user info. Shouldn't happen." << endreq;
    err() << "(Are you running with both the HistorianTrackAction AND the HistorianStepAction?)" << endreq;
    assert(info);    
  }
  
  return info->vertex();
}
bool HistorianStepAction::IsInterestingTrack ( const G4Track *  ) [virtual]

Definition at line 369 of file HistorianStepAction.cc.

{
  // All the magic is done by the RuleParser system, and the Rules it has procuded:
  return mTrackRule->select(this);
}
bool HistorianStepAction::IsInterestingVertex ( const G4Step *  ) [virtual]

Definition at line 363 of file HistorianStepAction.cc.

{
  // All the magic is done by the RuleParser system, and the Rules it has procuded:
  return mVertexRule->select(this);
}
StatusCode HistorianStepAction::CreateTrack ( const G4Track *  g4track,
DayaBay::SimTrack *&  track 
) [private]

Definition at line 268 of file HistorianStepAction.cc.

{
  SimTrackReference parentTrack;
  SimVertexReference parentVertex;
  int parentPdg = 0;
 
  const HepMC::GenParticle*  genPart = 0;

  if(g4track->GetParentID()==0){
    //info() << "Building primary particle" << endreq;
    // Primary track is a special case.
    parentTrack = SimTrackReference(0,0);
    parentVertex= SimVertexReference(0,0);

    const G4DhPrimaryParticleInformation* primaryinfo = 
        dynamic_cast<const G4DhPrimaryParticleInformation*>(
             g4track->GetDynamicParticle()->GetPrimaryParticle()->GetUserInformation() );
     if(primaryinfo) {
       genPart  = primaryinfo->GetHepParticle();      
     }

  } else {
    const G4HistoryUserTrackInfo* userinfo = getUserTrackInfo();
    if(userinfo) {
      parentTrack = userinfo->track();
      parentVertex = userinfo->vertex();
      parentPdg = userinfo->parentPdg();
    } else {
      err() << "Found track id " << mCurrentTrack->GetTrackID() << " with no user info. Shouldn't happen." << endreq;
      err() << "(Are you running with both the HistorianTrackAction AND the HistorianStepAction?)" << endreq;
      assert(userinfo);    
    }
    genPart = parentTrack.track()->primaryParticle();
  }

  int pdg = pdgFromG4Track(g4track);
  if(pdg == 0 ) {
    warning() << "Found particle with zero pdg code. Geant claims this particle is: " 
              << g4track->GetDefinition()->GetParticleName() << endreq;
  }
 
  track = new SimTrack();
  track->setTrackId(g4track->GetTrackID());
  track->setParentParticle(parentPdg);
  track->setAncestorTrack(parentTrack);
  track->setAncestorVertex(parentVertex);
  track->setPrimaryParticle(genPart);
  track->setParticle(pdg);
    
  // Create the start vertex.
  const G4ThreeVector& pos = mCurrentStepPoint->GetPosition();
  const G4ThreeVector mom = mCurrentStepPoint->GetMomentum();

  // Figure out the process.
  SimProcess proc;
  if(g4track->GetParentID()==0){
    // This is a primary.
    proc = SimProcess(SimProcess::kPrimaryVertex,"");
  } else {
    if(mCurrentTrack->GetCreatorProcess())
      proc = SimProcess(SimProcess::kParticleStart,mCurrentTrack->GetCreatorProcess()->GetProcessName());
    else
      proc = SimProcess(SimProcess::kParticleStart,"");
  }

  SimVertex* vertex = new SimVertex(
                       SimTrackReference(track,0),
                       proc,
                       mCurrentStepPoint->GetGlobalTime(),
                       Gaudi::XYZPoint(pos.x(),pos.y(),pos.z()),
                       mCurrentStepPoint->GetTotalEnergy(),
                       Gaudi::XYZVector(mom.x(),mom.y(),mom.z())
                      );
    
  track->addVertex(vertex);
  mCurrentHistory->addTrack(track);
  mCurrentHistory->addTrackRef(g4track->GetTrackID(),SimTrackReference(track,0));
  mCurrentHistory->addVertex(vertex);
  if(g4track->GetParentID()==0){ 
    // This is a primary.
    mCurrentHistory->addPrimaryTrack(track);
  } else {
    // This track is not a primary, so try to link from our parent vertex to us.
    // Invert the parent-vertex indirection to find out the parentVertex->this track reference.
    SimTrackReference downref(track,parentVertex.indirection());
    if(!parentVertex.vertex())
      err() << "SimTrack being created, but no Parent vertex is stored. This shouldn't happen." << endreq;
    else
      parentVertex.vertex()->addSecondary(downref);
  }  
  return StatusCode::SUCCESS;
}
StatusCode HistorianStepAction::CreateVertex ( DayaBay::SimVertex *&  v,
const DayaBay::SimTrackReference parent 
) [private]

Fill a DayaBay::SimVertex object with data from a G4Step. Only does something if vertex hasn't yet been filled. Note that a call here doesn't mean a vertex is saved; IsInteresting is allowed to call this, in case it wants to know something.

Definition at line 239 of file HistorianStepAction.cc.

{ 
     
  const G4ThreeVector& pos = mCurrentStepPoint->GetPosition();
  const G4ThreeVector mom = mCurrentStepPoint->GetMomentum();

  vertex = new SimVertex(
                      parent,
                      getProcess(),
                      mCurrentStepPoint->GetGlobalTime(),
                      Gaudi::XYZPoint(pos.x(),pos.y(),pos.z()),
                      mCurrentStepPoint->GetTotalEnergy(),
                      Gaudi::XYZVector(mom.x(),mom.y(),mom.z())
                     );

  return StatusCode::SUCCESS;
}
HistorianStepAction& HistorianStepAction::operator= ( const HistorianStepAction ) [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;
}
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

G4bool HistorianStepAction::idflag [private]

Definition at line 57 of file HistorianStepAction.h.

Definition at line 60 of file HistorianStepAction.h.

Definition at line 61 of file HistorianStepAction.h.

Definition at line 62 of file HistorianStepAction.h.

Definition at line 63 of file HistorianStepAction.h.

Definition at line 64 of file HistorianStepAction.h.

Definition at line 65 of file HistorianStepAction.h.

Definition at line 68 of file HistorianStepAction.h.

Definition at line 69 of file HistorianStepAction.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