#include <SG_StepNtuple.h>

Inheritance diagram for G4UA::SG_StepNtuple:

Collaboration diagram for G4UA::SG_StepNtuple:

Public Member Functions
	SG_StepNtuple (const std::vector< int > &)

virtual void	BeginOfRunAction (const G4Run *) override

virtual void	BeginOfEventAction (const G4Event *) override

virtual void	EndOfEventAction (const G4Event *) override

virtual void	UserSteppingAction (const G4Step *) override

bool	msgLvl (const MSG::Level lvl) const
	Test the output level. More...

MsgStream &	msg () const
	The standard message stream. More...

MsgStream &	msg (const MSG::Level lvl) const
	The standard message stream. More...

void	setLevel (MSG::Level lvl)
	Change the current logging level. More...

Private Member Functions
bool	isSUSYParticle (const int id) const

void	initMessaging () const
	Initialize our message level and MessageSvc. More...

Private Attributes
NTuple::Item< long >	m_nsteps

NTuple::Item< long >	m_evtid

NTuple::Array< int >	m_pdg

NTuple::Array< int >	m_charge

NTuple::Array< int >	m_baryon

NTuple::Array< float >	m_x1

NTuple::Array< float >	m_y1

NTuple::Array< float >	m_z1

NTuple::Array< float >	m_t1

NTuple::Array< float >	m_pt1

NTuple::Array< float >	m_x2

NTuple::Array< float >	m_y2

NTuple::Array< float >	m_z2

NTuple::Array< float >	m_t2

NTuple::Array< float >	m_pt2

NTuple::Array< float >	m_minA

NTuple::Array< float >	m_v2

NTuple::Array< float >	m_vthresh

NTuple::Array< float >	m_vbelowthresh

NTuple::Array< float >	m_dep

NTuple::Array< float >	m_mass

NTuple::Array< float >	m_ke1

NTuple::Array< float >	m_ke2

NTuple::Array< int >	m_rh

NTuple::Array< int >	m_rhid

NTuple::Array< int >	m_step

std::set< int >	m_rhs

long	m_nevents = 0

long	m_rhadronIndex = 0

std::string	m_nm
	Message source name. More...

boost::thread_specific_ptr< MsgStream >	m_msg_tls
	MsgStream instance (a std::cout like with print-out levels) More...

std::atomic< IMessageSvc * >	m_imsg { nullptr }
	MessageSvc pointer. More...

std::atomic< MSG::Level >	m_lvl { MSG::NIL }
	Current logging level. More...

std::atomic_flag m_initialized	ATLAS_THREAD_SAFE = ATOMIC_FLAG_INIT
	Messaging initialized (initMessaging) More...

Detailed Description

Definition at line 20 of file SG_StepNtuple.h.

Constructor & Destructor Documentation

◆ SG_StepNtuple()

G4UA::SG_StepNtuple::SG_StepNtuple ( const std::vector< int > & pdgids )

Definition at line 26 of file SG_StepNtuple.cxx.

     : AthMessaging(Gaudi::svcLocator()->service< IMessageSvc >( "MessageSvc" ),
                    "SG_StepNtuple")
   {
     // Load up the PDG IDs
     for (int i : pdgids ){
       m_rhs.insert(i);
     }
   }

Member Function Documentation

◆ BeginOfEventAction()

void G4UA::SG_StepNtuple::BeginOfEventAction ( const G4Event * )

overridevirtual

Definition at line 91 of file SG_StepNtuple.cxx.

   {
     m_nsteps=0;
     m_rhadronIndex=0;//the rhadron index (either the first or second rhadon, usually)
     m_nevents++;
     m_evtid=m_nevents;//since it gets cleared out after every fill...
   }

◆ BeginOfRunAction()

void G4UA::SG_StepNtuple::BeginOfRunAction ( const G4Run * )

overridevirtual

Definition at line 36 of file SG_StepNtuple.cxx.

   {
     NTupleFilePtr file1(ntupleSvc(), "/NTUPLES/FILE1");
  
     SmartDataPtr<NTuple::Directory>
       ntdir(ntupleSvc(),"/NTUPLES/FILE1/StepNtuple");
  
     if ( !ntdir ) ntdir = ntupleSvc()->createDirectory(file1,"StepNtuple");
     //if ( !ntdir ) log << MSG::ERROR << " failed to get ntuple directory" << endmsg;
     NTuplePtr nt(ntupleSvc(), "/NTUPLES/FILE1/StepNtuple/10");
     if ( !nt ) {    // Check if already booked
       nt = ntupleSvc()->book (ntdir.ptr(), 10,CLID_ColumnWiseTuple, "GEANT4 Step NTuple");
       if ( nt ) {
  
         ATH_MSG_INFO("booked step ntuple ");
  
         if( nt->addItem ("nstep", m_nsteps,0 ,50000)!=StatusCode::SUCCESS // WARNING!! Force limit to 50k tracks
             || nt->addItem ("pdg", m_nsteps, m_pdg)!=StatusCode::SUCCESS
             || nt->addItem ("charge", m_nsteps, m_charge)!=StatusCode::SUCCESS
             || nt->addItem ("mass", m_nsteps, m_mass)!=StatusCode::SUCCESS
             || nt->addItem ("baryon", m_nsteps, m_baryon)!=StatusCode::SUCCESS
             || nt->addItem ("x1", m_nsteps, m_x1)!=StatusCode::SUCCESS
             || nt->addItem ("y1", m_nsteps, m_y1)!=StatusCode::SUCCESS
             || nt->addItem ("z1", m_nsteps, m_z1)!=StatusCode::SUCCESS
             || nt->addItem ("t1", m_nsteps, m_t1)!=StatusCode::SUCCESS
             || nt->addItem ("x2", m_nsteps, m_x2)!=StatusCode::SUCCESS
             || nt->addItem ("y2", m_nsteps, m_y2)!=StatusCode::SUCCESS
             || nt->addItem ("z2", m_nsteps, m_z2)!=StatusCode::SUCCESS
             || nt->addItem ("t2", m_nsteps, m_t2)!=StatusCode::SUCCESS
             || nt->addItem ("dep", m_nsteps, m_dep)!=StatusCode::SUCCESS
             || nt->addItem ("ke1", m_nsteps, m_ke1)!=StatusCode::SUCCESS
             || nt->addItem ("ke2", m_nsteps, m_ke2)!=StatusCode::SUCCESS
             || nt->addItem ("rh", m_nsteps, m_rh)!=StatusCode::SUCCESS
             || nt->addItem ("rhid", m_nsteps, m_rhid)!=StatusCode::SUCCESS
             || nt->addItem ("step", m_nsteps, m_step)!=StatusCode::SUCCESS
             || nt->addItem ("pt1", m_nsteps, m_pt1)!=StatusCode::SUCCESS
             || nt->addItem ("pt2", m_nsteps, m_pt2)!=StatusCode::SUCCESS
             || nt->addItem ("minA",m_nsteps, m_minA)!=StatusCode::SUCCESS
             || nt->addItem ("v2",m_nsteps, m_v2)!=StatusCode::SUCCESS
             || nt->addItem ("vthresh",m_nsteps, m_vthresh)!=StatusCode::SUCCESS
             || nt->addItem ("vbelowthresh",m_nsteps, m_vbelowthresh)!=StatusCode::SUCCESS
             || nt->addItem ("evtid", m_evtid)!=StatusCode::SUCCESS)
  
           ATH_MSG_ERROR("Could not configure branches ");
  
       }
  
       else ATH_MSG_ERROR("Could not book step ntuple!! ");
     }
  
     //set initial values
     m_nevents=0;
  
   }

◆ EndOfEventAction()

void G4UA::SG_StepNtuple::EndOfEventAction ( const G4Event * )

overridevirtual

Definition at line 99 of file SG_StepNtuple.cxx.

   {
     if(! ntupleSvc()->writeRecord("/NTUPLES/FILE1/StepNtuple/10").isSuccess()) {
       ATH_MSG_ERROR( " failed to write record for this event" );
     }
  
     //this also seems to zero out all the arrays... so beware!
   }

◆ initMessaging()

void AthMessaging::initMessaging ( ) const

privateinherited

Initialize our message level and MessageSvc.

This method should only be called once.

Definition at line 39 of file AthMessaging.cxx.

 {
   m_imsg = Athena::getMessageSvc();
   m_lvl = m_imsg ?
     static_cast<MSG::Level>( m_imsg.load()->outputLevel(m_nm) ) :
     MSG::INFO;
 }

◆ isSUSYParticle()

bool G4UA::SG_StepNtuple::isSUSYParticle ( const int id ) const

private

Definition at line 208 of file SG_StepNtuple.cxx.

   {
     if (id==1000021 || id==1000005 || id==1000006 || id==1000512 || id==1000522 || id==1000991 || id==1000993 ||
         id==1000612 || id==1000622 || id==1000632 || id==1000642 || id==1000652 || id==1005211 ||
         id==1006113 || id==1006211 || id==1006213 || id==1006223 || id==1006311 ||
         id==1006313 || id==1006321 || id==1006323 || id==1006333 ||
         id==1009111 || id==1009113 || id==1009211 || id==1009213 || id==1009311 ||
         id==1009313 || id==1009321 || id==1009323 || id==1009223 || id==1009333 ||
         id==1092112 || id==1091114 || id==1092114 || id==1092212 || id==1092214 || id==1092224 ||
         id==1093114 || id==1093122 || id==1093214 || id==1093224 || id==1093314 || id==1093324 || id==1093334)
       return true;
     return false;
   }

◆ msg() [1/2]

MsgStream & AthMessaging::msg ( ) const

inlineinherited

The standard message stream.

Returns a reference to the default message stream May not be invoked before sysInitialize() has been invoked.

Definition at line 164 of file AthMessaging.h.

 {
   MsgStream* ms = m_msg_tls.get();
   if (!ms) {
     if (!m_initialized.test_and_set()) initMessaging();
     ms = new MsgStream(m_imsg,m_nm);
     m_msg_tls.reset( ms );
   }
  
   ms->setLevel (m_lvl);
   return *ms;
 }

◆ msg() [2/2]

MsgStream & AthMessaging::msg ( const MSG::Level lvl ) const

inlineinherited

The standard message stream.

Returns a reference to the default message stream May not be invoked before sysInitialize() has been invoked.

Definition at line 179 of file AthMessaging.h.

180 { return msg() << lvl; }

◆ msgLvl()

bool AthMessaging::msgLvl ( const MSG::Level lvl ) const

inlineinherited

Test the output level.

Parameters

lvl	The message level to test against

Returns: boolean Indicating if messages at given level will be printed

Return values

true	Messages at level "lvl" will be printed

Definition at line 151 of file AthMessaging.h.

 {
   if (!m_initialized.test_and_set()) initMessaging();
   if (m_lvl <= lvl) {
     msg() << lvl;
     return true;
   } else {
     return false;
   }
 }

◆ setLevel()

void AthMessaging::setLevel ( MSG::Level lvl )

inherited

Change the current logging level.

Use this rather than msg().setLevel() for proper operation with MT.

Definition at line 28 of file AthMessaging.cxx.

 {
   m_lvl = lvl;
 }

◆ UserSteppingAction()

void G4UA::SG_StepNtuple::UserSteppingAction ( const G4Step * aStep )

overridevirtual

Definition at line 108 of file SG_StepNtuple.cxx.

   {
     if(m_nsteps<50000){
       int pdg = aStep->GetTrack()->GetDefinition()->GetPDGEncoding();
       bool rhad=false;
       if (std::find(m_rhs.begin(),m_rhs.end(),std::abs(pdg))!=m_rhs.end()) {
         rhad=true;
       }
  
       //
       if (!rhad && isSUSYParticle(std::abs(pdg))){
         std::cout<<"ACH139: SG_StepNtuple: other code thinks "<<pdg<<" is an Rhadron!"<<std::endl;
         rhad=true;
       }
       //
  
       if (rhad){
  
         //
         G4Material * mat = aStep->GetTrack()->GetMaterial();
         double minA=1500000.;
         for (unsigned int i=0;i<mat->GetNumberOfElements();++i){
           if (mat->GetElement(i) && minA>mat->GetElement(i)->GetN()){
             minA=mat->GetElement(i)->GetN();
           }
         }
         //
  
         bool firstslow = aStep->GetPostStepPoint()->GetVelocity()<0.15*std::pow(minA,-2./3.)*CLHEP::c_light;
         //just save the first slow step for the rhadron
         for (int i=0; i<m_nsteps; ++i){
           if (m_rhid[i]==m_rhadronIndex && m_vbelowthresh[i]>0) firstslow=false;
         }
         if (firstslow || aStep->GetTrack()->GetCurrentStepNumber()<=1 || aStep->GetPostStepPoint()->GetKineticEnergy()==0.){
  
           //
           //int id = aStep->GetTrack()->GetDynamicParticle()->GetDefinition()->GetPDGEncoding();
           int id = std::abs(aStep->GetTrack()->GetDynamicParticle()->GetDefinition()->GetPDGEncoding());
           if (id>=1000000 && id<=1100000 && isSUSYParticle(id)){
             m_rh[m_nsteps] = 1;//other code agrees it's an Rhadron
           }
           else{
             m_rh[m_nsteps] = 0;//other code doesn't agree it's an Rhadron
           }
           //
  
           if (aStep->GetPreStepPoint()->GetGlobalTime()==0) m_rhadronIndex++;
           m_rhid[m_nsteps]=m_rhadronIndex;
  
           m_pdg[m_nsteps]=aStep->GetTrack()->GetDefinition()->GetPDGEncoding();
           m_charge[m_nsteps]=aStep->GetTrack()->GetDefinition()->GetPDGCharge();
           m_dep[m_nsteps]=aStep->GetTotalEnergyDeposit();
           m_mass[m_nsteps]=aStep->GetTrack()->GetDefinition()->GetPDGMass();
           m_baryon[m_nsteps]=aStep->GetTrack()->GetDefinition()->GetBaryonNumber();
           m_t1[m_nsteps]=aStep->GetPreStepPoint()->GetGlobalTime();
           m_t2[m_nsteps]=aStep->GetPostStepPoint()->GetGlobalTime();
           G4ThreeVector pos1=aStep->GetPreStepPoint()->GetPosition();
           m_x1[m_nsteps]=pos1.x();
           m_y1[m_nsteps]=pos1.y();
           m_z1[m_nsteps]=pos1.z();
           G4ThreeVector pos2=aStep->GetPostStepPoint()->GetPosition();
           m_x2[m_nsteps]=pos2.x();
           m_y2[m_nsteps]=pos2.y();
           m_z2[m_nsteps]=pos2.z();
           m_ke1[m_nsteps]=aStep->GetPreStepPoint()->GetKineticEnergy();
           m_ke2[m_nsteps]=aStep->GetPostStepPoint()->GetKineticEnergy();
           m_pt1[m_nsteps]=aStep->GetPreStepPoint()->GetMomentum().perp();
           m_pt2[m_nsteps]=aStep->GetPostStepPoint()->GetMomentum().perp();
  
           //
           m_minA[m_nsteps]=minA;
           m_v2[m_nsteps]=aStep->GetPostStepPoint()->GetVelocity();
           m_vthresh[m_nsteps]=0.15*std::pow(minA,-2./3.)*CLHEP::c_light;
           m_vbelowthresh[m_nsteps]=(firstslow?1:0);
           //
  
           m_step[m_nsteps]=m_nsteps;
           ++m_nsteps;
           //std::cout<<"stepping, size is "<<m_nsteps<<std::endl;
         } // writing the step because it stopped or is the start of the "track"
       } //rhad true
       else {
  
         //KILL the particles here, so we don't waste time in GEANT tracking what happens to it!
         if (std::abs(pdg)>1000000 && std::abs(pdg)<10000000){
           std::cout<<"ACH129: SG_StepNtuple: Killing non-rh track with pdg "<<pdg<<std::endl;
         }
         aStep->GetTrack()->SetTrackStatus(fKillTrackAndSecondaries);
         const G4TrackVector *tv = aStep->GetSecondary();
         //if ((*tv).size()>0) std::cout<<" ... and its "<<(*tv).size()<<" secondaries"<<std::endl;
         for (unsigned int i=0;i<tv->size();i++){
           G4Track *t = (*tv)[i];
           t->SetTrackStatus(fKillTrackAndSecondaries);
         }
  
       } // not an rhad
  
     } //m_nsteps<50000
   }

Member Data Documentation

◆ ATLAS_THREAD_SAFE

std::atomic_flag m_initialized AthMessaging::ATLAS_THREAD_SAFE = ATOMIC_FLAG_INIT

mutableprivateinherited

Messaging initialized (initMessaging)

Definition at line 141 of file AthMessaging.h.

◆ m_baryon

NTuple::Array<int> G4UA::SG_StepNtuple::m_baryon

private

Definition at line 38 of file SG_StepNtuple.h.

◆ m_charge

NTuple::Array<int> G4UA::SG_StepNtuple::m_charge

private

Definition at line 38 of file SG_StepNtuple.h.

◆ m_dep

NTuple::Array<float> G4UA::SG_StepNtuple::m_dep

private

Definition at line 42 of file SG_StepNtuple.h.

◆ m_evtid

NTuple::Item<long> G4UA::SG_StepNtuple::m_evtid

private

Definition at line 37 of file SG_StepNtuple.h.

◆ m_imsg

std::atomic<IMessageSvc*> AthMessaging::m_imsg { nullptr }

mutableprivateinherited

MessageSvc pointer.

Definition at line 135 of file AthMessaging.h.

◆ m_ke1

NTuple::Array<float> G4UA::SG_StepNtuple::m_ke1

private

Definition at line 43 of file SG_StepNtuple.h.

◆ m_ke2

NTuple::Array<float> G4UA::SG_StepNtuple::m_ke2

private

Definition at line 43 of file SG_StepNtuple.h.

◆ m_lvl

std::atomic<MSG::Level> AthMessaging::m_lvl { MSG::NIL }

mutableprivateinherited

Current logging level.

Definition at line 138 of file AthMessaging.h.

◆ m_mass

NTuple::Array<float> G4UA::SG_StepNtuple::m_mass

private

Definition at line 42 of file SG_StepNtuple.h.

◆ m_minA

NTuple::Array<float> G4UA::SG_StepNtuple::m_minA

private

Definition at line 41 of file SG_StepNtuple.h.

◆ m_msg_tls

boost::thread_specific_ptr<MsgStream> AthMessaging::m_msg_tls

mutableprivateinherited

MsgStream instance (a std::cout like with print-out levels)

Definition at line 132 of file AthMessaging.h.

◆ m_nevents

long G4UA::SG_StepNtuple::m_nevents = 0

private

Definition at line 46 of file SG_StepNtuple.h.

◆ m_nm

std::string AthMessaging::m_nm

privateinherited

Message source name.

Definition at line 129 of file AthMessaging.h.

◆ m_nsteps

NTuple::Item<long> G4UA::SG_StepNtuple::m_nsteps

private

Definition at line 37 of file SG_StepNtuple.h.

◆ m_pdg

NTuple::Array<int> G4UA::SG_StepNtuple::m_pdg

private

Definition at line 38 of file SG_StepNtuple.h.

◆ m_pt1

NTuple::Array<float> G4UA::SG_StepNtuple::m_pt1

private

Definition at line 39 of file SG_StepNtuple.h.

◆ m_pt2

NTuple::Array<float> G4UA::SG_StepNtuple::m_pt2

private

Definition at line 40 of file SG_StepNtuple.h.

◆ m_rh

NTuple::Array<int> G4UA::SG_StepNtuple::m_rh

private

Definition at line 44 of file SG_StepNtuple.h.

◆ m_rhadronIndex

long G4UA::SG_StepNtuple::m_rhadronIndex = 0

private

Definition at line 47 of file SG_StepNtuple.h.

◆ m_rhid

NTuple::Array<int> G4UA::SG_StepNtuple::m_rhid

private

Definition at line 44 of file SG_StepNtuple.h.

◆ m_rhs

std::set<int> G4UA::SG_StepNtuple::m_rhs

private

Definition at line 45 of file SG_StepNtuple.h.

◆ m_step

NTuple::Array<int> G4UA::SG_StepNtuple::m_step

private

Definition at line 44 of file SG_StepNtuple.h.

◆ m_t1

NTuple::Array<float> G4UA::SG_StepNtuple::m_t1

private

Definition at line 39 of file SG_StepNtuple.h.

◆ m_t2

NTuple::Array<float> G4UA::SG_StepNtuple::m_t2

private

Definition at line 40 of file SG_StepNtuple.h.

◆ m_v2

NTuple::Array<float> G4UA::SG_StepNtuple::m_v2

private

Definition at line 41 of file SG_StepNtuple.h.

◆ m_vbelowthresh

NTuple::Array<float> G4UA::SG_StepNtuple::m_vbelowthresh

private

Definition at line 41 of file SG_StepNtuple.h.

◆ m_vthresh

NTuple::Array<float> G4UA::SG_StepNtuple::m_vthresh

private

Definition at line 41 of file SG_StepNtuple.h.

◆ m_x1

NTuple::Array<float> G4UA::SG_StepNtuple::m_x1

private

Definition at line 39 of file SG_StepNtuple.h.

◆ m_x2

NTuple::Array<float> G4UA::SG_StepNtuple::m_x2

private

Definition at line 40 of file SG_StepNtuple.h.

◆ m_y1

NTuple::Array<float> G4UA::SG_StepNtuple::m_y1

private

Definition at line 39 of file SG_StepNtuple.h.

◆ m_y2

NTuple::Array<float> G4UA::SG_StepNtuple::m_y2

private

Definition at line 40 of file SG_StepNtuple.h.

◆ m_z1

NTuple::Array<float> G4UA::SG_StepNtuple::m_z1

private

Definition at line 39 of file SG_StepNtuple.h.

◆ m_z2

NTuple::Array<float> G4UA::SG_StepNtuple::m_z2

private

Definition at line 40 of file SG_StepNtuple.h.

The documentation for this class was generated from the following files:

Public Member Functions

Private Member Functions

Private Attributes

Detailed Description

Constructor & Destructor Documentation

◆ SG_StepNtuple()

Member Function Documentation

◆ BeginOfEventAction()

◆ BeginOfRunAction()

◆ EndOfEventAction()

◆ initMessaging()

◆ isSUSYParticle()

◆ msg() [1/2]

◆ msg() [2/2]

◆ msgLvl()

◆ setLevel()

◆ UserSteppingAction()

Member Data Documentation

◆ ATLAS_THREAD_SAFE

◆ m_baryon

◆ m_charge

◆ m_dep

◆ m_evtid

◆ m_imsg

◆ m_ke1

◆ m_ke2

◆ m_lvl

◆ m_mass

◆ m_minA

◆ m_msg_tls

◆ m_nevents

◆ m_nm

◆ m_nsteps

◆ m_pdg

◆ m_pt1

◆ m_pt2

◆ m_rh

◆ m_rhadronIndex

◆ m_rhid

◆ m_rhs

◆ m_step

◆ m_t1

◆ m_t2

◆ m_v2

◆ m_vbelowthresh

◆ m_vthresh

◆ m_x1

◆ m_x2

◆ m_y1

◆ m_y2

◆ m_z1

◆ m_z2