Inheritance diagram for filt::GenFilter:

Classes
struct	CounterSetPair

Public Member Functions
	GenFilter (fhicl::ParameterSet const &pset)

virtual	~GenFilter ()

virtual bool	filter (art::Event &e)

void	reconfigure (fhicl::ParameterSet const &pset)

void	beginJob ()

Public Member Functions inherited from art::EDFilter
	EDFilter (fhicl::ParameterSet const &pset)

template<typename Config >
	EDFilter (Table< Config > const &config)

std::string	workerType () const

Public Member Functions inherited from art::detail::Filter
virtual	~Filter () noexcept

	Filter (fhicl::ParameterSet const &)

	Filter (Filter const &)=delete

	Filter (Filter &&)=delete

Filter &	operator= (Filter const &)=delete

Filter &	operator= (Filter &&)=delete

void	doBeginJob (SharedResources const &resources)

void	doEndJob ()

void	doRespondToOpenInputFile (FileBlock const &fb)

void	doRespondToCloseInputFile (FileBlock const &fb)

void	doRespondToOpenOutputFiles (FileBlock const &fb)

void	doRespondToCloseOutputFiles (FileBlock const &fb)

bool	doBeginRun (RunPrincipal &rp, ModuleContext const &)

bool	doEndRun (RunPrincipal &rp, ModuleContext const &mc)

bool	doBeginSubRun (SubRunPrincipal &srp, ModuleContext const &mc)

bool	doEndSubRun (SubRunPrincipal &srp, ModuleContext const &mc)

bool	doEvent (EventPrincipal &ep, ModuleContext const &mc, std::atomic< std::size_t > &counts_run, std::atomic< std::size_t > &counts_passed, std::atomic< std::size_t > &counts_failed)

Public Member Functions inherited from art::Modifier
	~Modifier () noexcept

	Modifier ()

	Modifier (Modifier const &)=delete

	Modifier (Modifier &&)=delete

Modifier &	operator= (Modifier const &)=delete

Modifier &	operator= (Modifier &&)=delete

Public Member Functions inherited from art::ModuleBase
virtual	~ModuleBase () noexcept

	ModuleBase ()

ModuleDescription const &	moduleDescription () const

void	setModuleDescription (ModuleDescription const &)

std::array< std::vector< ProductInfo >, NumBranchTypes > const &	getConsumables () const

void	sortConsumables (std::string const &current_process_name)

template<typename T , BranchType BT>
ViewToken< T >	consumesView (InputTag const &tag)

template<typename T , BranchType BT>
ViewToken< T >	mayConsumeView (InputTag const &tag)

Private Member Functions
bool	IsInterestingParticle (const simb::MCParticle &particle)

bool	ParticleHitsCounterSetPairs (const simb::MCParticle &particle, const CounterSetPair &CSP)

bool	ParticleHitsCounterSet (const simb::MCParticle &particle, const std::vector< geo::AuxDetGeo const * > &counters, const TVector3 &counter_norm)

Private Attributes
std::vector< CounterSetPair >	fCounterSetPairs

bool	fUseEWCounterPair

bool	fUseNupSdownCounterPair

bool	fUseNdownSupCounterPair

std::vector< int >	fInterestingPDGs

double	fParticleMinEnergy

double	fParticleMaxEnergy

double	fCounterSizeScaleFactor

Additional Inherited Members
Public Types inherited from art::EDFilter
using	ModuleType = EDFilter

using	WorkerType = WorkerT< EDFilter >

Public Types inherited from art::detail::Filter
template<typename UserConfig >
using	Table = Modifier::Table< UserConfig >

Public Types inherited from art::Modifier
template<typename UserConfig , typename UserKeysToIgnore = void>
using	Table = ProducerTable< UserConfig, detail::ModuleConfig, UserKeysToIgnore >

Static Public Attributes inherited from art::detail::Filter
static constexpr bool	Pass {true}

static constexpr bool	Fail {false}

Protected Member Functions inherited from art::ModuleBase
ConsumesCollector &	consumesCollector ()

template<typename T , BranchType = InEvent>
ProductToken< T >	consumes (InputTag const &)

template<typename Element , BranchType = InEvent>
ViewToken< Element >	consumesView (InputTag const &)

template<typename T , BranchType = InEvent>
void	consumesMany ()

template<typename T , BranchType = InEvent>
ProductToken< T >	mayConsume (InputTag const &)

template<typename Element , BranchType = InEvent>
ViewToken< Element >	mayConsumeView (InputTag const &)

template<typename T , BranchType = InEvent>
void	mayConsumeMany ()

Detailed Description

Definition at line 18 of file GenExtCounterFilter_module.cc.

Constructor & Destructor Documentation

filt::GenFilter::GenFilter ( fhicl::ParameterSet const & pset )

explicit

virtual filt::GenFilter::~GenFilter ( )

inlinevirtual

Definition at line 21 of file GenExtCounterFilter_module.cc.

21 {};

Member Function Documentation

void filt::GenFilter::beginJob ( )

virtual

Reimplemented from art::EDFilter.

Definition at line 126 of file GenExtCounterFilter_module.cc.

                            {
     art::ServiceHandle<geo::Geometry> geom;
 
     CounterSetPair EWCounterSetPair;
     CounterSetPair NupSdownCounterSetPair;
     CounterSetPair NdownSupCounterSetPair;
 
     for (unsigned int i = 0; i < geom->NAuxDets(); i++){
       //The WE counter pairs
       geo::AuxDetGeo const* counter = &(geom->AuxDet(i));
       if (i >=6 && i <= 15) EWCounterSetPair.setA.push_back(counter);
       else if (i >= 28 && i <=37) EWCounterSetPair.setB.push_back(counter);
       //The N (up) S (down) counter pairs
       else if (i >= 22 && i <= 27) NupSdownCounterSetPair.setA.push_back(counter);
       else if (i <= 5) NupSdownCounterSetPair.setB.push_back(counter);
       //The N (down) S (up) counter pairs
       else if (i >= 16 && i <= 21) NdownSupCounterSetPair.setA.push_back(counter);
       else if (i >= 38 && i <= 43) NdownSupCounterSetPair.setB.push_back(counter);
     }
 
     EWCounterSetPair.normalVec[0] = 0.;
     EWCounterSetPair.normalVec[1] = 0.;
     EWCounterSetPair.normalVec[2] = 1.;
     EWCounterSetPair.isRequested = fUseEWCounterPair;
     fCounterSetPairs.push_back(EWCounterSetPair);
 
     NupSdownCounterSetPair.normalVec[0] = 1.;
     NupSdownCounterSetPair.normalVec[1] = 0.;
     NupSdownCounterSetPair.normalVec[2] = 0.;
     NupSdownCounterSetPair.isRequested = fUseNupSdownCounterPair;
     fCounterSetPairs.push_back(NupSdownCounterSetPair);
 
     NdownSupCounterSetPair.normalVec[0] = 1.;
     NdownSupCounterSetPair.normalVec[1] = 0.;
     NdownSupCounterSetPair.normalVec[2] = 0.;
     NdownSupCounterSetPair.isRequested = fUseNdownSupCounterPair;
     fCounterSetPairs.push_back(NdownSupCounterSetPair);
 
     for (unsigned int i = 0; i < fCounterSetPairs.size(); i++){
       std::cout<<"Counter set pair: "<<i<<std::endl;
       std::cout<<"--setA size: " << fCounterSetPairs[i].setA.size() << std::endl;
       std::cout<<"--setB size: " << fCounterSetPairs[i].setB.size() << std::endl;
 
     }
 
 
   }

bool filt::GenFilter::filter ( art::Event & e )

virtual

Implements art::EDFilter.

Definition at line 81 of file GenExtCounterFilter_module.cc.

                                     {
 
     auto const mclists = e.getMany<std::vector<simb::MCTruth>>();
     for (unsigned int i = 0; i < mclists.size() ; i++){
       for (unsigned int j = 0; j < mclists[i]->size(); j++){
         //Should have the truth record for the event now
         const art::Ptr<simb::MCTruth> mc_truth(mclists[i],j);
         for (int part = 0; part < mc_truth->NParticles(); part++){
           const simb::MCParticle particle = mc_truth->GetParticle(part);
           if (!IsInterestingParticle(particle)) continue;
           TVector3 particle_pos = particle.Position().Vect();
           TVector3 particle_dir = particle.Momentum().Vect().Unit();
           for (unsigned CSP_i = 0; CSP_i < fCounterSetPairs.size(); CSP_i++){
             if (!fCounterSetPairs[CSP_i].isRequested) continue;
             if (ParticleHitsCounterSetPairs(particle,fCounterSetPairs[i])){
               std::cout<<"HIT COUNTER SET"<<std::endl;
               return true;
             }
             /*
             TVector3 counter_norm(fCounterSetPairs[i].normalVec[0],fCounterSetPairs[i].normalVec[1],fCounterSetPairs[i].normalVec[2]);
             double counter_pos_array[3];
             fCounterSetPairs[CSP_i].setA.front()->GetCenter(counter_pos_array);
             TVector3 counter_pos(counter_pos_array[0], counter_pos_array[1], counter_pos_array[2]);
             double scale_factor = counter_norm.Dot(counter_pos - particle_pos)/(counter_norm.Dot(particle_dir));
             TVector3 particle_pos_in_plane = particle_pos + scale_factor * particle_dir;
             std::cout<<"Particle pos: " << "("<<particle_pos_in_plane.X()<<","<<particle_pos_in_plane.Y()<<","<<particle_pos_in_plane.Z()<<")"<<std::endl;
             std::cout<<"Counter pos: " << "("<<counter_pos.X()<<","<<counter_pos.Y()<<","<<counter_pos.Z()<<")"<<std::endl;
             */
 
           }
           /*
           geo::AuxDetGeo const*counter = fCounterSetPairs.front().setA.front();
           std::cout<<"Length: " << counter->Length() << std::endl;
           std::cout<<"HalfWidth1: " << counter->HalfWidth1() << std::endl;
           std::cout<<"HalfWidth2: " << counter->HalfWidth2() << std::endl;
 
           std::cout<<"HalfHeight: " << counter->HalfHeight() << std::endl;
           */
         }
       }
     }
 
     return false;
   }

bool filt::GenFilter::IsInterestingParticle ( const simb::MCParticle & particle )

private

Definition at line 174 of file GenExtCounterFilter_module.cc.

                                                                      {
 
     for (unsigned int i = 0; i < fInterestingPDGs.size(); i++){
       //Check if the particle under consideration has a requested PDG
       if (particle.PdgCode() == fInterestingPDGs[i]){
         //Got a requested PDG,  now check that the energy matches the requested range
         TLorentzVector mom4 = particle.Momentum();
         if (mom4.T() > fParticleMinEnergy && mom4.T() < fParticleMaxEnergy){
           //std::cout<<"FOUND INTERESTING PARTICLE"<<std::endl;
           return true;
         }
       }
     }
 
     return false;
   }

bool filt::GenFilter::ParticleHitsCounterSet	(	const simb::MCParticle &	particle,
		const std::vector< geo::AuxDetGeo const * > &	counters,
		const TVector3 &	counter_norm
	)

private

Definition at line 204 of file GenExtCounterFilter_module.cc.

                                                                                                                                                     {
 
     //Loop through the counters
     for (unsigned int i = 0; i < counters.size(); i++){
       //First step is to push the particle to the counter plane
       geo::AuxDetGeo const*counter = counters[i];
 
       double counter_pos_array[3];
       //fCounterSetPairs[CSP_i].setA.front()->GetCenter(counter_pos_array);
       counter->GetCenter(counter_pos_array);
       TVector3 counter_pos(counter_pos_array[0], counter_pos_array[1], counter_pos_array[2]);
       TVector3 particle_pos = particle.Position().Vect();
       TVector3 particle_dir = particle.Momentum().Vect().Unit();
 
       /*
           Length: 62.992
           HalfWidth1: 16.2814
           HalfWidth2: 13.5255
           HalfHeight: 0.475
           */
 
       double scale_factor = counter_norm.Dot(counter_pos - particle_pos)/(counter_norm.Dot(particle_dir));
 
       TVector3 particle_pos_in_plane = particle_pos + scale_factor * particle_dir;
 
       //We now have the particle position in the plane of the counter.  We now just need to calculate whether it sits inside the box
       //Create two TVector3s, each will hold the coordinates of opposing corners of the counter in question
       //We need to use the normals associated with the counter.  We already have two, one is a member of the C++ object and the other was passed to this function
       //Create the two TVector3s
       TVector3 pos_corner, neg_corner;
       //Lets start with the one passed to this function
       //The thin dimension of the counter is the one associated with normal passed to this function
       pos_corner += counter->HalfHeight()*counter_norm*fCounterSizeScaleFactor;
       neg_corner += -1.*counter->HalfHeight()*counter_norm*fCounterSizeScaleFactor; 
 
       //Now lets to the same for the normal stored in the C++ object (this "normal" actually points out the top of a counter)
       double counter_top_norm_array[3];
       counter->GetNormalVector(counter_top_norm_array);
       //Now package this up a TVector
       TVector3 counter_top_norm;
       counter_top_norm.SetX(counter_top_norm_array[0]);
       counter_top_norm.SetY(counter_top_norm_array[1]);
       counter_top_norm.SetZ(counter_top_norm_array[2]);
       //OK now we can add the dimension to the corner vectors.  The relevant dimension this time is Length/2
       pos_corner += counter->Length()*counter_top_norm*0.5*fCounterSizeScaleFactor;
       neg_corner += -1.*counter->Length()*counter_top_norm*0.5*fCounterSizeScaleFactor;
       //Now we need to the same for the vector pointing along the counter.  
       //The relevant dimension in this case in HalfWidth1 (going to assume the counter is a square and take the bigger width)
       //Because we have the other two vectors already, we can very easily get the final one by taking the cross product of them
       TVector3 counter_side_norm = counter_norm.Cross(counter_top_norm); 
       //now add the dimenions onto the corner vectors
       pos_corner += counter->HalfWidth1()*counter_side_norm*fCounterSizeScaleFactor;
       neg_corner += -1.*counter->HalfWidth1()*counter_side_norm*fCounterSizeScaleFactor;
 
       //Almost ready
       //We have calculated the corners assuming the centre of the counter is the origin.  Two choices, either translate the corners OR translate the propagated particle position
       //The latter is less lines of code so lets to that
       particle_pos_in_plane += -1.*counter_pos;
 
       //We are now ready to check if the particle sits in the counter
       //We don't know by default which way round the corners are oriented, but that doesn't matter as we can take abs, max and min
       //We are going to take abs of the particle position, the pos corner and the neg corner first
       /*
       pos_corner.SetXYZ(std::abs(pos_corner.X()),std::abs(pos_corner.Y()),std::abs(pos_corner.Z()));
       neg_corner.SetXYZ(std::abs(neg_corner.X()),std::abs(neg_corner.Y()),std::abs(neg_corner.Z()));
       particle_pos_in_plane.SetXYZ(std::abs(particle_pos_in_plane.X()),std::abs(particle_pos_in_plane.Y()),std::abs(particle_pos_in_plane.Z()));
       */
 
       /*
       //Dump the corners and particle pos
       std::cout<<"pos_corner: " << pos_corner.X()<<","<<pos_corner.Y()<<","<<pos_corner.Z()<<std::endl;
       std::cout<<"neg_corner: " << neg_corner.X()<<","<<neg_corner.Y()<<","<<neg_corner.Z()<<std::endl;
       std::cout<<"particle_pos_in_plane: " << particle_pos_in_plane.X()<<","<<particle_pos_in_plane.Y()<<","<<particle_pos_in_plane.Z()<<std::endl;
       */
 
 
 
       //Now we can check
       //This is going to be one huge if statement...
       if (particle_pos_in_plane.X() > std::min(pos_corner.X(), neg_corner.X()) && particle_pos_in_plane.X() < std::max(pos_corner.X(),neg_corner.X()) && particle_pos_in_plane.Y() > std::min(pos_corner.Y(), neg_corner.Y()) && particle_pos_in_plane.Y() < std::max(pos_corner.Y(),neg_corner.Y()) && particle_pos_in_plane.Z() > std::min(pos_corner.Z(), neg_corner.Z()) && particle_pos_in_plane.Z() < std::max(pos_corner.Z(),neg_corner.Z())){
         std::cout<<"Particle uses counter in set"<<std::endl;
         return true;
       }
 
     }
 
     return false;
   }

bool filt::GenFilter::ParticleHitsCounterSetPairs	(	const simb::MCParticle &	particle,
		const CounterSetPair &	CSP
	)

private

Definition at line 191 of file GenExtCounterFilter_module.cc.

                                                                                                       {
     //Need the normal to the counters
     TVector3 counter_norm(CSP.normalVec[0],CSP.normalVec[1],CSP.normalVec[2]);
 
     //Loop through one of the counter sets
     if (ParticleHitsCounterSet(particle, CSP.setA, counter_norm)){
       if (ParticleHitsCounterSet(particle, CSP.setB, counter_norm)){
         return true;
       }
     }
     return false;
   }

void filt::GenFilter::reconfigure ( fhicl::ParameterSet const & pset )

Definition at line 59 of file GenExtCounterFilter_module.cc.

                                                           {
     fUseEWCounterPair = pset.get<bool>("UseEWCounterPair",1);
     std::cout<<"Use EW counter pair: " << fUseEWCounterPair<<std::endl;
     fUseNupSdownCounterPair = pset.get<bool>("UseNupSdownCounterPair",1);
     std::cout<<"Use N (up) S (down) counter pair: " << fUseNupSdownCounterPair << std::endl;
     fUseNdownSupCounterPair = pset.get<bool>("UseNdownSupCounterPair",1);
     std::cout<<"Use N (down) S (up) counter pair: " << fUseNdownSupCounterPair << std::endl;
     fInterestingPDGs = pset.get<std::vector<int> >("InterestingPDGs");
     std::cout<<"NInteresting PDGs: " << fInterestingPDGs.size() << std::endl;
     for (unsigned int i = 0; i < fInterestingPDGs.size(); i++){
       std::cout<<"-- PDG: " << fInterestingPDGs[i] << std::endl; 
     }
     fParticleMinEnergy = pset.get<double>("ParticleMinEnergy",0);
     std::cout<<"Particle min energy: " << fParticleMinEnergy << std::endl;
     fParticleMaxEnergy = pset.get<double>("ParticleMaxEnergy",999999999);
     std::cout<<"Particle max energy: " << fParticleMaxEnergy << std::endl;
 
     fCounterSizeScaleFactor = pset.get<double>("CounterSizeScaleFactor",1.);
     std::cout<<"Counter size scale factor: " << fCounterSizeScaleFactor << std::endl;
   }