| Message: Re: Finding (grand^n)mother particle of a track [solved] | Not Logged In (login) |
|
Thank you!
This seems kindof similar to my current solution, which is using UserTrackingAction to store the data internaly in the UTA class. However you avoids having making the filter having to access the UTA-class. The solution I eventually came up with was: //Main filters:
#ifndef FilterMotherIon_hh #define FilterMotherIon_hh
#include "G4VSDFilter.hh" #include "globals.hh"
/* * Filter for a simple detector, in order to filter * hits depending on the originally decaying ion. * * True if hit comes from MotherIon or a child of MotherIon, else false. */ class EdepSpectraTrackingAction;
class FilterMotherIon : public G4VSDFilter {
public:
FilterMotherIon(G4String name, G4int z, G4int a) :
G4VSDFilter(name), Z(z), A(a) {};
G4bool Accept(const G4Step* origStep) const;
private: G4int Z; G4int A; }; #endif
#include "FilterMotherIon.hh" #include "G4Track.hh" #include "G4ParticleDefinition.hh"
#include "EdepSpectraTrackingAction.hh" #include "G4RunManager.hh"
G4bool FilterMotherIon::Accept(const G4Step* origStep) const {
//We need the trackingAction to do the digging
EdepSpectraTrackingAction* trkAct = (EdepSpectraTrackingAction*) G4RunManager::GetRunManager()->GetUserTrackingAction();
/* * Run down the particle tree until either condition is met: * - Correct particle is found (accept) * - Pre-primary particle is hit (not found - reject) */
//Initialize on current particle
G4int thisID = origStep->GetTrack()->GetTrackID();
G4int parentID = origStep->GetTrack()->GetParentID();
G4ParticleDefinition* thisDef = origStep->GetTrack()->GetDefinition();
//G4cout << "StartSearch, Z = " << Z << ", A = " << A << G4endl;
//Search!
while (true) {
if (thisDef->GetAtomicNumber() == Z && thisDef->GetAtomicMass() == A) {
//Gold!
//G4cout << "Gold! ID = " << thisID << G4endl;
return true;
}
else {
//No luck. Try to dig deeper
if (parentID == 0) {
//G4cout << "Rock" << G4endl;
return false; //Rock bottom, no gold.
}
//Move the definition one step deeper
thisID = parentID;
parentID = trkAct->GetParentID(thisID);
thisDef = trkAct->GetDefByID(thisID);
//G4cout << "diggin'..." << G4endl;
}
}
}
//Some simple boolean filters (I have only tested and used the Not filter: #ifndef FilterBools_hh #define FilterBools_hh
#include "G4VSDFilter.hh" #include "globals.hh"
/* * Set of boolean operations between two filters */
class FilterAnd : public G4VSDFilter {
public:
FilterAnd(G4String name, G4VSDFilter* filter1, G4VSDFilter* filter2) :
G4VSDFilter(name), Filter1(filter1), Filter2(filter2) {};
inline G4bool Accept(const G4Step* origStep) const {
return (Filter1->Accept(origStep) && Filter2->Accept(origStep));
}
private: G4VSDFilter* Filter1; G4VSDFilter* Filter2;
};
class FilterOr : public G4VSDFilter {
public:
FilterOr(G4String name, G4VSDFilter* filter1, G4VSDFilter* filter2) :
G4VSDFilter(name), Filter1(filter1), Filter2(filter2) {};
inline G4bool Accept(const G4Step* origStep) const {
return (Filter1->Accept(origStep) || Filter2->Accept(origStep));
}
private: G4VSDFilter* Filter1; G4VSDFilter* Filter2; };
class FilterNot : public G4VSDFilter {
public:
FilterNot(G4String name, G4VSDFilter* filter) :
G4VSDFilter(name), Filter(filter) {};
inline G4bool Accept(const G4Step* origStep) const {
return (not Filter->Accept(origStep));
}
private: G4VSDFilter* Filter; }; #endif //UserTrackingAction class: ifndef EdepSpectraTrackingAction_hh #define EdepSpectraTrackingAction_hh
#include "G4UserTrackingAction.hh" #include "G4ParticleDefinition.hh" #include <map>
class EdepSpectraTrackingAction : public G4UserTrackingAction{
public:
EdepSpectraTrackingAction();
void PreUserTrackingAction(const G4Track*);
//Clear the maps for a new event void ResetMaps();
inline G4int GetParentID(G4int child) {
return parentMap[child];
}
inline G4ParticleDefinition* GetDefByID(G4int ID) {
return defMap[ID];
}
private:
//What trackID <value> is parent to trackID <key> ? std::map<G4int, G4int> parentMap; // trackID <key>, particleDefinition<value> std::map<G4int, G4ParticleDefinition*> defMap; }; #endif include "EdepSpectraTrackingAction.hh" #include "G4Track.hh" EdepSpectraTrackingAction::EdepSpectraTrackingAction() { }
void EdepSpectraTrackingAction::PreUserTrackingAction(const G4Track* track) {
//Method only called once pr. track - no need to check if already added
G4int trackID = track->GetTrackID();
parentMap[trackID] = track->GetParentID(); defMap[trackID] = track->GetDefinition(); }
void EdepSpectraTrackingAction::ResetMaps() {
//DEBUG: print the maps
/*
std::map<G4int,G4int>::iterator it;
for ( it=parentMap.begin() ; it != parentMap.end(); it++ ) {
G4cout << (*it).first << " => " << (*it).second << G4endl;
}
*/
parentMap.clear(); defMap.clear(); } This is a somewhat simpler solution, and it does what I need it to do. Note that I need to calle EdepSpectraTrackingAction::ResetMaps() before each event -- i do this in UserPrimaryGeneratorAction::GeneratePrimaries() just before actually generating any primaries. --- Kyrre
|
|
to: |