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

A mix-in virtual class used with the HistorianStepAction. More...

#include <QueriableStepAction.h>

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

List of all members.

Public Member Functions

 QueriableStepAction (const std::string &type, const std::string &name, const IInterface *parent)
virtual ~QueriableStepAction ()
virtual StatusCode initialize ()
virtual StatusCode finalize ()
virtual void UserSteppingAction (const G4Step *)=0
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

 QueriableStepAction ()
 no default constructor
 QueriableStepAction (const QueriableStepAction &)
 no copy
QueriableStepActionoperator= (const QueriableStepAction &)
 no =

Private Attributes

INeutronCaptureInfom_capinfo
double m_BirksConstant1
double m_BirksConstant2

Detailed Description

A mix-in virtual class used with the HistorianStepAction.

Author:
Nathaniel Tagg (tagg@minos.phy.tufts.edu)

Definition at line 30 of file QueriableStepAction.h.


Member Typedef Documentation

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

Definition at line 88 of file QueriableStepAction.h.


Member Enumeration Documentation

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

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

Definition at line 34 of file QueriableStepAction.cc.

  : GiGaStepActionBase( type , name , parent )
  , m_DetectorElementSearchPath(1,"/dd/Structure")
  , mCurrentStep(0)
  , mCurrentTrack(0)
  , mCurrentStepPoint(0)
  , mTouchToDetElem(0)
{ 
  declareProperty("TouchableToDetelem",m_TouchableToDetelem_name = "TH2DE");
  declareProperty("DetectorElementSearchPath",m_DetectorElementSearchPath);

  m_IdParameterNames.push_back("PmtID");
  m_IdParameterNames.push_back("DetectorID");
  m_IdParameterNames.push_back("SiteID");
  declareProperty("IdParameterNames",m_IdParameterNames,
                  "Names of the user parameters attached to detector elements used to give ID numbers.");
  
  // Liang Zhan, Jan 18, 2009.
  declareProperty("BirksConstant1", m_BirksConstant1 = 6.5e-3*g/cm2/MeV, "Birks constant C1");
  declareProperty("BirksConstant2", m_BirksConstant2 = 2.1e-6*(g/cm2/MeV)*(g/cm2/MeV), "Birks constant C2");
}
virtual QueriableStepAction::~QueriableStepAction ( ) [inline, virtual]

Definition at line 36 of file QueriableStepAction.h.

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

no default constructor

QueriableStepAction::QueriableStepAction ( const QueriableStepAction ) [private]

no copy


Member Function Documentation

StatusCode QueriableStepAction::initialize ( ) [virtual]

Reimplemented in HistorianStepAction, and UnObserverStepAction.

Definition at line 59 of file QueriableStepAction.cc.

{
  StatusCode sc =  GiGaStepActionBase::initialize();
  if(sc.isFailure()) return sc;

  std::string t2de_name = std::string("t2de_") + name();
  mTouchToDetElem = tool<ITouchableToDetectorElement>(
     m_TouchableToDetelem_name,
     name(),
     this);

  // Get the Neutron Capture tool which is updated in DsPhysConHadron.cc
  // to gain access to the neutron capture information.
  // --- Wei, Aug 15, 2008
  m_capinfo = tool<INeutronCaptureInfo>("G4DhNeutronCaptureInfoTool");

  verbose() << " QueriableStepAction::initialize() done" << endreq;

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

Reimplemented in HistorianStepAction, and UnObserverStepAction.

Definition at line 80 of file QueriableStepAction.cc.

virtual void QueriableStepAction::UserSteppingAction ( const G4Step *  ) [pure virtual]
StatusCode QueriableStepAction::GetTrackParameterList ( RuleParser::ParameterList tplist) [virtual]

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]

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 
)

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]

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]

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

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

Definition at line 557 of file QueriableStepAction.cc.

const ILVolume * QueriableStepAction::getLogicalVolume ( ) const

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

Definition at line 572 of file QueriableStepAction.cc.

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

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

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

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

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

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;
}
QueriableStepAction& QueriableStepAction::operator= ( const QueriableStepAction ) [private]

no =


Member Data Documentation

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

Definition at line 68 of file QueriableStepAction.h.

Definition at line 69 of file QueriableStepAction.h.

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

Definition at line 70 of file QueriableStepAction.h.

const G4Step* QueriableStepAction::mCurrentStep [protected]

Definition at line 73 of file QueriableStepAction.h.

const G4Track* QueriableStepAction::mCurrentTrack [protected]

Definition at line 74 of file QueriableStepAction.h.

const G4StepPoint* QueriableStepAction::mCurrentStepPoint [protected]

Definition at line 75 of file QueriableStepAction.h.

Definition at line 76 of file QueriableStepAction.h.

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

Definition at line 80 of file QueriableStepAction.h.

int QueriableStepAction::mDetElementMatch [mutable, protected]

Definition at line 81 of file QueriableStepAction.h.

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

Definition at line 82 of file QueriableStepAction.h.

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

Definition at line 83 of file QueriableStepAction.h.

Definition at line 84 of file QueriableStepAction.h.

double QueriableStepAction::mQuenchedEnergy [mutable, protected]

Definition at line 85 of file QueriableStepAction.h.

Definition at line 86 of file QueriableStepAction.h.

unsigned int QueriableStepAction::mDetectorId [mutable, protected]

Definition at line 87 of file QueriableStepAction.h.

Definition at line 89 of file QueriableStepAction.h.

Definition at line 174 of file QueriableStepAction.h.

Definition at line 177 of file QueriableStepAction.h.

Definition at line 178 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