FairRoot/PandaRoot
Public Types | Public Member Functions | Protected Attributes | Static Protected Attributes | Private Member Functions | Private Attributes | Static Private Attributes | List of all members
FairEvtFilterOnSingleParticleCounts Class Reference

#include <FairEvtFilterOnSingleParticleCounts.h>

Inheritance diagram for FairEvtFilterOnSingleParticleCounts:
FairEvtFilter

Public Types

enum  ChargeState {
  kNeutral =0, kPlus, kMinus, kCharged,
  kAll, kChargeLastElement
}
 
enum  MomState { kMomTotal =0, kMomTrans, kMomZ, kMomLastElement }
 
enum  GeomState {
  kTheta =0, kPhi, kVertexZ, kVertexRho,
  kVertexRadius, kGeomLastElement
}
 
enum  LogicOp { kAnd =0, kOr, kLogicOpLastElement }
 

Public Member Functions

 FairEvtFilterOnSingleParticleCounts ()
 
 FairEvtFilterOnSingleParticleCounts (const char *name, const char *title="FairEvtFilterOnSingleParticleCounts")
 
virtual ~FairEvtFilterOnSingleParticleCounts ()
 
Bool_t AndMinMaxPdgCodes (Int_t min, Int_t max, Int_t pdgCode1, Int_t pdgCode2=kInvalidPdgCode, Int_t pdgCode3=kInvalidPdgCode, Int_t pdgCode4=kInvalidPdgCode, Int_t pdgCode5=kInvalidPdgCode, Int_t pdgCode6=kInvalidPdgCode, Int_t pdgCode7=kInvalidPdgCode, Int_t pdgCode8=kInvalidPdgCode)
 
Bool_t AndMinMaxPdgCodes (Int_t min, Int_t max, std::vector< Int_t > &pdgCodes)
 
Bool_t AndMinPdgCodes (Int_t min, Int_t pdgCode1, Int_t pdgCode2=kInvalidPdgCode, Int_t pdgCode3=kInvalidPdgCode, Int_t pdgCode4=kInvalidPdgCode, Int_t pdgCode5=kInvalidPdgCode, Int_t pdgCode6=kInvalidPdgCode, Int_t pdgCode7=kInvalidPdgCode, Int_t pdgCode8=kInvalidPdgCode)
 
Bool_t AndMinPdgCodes (Int_t min, std::vector< Int_t > &pdgCodes)
 
Bool_t AndMaxPdgCodes (Int_t max, std::vector< Int_t > &pdgCodes)
 
Bool_t AndMaxPdgCodes (Int_t max, Int_t pdgCode1, Int_t pdgCode2=kInvalidPdgCode, Int_t pdgCode3=kInvalidPdgCode, Int_t pdgCode4=kInvalidPdgCode, Int_t pdgCode5=kInvalidPdgCode, Int_t pdgCode6=kInvalidPdgCode, Int_t pdgCode7=kInvalidPdgCode, Int_t pdgCode8=kInvalidPdgCode)
 
Bool_t AndMinMaxCharge (Int_t min, Int_t max, ChargeState charge)
 
Bool_t AndMinCharge (Int_t min, ChargeState charge)
 
Bool_t AndMaxCharge (Int_t max, ChargeState charge)
 
Bool_t AndMinMaxAllParticles (Int_t min, Int_t max)
 
Bool_t AndMinAllParticles (Int_t min)
 
Bool_t AndMaxAllParticles (Int_t max)
 
Bool_t AndPRange (Double_t min, Double_t max)
 
Bool_t AndPtRange (Double_t min, Double_t max)
 
Bool_t AndPzRange (Double_t min, Double_t max)
 
Bool_t AndThetaRange (Double_t min, Double_t max)
 
Bool_t AndPhiRange (Double_t min, Double_t max)
 
Bool_t AndVzRange (Double_t min, Double_t max)
 
Bool_t AndRhoRange (Double_t min, Double_t max)
 
Bool_t AndRadiusRange (Double_t min, Double_t max)
 
Bool_t EventMatches (Int_t evtNr)
 
Bool_t FilterActive ()
 
Bool_t SetListOfParticles (TClonesArray *ParticleList)
 
void PrintAllTParticleInEvent ()
 
Bool_t Init ()
 
void SetVerbose (Int_t verbose=12)
 
void ShowAcceptedEventNumbers ()
 
void ShowEvtNrsToAccept ()
 
void SetTestMode (Int_t *arrayPtr, Int_t nEntries)
 
Bool_t TestPassed ()
 
Bool_t GetCharge (Int_t inPdgCode, Double_t *pdgCodeCharge)
 

Protected Attributes

TDatabasePDG * fdbPdg
 
TClonesArray * fParticleList
 
Int_t fVerbose
 
Bool_t fTestMode
 
std::set< Int_t > fAcceptedEventNumbers
 
std::set< Int_t > fEvtNrsToAccept
 
Int_t fEventNr
 

Static Protected Attributes

static const Double_t kNoChargeSpecified = -999.9
 

Private Member Functions

Bool_t AndMinMaxMom (Double_t min, Double_t max, MomState mom)
 
Bool_t AndMinMaxGeom (Double_t min, Double_t max, GeomState geom)
 
void InitCounters ()
 
void SetDefaultBoundaries ()
 
void CountCharge (TParticle *particle)
 
void CountPdg (TParticle *particle)
 
Bool_t AcceptMomentum (TParticle *particle)
 
Bool_t AcceptGeometry (TParticle *particle)
 
Bool_t AcceptPdgCounter ()
 
Bool_t AcceptChargeCounter ()
 
 ClassDef (FairEvtFilterOnSingleParticleCounts, 1)
 

Private Attributes

PdgGroupId fPdgGroupId
 
std::vector< Int_t > fCountGroupId
 
std::vector< Int_t > fCountCharge
 
std::vector< std::pair< Int_t,
Int_t > > 
fGroupIdCountsMinMax
 
std::vector< std::pair< Int_t,
Int_t > > 
fChargeCountsMinMax
 
std::vector< std::pair
< Double_t, Double_t > > 
fMomMinMax
 
std::vector< std::pair
< Double_t, Double_t > > 
fGeomMinMax
 
Bool_t fFilterPdg
 is kTRUE if any pdg / group ID filter is set More...
 
Bool_t fFilterCharge
 is kTRUE if any filter on electrical charge is set (also kTRUE if filter on the total number of particles is set - charge state kAll) More...
 
Bool_t fFilterMom
 is kTRUE if any momentum filter is set More...
 
Bool_t fFilterGeom
 is kTRUE if any geometry filter is set More...
 

Static Private Attributes

static const Int_t kInvalidPdgCode = 0
 constant holding an integer number which is not used as a pdg code this serves as a place holder when constructing a vector from integers given by the user More...
 

Detailed Description

Definition at line 56 of file FairEvtFilterOnSingleParticleCounts.h.

Member Enumeration Documentation

enum FairEvtFilter::GeomState
inherited
enum FairEvtFilter::LogicOp
inherited
Enumerator
kAnd 
kOr 
kLogicOpLastElement 

Definition at line 42 of file FairEvtFilter.h.

enum FairEvtFilter::MomState
inherited
Enumerator
kMomTotal 
kMomTrans 
kMomZ 
kMomLastElement 

Definition at line 40 of file FairEvtFilter.h.

Constructor & Destructor Documentation

FairEvtFilterOnSingleParticleCounts::FairEvtFilterOnSingleParticleCounts ( )

Definition at line 88 of file FairEvtFilterOnSingleParticleCounts.cxx.

References SetDefaultBoundaries().

88  :
89  FairEvtFilter(" ", "FairEvtFilterOnSingleParticleCounts"), fFilterPdg(kFALSE),fFilterCharge(kFALSE),
90  fFilterMom(kFALSE),fFilterGeom(kFALSE)
91 {
93 }
Bool_t fFilterMom
is kTRUE if any momentum filter is set
Bool_t fFilterGeom
is kTRUE if any geometry filter is set
Bool_t fFilterCharge
is kTRUE if any filter on electrical charge is set (also kTRUE if filter on the total number of parti...
Bool_t fFilterPdg
is kTRUE if any pdg / group ID filter is set
FairEvtFilterOnSingleParticleCounts::FairEvtFilterOnSingleParticleCounts ( const char *  name,
const char *  title = "FairEvtFilterOnSingleParticleCounts" 
)

Definition at line 97 of file FairEvtFilterOnSingleParticleCounts.cxx.

References SetDefaultBoundaries().

98 : FairEvtFilter(name, title),fFilterPdg(kFALSE), fFilterCharge(kFALSE),
99  fFilterMom(kFALSE),fFilterGeom(kFALSE)
100 {
102 }
Bool_t fFilterMom
is kTRUE if any momentum filter is set
Bool_t fFilterGeom
is kTRUE if any geometry filter is set
Bool_t fFilterCharge
is kTRUE if any filter on electrical charge is set (also kTRUE if filter on the total number of parti...
Bool_t fFilterPdg
is kTRUE if any pdg / group ID filter is set
TString name
FairEvtFilterOnSingleParticleCounts::~FairEvtFilterOnSingleParticleCounts ( )
virtual

Definition at line 123 of file FairEvtFilterOnSingleParticleCounts.cxx.

123 {}

Member Function Documentation

Bool_t FairEvtFilterOnSingleParticleCounts::AcceptChargeCounter ( )
private

Definition at line 466 of file FairEvtFilterOnSingleParticleCounts.cxx.

References fChargeCountsMinMax, fCountCharge, FairEvtFilter::fVerbose, and FairEvtFilter::kChargeLastElement.

Referenced by EventMatches().

467 { // check if event matches the filter
468  // event matches if
469  // fCountCharge[icountPdg] >= fChargeCountsMinMax[icountPdg].first
470  // and
471  // fCountCharge[icountPdg] <= fChargeCountsMinMax[icountPdg].second
472  for (UInt_t icountCharge = 0; icountCharge < FairEvtFilter::kChargeLastElement; ++icountCharge){
473  if ( fCountCharge[icountCharge] < fChargeCountsMinMax[icountCharge].first){
474  if (fVerbose >9){
475  std::cout << "Event does not match filter because of fCountCharge[icountChart] < fChargeCountMinMax[icountCharge].first for icountCharge == " << icountCharge << "\n";
476  }
477  return kFALSE;
478  }
479  if ( fCountCharge[icountCharge] > fChargeCountsMinMax[icountCharge].second ){
480  if (fVerbose >9){
481  std::cout << "Event does not match filter because of fCountCharge[icountCharge] > fChargeCountsMinMax[icountCharge].second for icountCharge == " << icountCharge << "\n";
482  }
483  return kFALSE;
484  }
485  }
486  return kTRUE;
487 }
std::vector< std::pair< Int_t, Int_t > > fChargeCountsMinMax
Bool_t FairEvtFilterOnSingleParticleCounts::AcceptGeometry ( TParticle *  particle)
private

Definition at line 390 of file FairEvtFilterOnSingleParticleCounts.cxx.

References Double_t, fGeomMinMax, FairEvtFilter::fVerbose, geom(), FairEvtFilter::kGeomLastElement, FairEvtFilter::kPhi, FairEvtFilter::kTheta, FairEvtFilter::kVertexRadius, FairEvtFilter::kVertexRho, FairEvtFilter::kVertexZ, and Pi.

Referenced by EventMatches().

391 {
392  if (fVerbose > 9){
393  std::cout << "{ Theta , Phi , Vz, VRho, VRadius }: { " << particle->Theta()*180./TMath::Pi() << ", " << particle->Phi()*180./TMath::Pi();
394  std::cout << ", " << particle->Vz() << ", " << particle->Rho() << ", " << particle->R()<< " }\n";
395  }
396 
397  Double_t geom;
398  std::string type;
399  for(UInt_t iGeom=0; iGeom < FairEvtFilter::kGeomLastElement; iGeom++){//check every geometry constraint for the considered particle
400  switch(iGeom){
401  case kTheta:
402  geom=particle->Theta();
403  type="Theta";
404  break;
405  case kPhi:
406  geom=particle->Phi();
407  type="Phi";
408  break;
409  case kVertexZ:
410  geom=particle->Vz();
411  type="Vz";
412  break;
413  case kVertexRho:
414  geom=particle->Rho();
415  type="VRho";
416  break;
417  case kVertexRadius:
418  geom=particle->R();
419  type="VRadius";
420  break;
421  default:
422  std::cout << "\n\n\n -WARNING from FairEvtFilterOnSingleParticleCounts: Check your filter calls!\n\n\n";
423  return kFALSE;
424  }
425  if ( geom < fGeomMinMax[iGeom].first){
426  if (fVerbose >9){
427  std::cout << "Set GeomOk to kFALSE because of " << type << " < And" << type << "Range() minimum." << "\n";
428  }
429  return kFALSE;
430  }
431  if ( geom > fGeomMinMax[iGeom].second ){
432  if (fVerbose >9){
433  std::cout << "Set GeomOk to kFALSE because of " << type << " > And" << type << "Range() maximum." << "\n";
434  }
435  return kFALSE;
436  }
437  }
438  return kTRUE;
439 }
const int particle
std::vector< std::pair< Double_t, Double_t > > fGeomMinMax
Double_t
Double_t Pi
Bool_t FairEvtFilterOnSingleParticleCounts::AcceptMomentum ( TParticle *  particle)
private

Definition at line 347 of file FairEvtFilterOnSingleParticleCounts.cxx.

References Double_t, fMomMinMax, FairEvtFilter::fVerbose, FairEvtFilter::kMomLastElement, FairEvtFilter::kMomTotal, FairEvtFilter::kMomTrans, and FairEvtFilter::kMomZ.

Referenced by EventMatches().

348 {
349  if (fVerbose > 9){
350  std::cout << "{ P , Pt , Pz }: { " << particle->P() << ", " << particle->Pt() << ", " << particle->Pz() << " }\n";
351  }
352 
353  Double_t momentum;
354  std::string type;
355  for(UInt_t iMom=0; iMom < FairEvtFilter::kMomLastElement; iMom++){//check every momentum constraint for the considered particle
356  switch(iMom){
357  case kMomTotal:
358  momentum=particle->P();
359  type="P";
360  break;
361  case kMomTrans:
362  momentum=particle->Pt();
363  type="Pt";
364  break;
365  case kMomZ:
366  momentum=particle->Pz();
367  type="Pz";
368  break;
369  default:
370  std::cout << "\n\n\n -WARNING from FairEvtFilterOnSingleParticleCounts: Check your filter calls!\n\n\n";
371  return kFALSE;
372  if ( momentum < fMomMinMax[iMom].first){
373  if (fVerbose >9){
374  std::cout << "Set MomOk to kFALSE because of " << type << "< And" << type << "Range() minimum." << "\n";
375  }
376  return kFALSE;
377  }
378  if ( momentum > fMomMinMax[iMom].second ){
379  if (fVerbose >9){
380  std::cout << "Set MomOk to kFALSE because of " << type << " > And" << type << "Range() maximum." << "\n";
381  }
382  return kFALSE;
383  }
384  }
385  }
386  return kTRUE;
387 }
const int particle
Double_t
std::vector< std::pair< Double_t, Double_t > > fMomMinMax
Bool_t FairEvtFilterOnSingleParticleCounts::AcceptPdgCounter ( )
private

Definition at line 442 of file FairEvtFilterOnSingleParticleCounts.cxx.

References fCountGroupId, fGroupIdCountsMinMax, and FairEvtFilter::fVerbose.

Referenced by EventMatches().

443 { // check if event matches the filter
444  // event matches if
445  // fCountGroupId[icountPdg] >= fGroupIdCountsMinMax[icountPdg].first
446  // and
447  // fCountGroupId[icountPdg] <= fGroupIdCountsMinMax[icountPdg].second
448  for (UInt_t icountPdg = 0; icountPdg < fCountGroupId.size(); ++icountPdg){
449  if ( fCountGroupId[icountPdg] < fGroupIdCountsMinMax[icountPdg].first){
450  if (fVerbose >9){
451  std::cout << "Event does NOT match filter because of fCountGroupId[icountPdg] < fGroupIdCountMinMax[icountPdg].first for icountPdg == " << icountPdg << "\n";
452  }
453  return kFALSE;
454  }
455  if ( fCountGroupId[icountPdg] > fGroupIdCountsMinMax[icountPdg].second ){
456  if (fVerbose >9){
457  std::cout << "Event does NOT match filter because of fCountGroupId[icountPdg] > fGroupIdCountsMinMax[icountPdg].second for icountPdg == " << icountPdg << "\n";
458  }
459  return kFALSE;
460  }
461  }
462  return kTRUE;
463 }
std::vector< std::pair< Int_t, Int_t > > fGroupIdCountsMinMax
Bool_t FairEvtFilterOnSingleParticleCounts::AndMaxAllParticles ( Int_t  max)
inline

Definition at line 144 of file FairEvtFilterOnSingleParticleCounts.h.

References AndMinMaxCharge(), and FairEvtFilter::kAll.

144 { return AndMinMaxCharge( 0, max, kAll); };
friend F32vec4 max(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:26
Bool_t AndMinMaxCharge(Int_t min, Int_t max, ChargeState charge)
Bool_t FairEvtFilterOnSingleParticleCounts::AndMaxCharge ( Int_t  max,
ChargeState  charge 
)
inline

Definition at line 131 of file FairEvtFilterOnSingleParticleCounts.h.

References AndMinMaxCharge().

131 {return AndMinMaxCharge( 0, max, charge);};
friend F32vec4 max(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:26
Bool_t AndMinMaxCharge(Int_t min, Int_t max, ChargeState charge)
Bool_t FairEvtFilterOnSingleParticleCounts::AndMaxPdgCodes ( Int_t  max,
std::vector< Int_t > &  pdgCodes 
)
inline

Definition at line 106 of file FairEvtFilterOnSingleParticleCounts.h.

References AndMinMaxPdgCodes().

Referenced by prod_fsim().

106  {
107  return AndMinMaxPdgCodes( 0, max, pdgCodes );
108  };
Bool_t AndMinMaxPdgCodes(Int_t min, Int_t max, Int_t pdgCode1, Int_t pdgCode2=kInvalidPdgCode, Int_t pdgCode3=kInvalidPdgCode, Int_t pdgCode4=kInvalidPdgCode, Int_t pdgCode5=kInvalidPdgCode, Int_t pdgCode6=kInvalidPdgCode, Int_t pdgCode7=kInvalidPdgCode, Int_t pdgCode8=kInvalidPdgCode)
friend F32vec4 max(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:26
Bool_t FairEvtFilterOnSingleParticleCounts::AndMaxPdgCodes ( Int_t  max,
Int_t  pdgCode1,
Int_t  pdgCode2 = kInvalidPdgCode,
Int_t  pdgCode3 = kInvalidPdgCode,
Int_t  pdgCode4 = kInvalidPdgCode,
Int_t  pdgCode5 = kInvalidPdgCode,
Int_t  pdgCode6 = kInvalidPdgCode,
Int_t  pdgCode7 = kInvalidPdgCode,
Int_t  pdgCode8 = kInvalidPdgCode 
)
inline

Definition at line 109 of file FairEvtFilterOnSingleParticleCounts.h.

References AndMinMaxPdgCodes().

109  {
110  return AndMinMaxPdgCodes( 0, max, pdgCode1, pdgCode2, pdgCode3, pdgCode4, pdgCode5, pdgCode6, pdgCode7, pdgCode8 );
111  };
Bool_t AndMinMaxPdgCodes(Int_t min, Int_t max, Int_t pdgCode1, Int_t pdgCode2=kInvalidPdgCode, Int_t pdgCode3=kInvalidPdgCode, Int_t pdgCode4=kInvalidPdgCode, Int_t pdgCode5=kInvalidPdgCode, Int_t pdgCode6=kInvalidPdgCode, Int_t pdgCode7=kInvalidPdgCode, Int_t pdgCode8=kInvalidPdgCode)
friend F32vec4 max(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:26
Bool_t FairEvtFilterOnSingleParticleCounts::AndMinAllParticles ( Int_t  min)
inline

Definition at line 143 of file FairEvtFilterOnSingleParticleCounts.h.

References AndMinMaxCharge(), and FairEvtFilter::kAll.

143 { return AndMinMaxCharge( min, 9999, kAll); };
Bool_t AndMinMaxCharge(Int_t min, Int_t max, ChargeState charge)
friend F32vec4 min(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:25
Bool_t FairEvtFilterOnSingleParticleCounts::AndMinCharge ( Int_t  min,
ChargeState  charge 
)
inline

Definition at line 130 of file FairEvtFilterOnSingleParticleCounts.h.

References AndMinMaxCharge().

Referenced by prod_fsim().

130 {return AndMinMaxCharge( min, 9999, charge);};
Bool_t AndMinMaxCharge(Int_t min, Int_t max, ChargeState charge)
friend F32vec4 min(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:25
Bool_t FairEvtFilterOnSingleParticleCounts::AndMinMaxAllParticles ( Int_t  min,
Int_t  max 
)
inline

Definition at line 139 of file FairEvtFilterOnSingleParticleCounts.h.

References AndMinMaxCharge(), and FairEvtFilter::kAll.

Referenced by sim_filter_ex2().

139 { return AndMinMaxCharge(min, max, kAll); };
friend F32vec4 max(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:26
Bool_t AndMinMaxCharge(Int_t min, Int_t max, ChargeState charge)
friend F32vec4 min(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:25
Bool_t FairEvtFilterOnSingleParticleCounts::AndMinMaxCharge ( Int_t  min,
Int_t  max,
ChargeState  charge 
)

Definition at line 218 of file FairEvtFilterOnSingleParticleCounts.cxx.

References fChargeCountsMinMax, fFilterCharge, FairEvtFilter::fVerbose, max(), and min().

Referenced by AndMaxAllParticles(), AndMaxCharge(), AndMinAllParticles(), AndMinCharge(), and AndMinMaxAllParticles().

219 {
220 
221  if (fVerbose > 0){
222  std::cout << this->GetTitle() << ": " << this->GetName() << "\n";
223  std::cout << "fChargeCountsMinMax: " << fChargeCountsMinMax << "\n\n";
224  }
225 
226  // check if filter can be applied (make sure that min and max make sense)
227  if ( min < 0 ){
228  std::cout << "\n\n\n -WARNING from FairEvtFilterOnSingleParticleCounts: Filter could not be added. You are trying to request that your events should have at least a number <= 0 of some particles. That makes no sense. Check your AndMinMaxCharge calls!\n\n\n";
229  return kFALSE;
230  }
231  if ( max < min ){
232  std::cout << "\n\n\n -WARNING from FairEvtFilterOnSingleParticleCounts: Filter could not be added. You are trying to request that your events should have at most a number of some particles which is less than the minimum number that you request. That makes no sense. Check your AndMinMaxCharge calls!\n\n\n";
233  return kFALSE;
234  }
235 
236  // saving the lower and upper limit of particle multiplicities for the corresponding ChargeState in fChargeCountsMinMax
237  fChargeCountsMinMax[charge]=std::pair<Int_t, Int_t> (min,max);
238 
239  // turn filter on
240  fFilterCharge=kTRUE;
241 
242  if (fVerbose > 0){
243  std::cout << "FairEvtFilterOnSingleParticleCounts: After adding ChargeFilter:\n";
244  std::cout << "fGroupIdCountsMinMax: " << fChargeCountsMinMax << "\n\n";
245  }
246  return kTRUE;
247 }
std::vector< std::pair< Int_t, Int_t > > fChargeCountsMinMax
Bool_t fFilterCharge
is kTRUE if any filter on electrical charge is set (also kTRUE if filter on the total number of parti...
friend F32vec4 max(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:26
friend F32vec4 min(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:25
Bool_t FairEvtFilterOnSingleParticleCounts::AndMinMaxGeom ( Double_t  min,
Double_t  max,
FairEvtFilter::GeomState  geom 
)
private

Definition at line 282 of file FairEvtFilterOnSingleParticleCounts.cxx.

References fFilterGeom, fGeomMinMax, FairEvtFilter::fVerbose, geom(), max(), min(), and Pi.

Referenced by AndPhiRange(), AndRadiusRange(), AndRhoRange(), AndThetaRange(), and AndVzRange().

283 {
284 
285  if (fVerbose > 0){
286  std::cout << this->GetTitle() << ": " << this->GetName() << "\n";
287  std::cout << "fGeomMinMax { Theta , Phi, VertexZ, VertexRho, VertexRadius }: " << fGeomMinMax << "\n\n";
288  }
289 
290  // check if filter can be applied (make sure that min and max make sense)
291  if ( min < 0 ){
292  std::cout << "\n\n\n -WARNING from FairEvtFilterOnSingleParticleCounts: Filter could not be added. You are trying to request that your events should have at least a number <= 0 of some particles. Check your AndMinMaxGeom calls.\n\n\n";
293  return kFALSE;
294  }
295  if ( max < min ){
296  std::cout << "\n\n\n -WARNING from FairEvtFilterOnSingleParticleCounts: Filter could not be added. You are trying to request that your events should have at most a number of some particles which is less than the minimum number that you request. Check your AndMinMaxGeom calls.\n\n\n";
297  return kFALSE;
298  }
299 
300  if(geom==0 || geom==1){
301  min=min*2*(TMath::Pi())/360;
302  max=max*2*(TMath::Pi())/360;
303  }
304 
305  // saving the lower and upper limit for the corresponding GeomState in fGeomMinMax
306  fGeomMinMax[geom]=std::pair<Double_t, Double_t> (min,max);
307 
308 
309  // turn filter on
310  fFilterGeom=kTRUE;
311 
312  if (fVerbose > 0){
313  std::cout << "FairEvtFilterOnSingleParticleCounts: After adding GeomFilter:\n";
314  std::cout << "fGeomMinMax { Theta , Phi, VertexZ, VertexRho, VertexRadius }: " << fGeomMinMax << "\n\n";
315  }
316 
317  return kTRUE;
318 }
Bool_t fFilterGeom
is kTRUE if any geometry filter is set
friend F32vec4 max(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:26
std::vector< std::pair< Double_t, Double_t > > fGeomMinMax
friend F32vec4 min(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:25
Double_t Pi
Bool_t FairEvtFilterOnSingleParticleCounts::AndMinMaxMom ( Double_t  min,
Double_t  max,
FairEvtFilter::MomState  mom 
)
private

Definition at line 250 of file FairEvtFilterOnSingleParticleCounts.cxx.

References fFilterMom, fMomMinMax, FairEvtFilter::fVerbose, max(), min(), and mom.

Referenced by AndPRange(), AndPtRange(), and AndPzRange().

251 {
252 
253  if (fVerbose > 0){
254  std::cout << this->GetTitle() << ": " << this->GetName() << "\n";
255  std::cout << "fMomMinMax { P , Pt , Pz }: " << fMomMinMax << "\n\n";
256  }
257 
258  // check if filter can be applied (make sure that min and max make sense)
259  if ( min < 0 ){
260  std::cout << "\n\n\n -WARNING from FairEvtFilterOnSingleParticleCounts: Filter could not be added. You are trying to request that your events should have at least a number <= 0 of some particles. Check your AndMinMaxMom calls.\n\n\n";
261  return kFALSE;
262  }
263  if ( max < min ){
264  std::cout << "\n\n\n -WARNING from FairEvtFilterOnSingleParticleCounts: Filter could not be added. You are trying to request that your events should have at most a number of some particles which is less than the minimum number that you request. Check your AndMinMaxMom calls.\n\n\n";
265  return kFALSE;
266  }
267 
268  // saving the lower and upper limit for the corresponding MomState in fMomMinMax
269  fMomMinMax[mom]=std::pair<Double_t, Double_t> (min,max);
270 
271  // turn filter on
272  fFilterMom=kTRUE;
273 
274  if (fVerbose > 0){
275  std::cout << "FairEvtFilterOnSingleParticleCounts: After adding MomFilter:\n";
276  std::cout << "fMomMinMax { P , Pt , Pz }: " << fMomMinMax << "\n\n";
277  }
278  return kTRUE;
279 }
Bool_t fFilterMom
is kTRUE if any momentum filter is set
Double_t mom
Definition: plot_dirc.C:14
friend F32vec4 max(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:26
friend F32vec4 min(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:25
std::vector< std::pair< Double_t, Double_t > > fMomMinMax
Bool_t FairEvtFilterOnSingleParticleCounts::AndMinMaxPdgCodes ( Int_t  min,
Int_t  max,
Int_t  pdgCode1,
Int_t  pdgCode2 = kInvalidPdgCode,
Int_t  pdgCode3 = kInvalidPdgCode,
Int_t  pdgCode4 = kInvalidPdgCode,
Int_t  pdgCode5 = kInvalidPdgCode,
Int_t  pdgCode6 = kInvalidPdgCode,
Int_t  pdgCode7 = kInvalidPdgCode,
Int_t  pdgCode8 = kInvalidPdgCode 
)

Definition at line 126 of file FairEvtFilterOnSingleParticleCounts.cxx.

Referenced by AndMaxPdgCodes(), AndMinPdgCodes(), and sim_filter_ex2().

127 {
128  // construct a vector pdgCodes containing all pdg codes given by caller
129  // and call the general AndMinMaxPdgCodes method
130  std::vector<Int_t> pdgCodes;
131  pdgCodes.push_back(pdgCode1);
132  pdgCodes.push_back(pdgCode2);
133  pdgCodes.push_back(pdgCode3);
134  pdgCodes.push_back(pdgCode4);
135  pdgCodes.push_back(pdgCode5);
136  pdgCodes.push_back(pdgCode6);
137  pdgCodes.push_back(pdgCode7);
138  pdgCodes.push_back(pdgCode8);
139 
140  return AndMinMaxPdgCodes( min, max, pdgCodes );
141 }
Bool_t AndMinMaxPdgCodes(Int_t min, Int_t max, Int_t pdgCode1, Int_t pdgCode2=kInvalidPdgCode, Int_t pdgCode3=kInvalidPdgCode, Int_t pdgCode4=kInvalidPdgCode, Int_t pdgCode5=kInvalidPdgCode, Int_t pdgCode6=kInvalidPdgCode, Int_t pdgCode7=kInvalidPdgCode, Int_t pdgCode8=kInvalidPdgCode)
friend F32vec4 max(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:26
friend F32vec4 min(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:25
Bool_t FairEvtFilterOnSingleParticleCounts::AndMinMaxPdgCodes ( Int_t  min,
Int_t  max,
std::vector< Int_t > &  pdgCodes 
)

Definition at line 144 of file FairEvtFilterOnSingleParticleCounts.cxx.

References Bool_t, fFilterPdg, fGroupIdCountsMinMax, finder, fPdgGroupId, FairEvtFilter::fVerbose, and kInvalidPdgCode.

145 {
146  // was adding of the filter successful
147  Bool_t filterAdded = kFALSE;
148 
149  if (fVerbose > 0){
150  std::cout << this->GetTitle() << ": " << this->GetName() << "\n";
151  std::cout << "fPdgGroupId: " << fPdgGroupId << "\n";
152  std::cout << "fGroupIdCountsMinMax: " << fGroupIdCountsMinMax << "\n\n";
153  }
154 
155 
156  // check if filter can be applied (make sure that min and max make sense)
157  if ( min < 0 ){
158  std::cout << "\n\n\n -WARNING from FairEvtFilterOnSingleParticleCounts: Filter could not be added. You are trying to request that your events should have at least a number <= 0 of some particles. That makes no sense. Check your AndMinMaxPdgCodes calls!\n\n\n";
159  return kFALSE;
160  }
161  if ( max < min ){
162  std::cout << "\n\n\n -WARNING from FairEvtFilterOnSingleParticleCounts: Filter could not be added. You are trying to request that your events should have at most a number of some particles which is less than the minimum number that you request. That makes no sense. Check your AndMinMaxPdgCodes calls!\n\n\n";
163  return kFALSE;
164  }
165 
166  // find the groupId for the pdgCodes
167  // the current pdgCodes are appended to the pre-existing vector
168  Int_t newGroupId = fGroupIdCountsMinMax.size();
169 
170  // for every entry in pdgCodes
171  // -- save a new entry in the map fPdgGroupId if the entry is not already in map
172  // -- add a new groupId to the vector in fPdgGroupId if the entry already exists in the map
173 
174  for (UInt_t iPdgCodes = 0; iPdgCodes < pdgCodes.size(); ++iPdgCodes){
175  // skip kInvalidPdgCode entries in pdgCodes
176  if ( kInvalidPdgCode != pdgCodes[iPdgCodes] ) {
177  // check if map has already an entry for the pdgCode
178  PdgGroupIdIterator finder = fPdgGroupId.find(pdgCodes[iPdgCodes]);
179  if ( fPdgGroupId.end() != finder){
180  // pdgCode is already in map, add groupId to vector
181 
182  std::vector<Int_t> & groupIdVector = fPdgGroupId[pdgCodes[iPdgCodes]];
183  groupIdVector.push_back(newGroupId);
184  }
185  else {
186  // pdgCode is not in map yet, create new entry of Int_t and vector<Int_t>
187 
188  std::vector<Int_t> groupIdVector;
189  groupIdVector.push_back(newGroupId);
190  fPdgGroupId.insert(PdgGroupIdPair(pdgCodes[iPdgCodes],groupIdVector));
191  }
192 
193  filterAdded = kTRUE;
194  }
195  }
196 
197  if ( kTRUE == filterAdded ){
198  // set the minimum and maximum number of particles for the current groupId
199  fGroupIdCountsMinMax.push_back(std::pair<Int_t, Int_t> (min,max));
200 
201  } else {
202  std::cout << "\n\n\n -WARNING from FairEvtFilterOnSingleParticleCounts: Tried to add filter, but that was not possible as all pdg codes were kInvalidPdgCode!\n\n\n";
203  return kFALSE;
204  }
205 
206  // turn filter on
207  fFilterPdg = kTRUE;
208 
209  if (fVerbose > 0){
210  std::cout << "FairEvtFilterOnSingleParticleCounts: After adding PdgFilter:\n";
211  std::cout << "fPdgGroupId: " << fPdgGroupId << "\n";
212  std::cout << "fGroupIdCountsMinMax: " << fGroupIdCountsMinMax << "\n\n";
213  }
214  return kTRUE;
215 }
static const Int_t kInvalidPdgCode
constant holding an integer number which is not used as a pdg code this serves as a place holder when...
std::map< Int_t, std::vector< Int_t > >::iterator PdgGroupIdIterator
friend F32vec4 max(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:26
PndStraightLineTrackFinderTask * finder
std::vector< std::pair< Int_t, Int_t > > fGroupIdCountsMinMax
Bool_t fFilterPdg
is kTRUE if any pdg / group ID filter is set
friend F32vec4 min(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:25
std::pair< Int_t, std::vector< Int_t > > PdgGroupIdPair
Bool_t FairEvtFilterOnSingleParticleCounts::AndMinPdgCodes ( Int_t  min,
Int_t  pdgCode1,
Int_t  pdgCode2 = kInvalidPdgCode,
Int_t  pdgCode3 = kInvalidPdgCode,
Int_t  pdgCode4 = kInvalidPdgCode,
Int_t  pdgCode5 = kInvalidPdgCode,
Int_t  pdgCode6 = kInvalidPdgCode,
Int_t  pdgCode7 = kInvalidPdgCode,
Int_t  pdgCode8 = kInvalidPdgCode 
)
inline

Definition at line 95 of file FairEvtFilterOnSingleParticleCounts.h.

References AndMinMaxPdgCodes().

Referenced by quickfsimana().

95  {
96  return AndMinMaxPdgCodes( min, 9999, pdgCode1, pdgCode2, pdgCode3, pdgCode4, pdgCode5, pdgCode6, pdgCode7, pdgCode8 );
97  };
Bool_t AndMinMaxPdgCodes(Int_t min, Int_t max, Int_t pdgCode1, Int_t pdgCode2=kInvalidPdgCode, Int_t pdgCode3=kInvalidPdgCode, Int_t pdgCode4=kInvalidPdgCode, Int_t pdgCode5=kInvalidPdgCode, Int_t pdgCode6=kInvalidPdgCode, Int_t pdgCode7=kInvalidPdgCode, Int_t pdgCode8=kInvalidPdgCode)
friend F32vec4 min(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:25
Bool_t FairEvtFilterOnSingleParticleCounts::AndMinPdgCodes ( Int_t  min,
std::vector< Int_t > &  pdgCodes 
)
inline

Definition at line 98 of file FairEvtFilterOnSingleParticleCounts.h.

References AndMinMaxPdgCodes().

98  {
99  return AndMinMaxPdgCodes( min, 9999, pdgCodes );
100  };
Bool_t AndMinMaxPdgCodes(Int_t min, Int_t max, Int_t pdgCode1, Int_t pdgCode2=kInvalidPdgCode, Int_t pdgCode3=kInvalidPdgCode, Int_t pdgCode4=kInvalidPdgCode, Int_t pdgCode5=kInvalidPdgCode, Int_t pdgCode6=kInvalidPdgCode, Int_t pdgCode7=kInvalidPdgCode, Int_t pdgCode8=kInvalidPdgCode)
friend F32vec4 min(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:25
Bool_t FairEvtFilterOnSingleParticleCounts::AndPhiRange ( Double_t  min,
Double_t  max 
)
inline

Definition at line 163 of file FairEvtFilterOnSingleParticleCounts.h.

References AndMinMaxGeom(), and FairEvtFilter::kPhi.

163 {return AndMinMaxGeom( min, max, kPhi);};
friend F32vec4 max(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:26
Bool_t AndMinMaxGeom(Double_t min, Double_t max, GeomState geom)
friend F32vec4 min(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:25
Bool_t FairEvtFilterOnSingleParticleCounts::AndPRange ( Double_t  min,
Double_t  max 
)
inline

Definition at line 152 of file FairEvtFilterOnSingleParticleCounts.h.

References AndMinMaxMom(), and FairEvtFilter::kMomTotal.

152 {return AndMinMaxMom( min, max, kMomTotal);};
friend F32vec4 max(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:26
Bool_t AndMinMaxMom(Double_t min, Double_t max, MomState mom)
friend F32vec4 min(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:25
Bool_t FairEvtFilterOnSingleParticleCounts::AndPtRange ( Double_t  min,
Double_t  max 
)
inline

Definition at line 153 of file FairEvtFilterOnSingleParticleCounts.h.

References AndMinMaxMom(), and FairEvtFilter::kMomTrans.

153 {return AndMinMaxMom( min, max, kMomTrans);};
friend F32vec4 max(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:26
Bool_t AndMinMaxMom(Double_t min, Double_t max, MomState mom)
friend F32vec4 min(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:25
Bool_t FairEvtFilterOnSingleParticleCounts::AndPzRange ( Double_t  min,
Double_t  max 
)
inline

Definition at line 154 of file FairEvtFilterOnSingleParticleCounts.h.

References AndMinMaxMom(), and FairEvtFilter::kMomZ.

154 {return AndMinMaxMom( min, max, kMomZ);};
friend F32vec4 max(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:26
Bool_t AndMinMaxMom(Double_t min, Double_t max, MomState mom)
friend F32vec4 min(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:25
Bool_t FairEvtFilterOnSingleParticleCounts::AndRadiusRange ( Double_t  min,
Double_t  max 
)
inline

Definition at line 173 of file FairEvtFilterOnSingleParticleCounts.h.

References AndMinMaxGeom(), and FairEvtFilter::kVertexRadius.

173 {return AndMinMaxGeom( min, max, kVertexRadius);};
friend F32vec4 max(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:26
Bool_t AndMinMaxGeom(Double_t min, Double_t max, GeomState geom)
friend F32vec4 min(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:25
Bool_t FairEvtFilterOnSingleParticleCounts::AndRhoRange ( Double_t  min,
Double_t  max 
)
inline

Definition at line 172 of file FairEvtFilterOnSingleParticleCounts.h.

References AndMinMaxGeom(), and FairEvtFilter::kVertexRho.

172 {return AndMinMaxGeom( min, max, kVertexRho);};
friend F32vec4 max(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:26
Bool_t AndMinMaxGeom(Double_t min, Double_t max, GeomState geom)
friend F32vec4 min(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:25
Bool_t FairEvtFilterOnSingleParticleCounts::AndThetaRange ( Double_t  min,
Double_t  max 
)
inline

Definition at line 162 of file FairEvtFilterOnSingleParticleCounts.h.

References AndMinMaxGeom(), and FairEvtFilter::kTheta.

162 {return AndMinMaxGeom( min, max, kTheta);};
friend F32vec4 max(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:26
Bool_t AndMinMaxGeom(Double_t min, Double_t max, GeomState geom)
friend F32vec4 min(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:25
Bool_t FairEvtFilterOnSingleParticleCounts::AndVzRange ( Double_t  min,
Double_t  max 
)
inline

Definition at line 171 of file FairEvtFilterOnSingleParticleCounts.h.

References AndMinMaxGeom(), and FairEvtFilter::kVertexZ.

171 {return AndMinMaxGeom( min, max, kVertexZ);};
friend F32vec4 max(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:26
Bool_t AndMinMaxGeom(Double_t min, Double_t max, GeomState geom)
friend F32vec4 min(const F32vec4 &a, const F32vec4 &b)
Definition: P4_F32vec4.h:25
FairEvtFilterOnSingleParticleCounts::ClassDef ( FairEvtFilterOnSingleParticleCounts  ,
 
)
private
void FairEvtFilterOnSingleParticleCounts::CountCharge ( TParticle *  particle)
private

Definition at line 517 of file FairEvtFilterOnSingleParticleCounts.cxx.

References Double_t, fCountCharge, FairEvtFilter::fVerbose, FairEvtFilter::GetCharge(), FairEvtFilter::kAll, FairEvtFilter::kCharged, FairEvtFilter::kMinus, FairEvtFilter::kNeutral, and FairEvtFilter::kPlus.

Referenced by EventMatches().

517  {
518  // get pdg code for particle
519  Int_t partPdg = particle->GetPdgCode();
520 
521 
522  // get charge for pdg code
523  Double_t pdgCodeCharge = 0;
524  if ( kFALSE == GetCharge( partPdg, &pdgCodeCharge ) ) {
525  // skip particles of unknown charge
526  return;
527  }
528 
529 
530  // Based on charge increment the according entry in fChargeCountsMinMax
531  if(pdgCodeCharge!=0){
533  if(pdgCodeCharge>0){
535  }else{
537  }
538  }else{
540  }
542 
543  if (fVerbose >9){
544  std::cout << "Charge:" << pdgCodeCharge << " -> " << "fCountCharge: " << fCountCharge << "\n";
545  }
546 }
Bool_t GetCharge(Int_t inPdgCode, Double_t *pdgCodeCharge)
const int particle
Double_t
void FairEvtFilterOnSingleParticleCounts::CountPdg ( TParticle *  particle)
private

Definition at line 490 of file FairEvtFilterOnSingleParticleCounts.cxx.

References fCountGroupId, fPdgGroupId, and FairEvtFilter::fVerbose.

Referenced by EventMatches().

490  {
491  // every restriction that is placed on a pdg code is described by a groupId
492  // one groupId can belong to several particles assuming that they are indistinguishable
493  // fPdgGroupId contains all the user defined constraints on particles with certain pdg codes organized in a map
494  // all elements in this map have the pdg code as first entry and a vector of groupId's as second entry
495  // if there exists an entry for the considered particle in the map the appropriate entry of fCountGroupId is counted up
496  // whereat the associated groupId's are exactly representing the index of the incremented fCountGroupId vector
497 
498  // get pdg code for particle
499  Int_t partPdg = particle->GetPdgCode();
500 
501  // check if the pdg code is contained in some groupId
502  PdgGroupIdIterator it = fPdgGroupId.find(partPdg);
503  if ( fPdgGroupId.end() != it ){
504  // pdg code was found
505  std::vector<Int_t> groupIdVector = it->second;
506  for (size_t iGroupId = 0; iGroupId < groupIdVector.size(); ++iGroupId){
507  ++fCountGroupId[groupIdVector[iGroupId]];
508  }
509  }
510 
511  if (fVerbose >9){
512  std::cout << "partPdg:" << partPdg << " -> " << "fCountGroupId: " << fCountGroupId << "\n";
513  }
514 }
std::map< Int_t, std::vector< Int_t > >::iterator PdgGroupIdIterator
const int particle
Bool_t FairEvtFilterOnSingleParticleCounts::EventMatches ( Int_t  evtNr)
virtual

Implements FairEvtFilter.

Definition at line 549 of file FairEvtFilterOnSingleParticleCounts.cxx.

References AcceptChargeCounter(), AcceptGeometry(), AcceptMomentum(), AcceptPdgCounter(), Bool_t, CountCharge(), CountPdg(), FairEvtFilter::fAcceptedEventNumbers, fChargeCountsMinMax, fCountCharge, fCountGroupId, fFilterCharge, fFilterGeom, fFilterMom, fFilterPdg, fGeomMinMax, fGroupIdCountsMinMax, fMomMinMax, FairEvtFilter::fParticleList, fPdgGroupId, FairEvtFilter::fVerbose, InitCounters(), particle, and FairEvtFilter::PrintAllTParticleInEvent().

550 {
551  if (fVerbose > 3){
552  std::cout << "\n\n";
553  std::cout << "Generated event: " << evtNr << "\n";
554  std::cout << "FairEvtFilterOnSingleParticleCounts: " << this->GetTitle() << ": " << this->GetName() << " Beginning of EventMatches\n";
555  std::cout << "Nr. of simulated particles: " << fParticleList->GetEntries()<<"\n";
557  std::cout << "fPdgGroupId: " << fPdgGroupId << "\n";
558  std::cout << "fGroupIdCountsMinMax: " << fGroupIdCountsMinMax << "\n";
559  std::cout << "fChargeCountsMinMax: " << fChargeCountsMinMax << "\n";
560  std::cout << "fMomMinMax: " << fMomMinMax << "\n";
561  std::cout << "fGeomMinMax: " << fGeomMinMax << "\n";
562  }
563 
564 
565  // Not necessary as FairFilteredPrimaryGenerator already checks this
566  // // event matches if filter is not set
567  // if ( kFALSE == FilterActive() ) {
568  // return kTRUE;
569  // }
570 
571 
572  // sanity checks
573  if (0==fParticleList){
574  std::cout << "\n\n\n FairEvtFilterOnSingleParticleCounts: FATAL ERROR: No particle list! Event does not match.\n\n\n";
575  return kFALSE;
576  }
577 
578  if (0==fParticleList->GetEntriesFast()){
579  std::cout << "\n\n\n FairEvtFilterOnSingleParticleCounts: Event contains 0 particles. Event does not match.\n\n\n";
580  return kFALSE;
581  }
582 
583  // reset all counters for each new event
584  InitCounters();
585 
586 
587  // loop the produced particles and count up
588  // the specified "groupIDs"
589  // and charge states
590  // for particles which satisfy the momentum, angular and vertex constraints
591  for (Int_t iPart= 0; iPart< fParticleList->GetEntries(); iPart++) {
592 
593  TParticle* particle=(TParticle*)fParticleList->At(iPart);
594  if (0 == particle) { continue; }
595 
596  // check whether the particle satisfies the momentum constraints
597  if(fFilterMom){
598  if(!AcceptMomentum(particle)){ continue; } // if the particle doesn't suit the momentum filter settings go ahead to the next particle
599  }
600 
601  // check whether the particle satisfies the geometric constraints
602  if(fFilterGeom){
603  if(!AcceptGeometry(particle)){ continue; } // if the particle doesn't suit the geometric filter settings go ahead to the next particle
604  }
605 
606  // count up the appropriate element of fCountGroupId and fCountCharge
607  // important: Do not count before all of the above tests were passed!
608  if(fFilterPdg){
609  CountPdg(particle);
610  }
611  if(fFilterCharge){
612  CountCharge(particle);
613  }
614 
615  }// end of the loop over all particles
616 
617  // reset PdgOk and ChargeOk to kTRUE for each new event
618  // for checking if event matches
619  Bool_t PdgOk= kTRUE;
620  Bool_t ChargeOk = kTRUE;
621 
622  if (fVerbose >9){
623  std::cout << "\n Check if event matches the filter \n";
624  std::cout << "fCountGroupId: " << fCountGroupId << "\n";
625  std::cout << "fGroupIdCountsMinMax: " << fGroupIdCountsMinMax << "\n";
626  std::cout << "fCountCharge: " << fCountCharge << "\n";
627  std::cout << "fChargeCountsMinMax: " << fChargeCountsMinMax << "\n";
628  }
629 
630 
631  if(fFilterPdg){
632  PdgOk=AcceptPdgCounter();
633  }
634 
635  if(fFilterCharge){
636  ChargeOk=AcceptChargeCounter();
637  }
638 
639  // event matches if pdg and charges are ok
640  Bool_t evtOk=PdgOk && ChargeOk;
641 
642  if ( kTRUE == evtOk ) {
643  fAcceptedEventNumbers.insert(evtNr); // for QA
644  }
645 
646  if (fVerbose >5){
647  if ( kTRUE == evtOk ) {
648  std::cout << "\n Event matches " << this->GetTitle() << ": " << this->GetName() << "\n\n";
649  } else {
650  std::cout << "\n Event does NOT match " << this->GetTitle() << ": " << this->GetName() << "\n\n";
651  }
652  }
653  return evtOk;
654 }
Bool_t fFilterMom
is kTRUE if any momentum filter is set
void PrintAllTParticleInEvent()
std::vector< std::pair< Int_t, Int_t > > fChargeCountsMinMax
Bool_t fFilterGeom
is kTRUE if any geometry filter is set
Bool_t fFilterCharge
is kTRUE if any filter on electrical charge is set (also kTRUE if filter on the total number of parti...
const int particle
std::vector< std::pair< Double_t, Double_t > > fGeomMinMax
std::vector< std::pair< Int_t, Int_t > > fGroupIdCountsMinMax
std::set< Int_t > fAcceptedEventNumbers
TClonesArray * fParticleList
Bool_t fFilterPdg
is kTRUE if any pdg / group ID filter is set
std::vector< std::pair< Double_t, Double_t > > fMomMinMax
Bool_t FairEvtFilterOnSingleParticleCounts::FilterActive ( )
inlinevirtual

Implements FairEvtFilter.

Definition at line 183 of file FairEvtFilterOnSingleParticleCounts.h.

References fFilterCharge, and fFilterPdg.

183  {
184  return ( fFilterPdg || fFilterCharge );
185  };
Bool_t fFilterCharge
is kTRUE if any filter on electrical charge is set (also kTRUE if filter on the total number of parti...
Bool_t fFilterPdg
is kTRUE if any pdg / group ID filter is set
Bool_t FairEvtFilter::GetCharge ( Int_t  inPdgCode,
Double_t pdgCodeCharge 
)
inherited

Definition at line 62 of file FairEvtFilter.cxx.

References FairEvtFilter::fdbPdg, FairEvtFilter::fVerbose, and FairEvtFilter::kNoChargeSpecified.

Referenced by CountCharge(), PndEvtFilter::FillList(), and PndEvtFilterOnInvMassCounts::SetPdgCodesToCombine().

63 {
64  // Try to find the pdg code
65  TParticlePDG *ptrToPdg = fdbPdg->GetParticle(inPdgCode);
66  if ( 0 == ptrToPdg) {
67  // ignore particles with unknown charges
68  std::cout << "WARNING from FairEvtFilter::GetCharge Charge of pdgCode " << inPdgCode << " is unknown and will be ignored!\n";
69  *pdgCodeCharge = kNoChargeSpecified;
70  return kFALSE;
71  }
72  *pdgCodeCharge = ptrToPdg->Charge()/3.; // TParticlePDG contains charge in units of |e|/3
73  if ( fVerbose > 1 ) std::cout << "Found pdgCodeCharge = " << *pdgCodeCharge << " for inPdgCode " << inPdgCode << '\n';
74  return kTRUE;
75 }
static const Double_t kNoChargeSpecified
TDatabasePDG * fdbPdg
Bool_t FairEvtFilter::Init ( )
inlineinherited

Definition at line 64 of file FairEvtFilter.h.

Referenced by FairFilteredPrimaryGenerator::Init().

64 { return kTRUE;}
void FairEvtFilterOnSingleParticleCounts::InitCounters ( )
private

Definition at line 321 of file FairEvtFilterOnSingleParticleCounts.cxx.

References fCountCharge, fCountGroupId, fGroupIdCountsMinMax, FairEvtFilter::fVerbose, and FairEvtFilter::kChargeLastElement.

Referenced by EventMatches().

322 {
323  // initialize a 0 fCountGroupId of same size as fGroupIdCountsMinMax
324  fCountGroupId.clear();
325  for (UInt_t icountPdg = 0; icountPdg < fGroupIdCountsMinMax.size(); ++icountPdg){
326  fCountGroupId.push_back(0);
327  }
328 
329  // fCountGroupId counts up, event matches if all entries in fCountGroupId are between the according entries in fGroupIdCountsMinMax
330  if (fVerbose >11){
331  std::cout << "fCountGroupId after initialization: " << fCountGroupId << "\n";
332  }
333 
334  // initialize a 0 fCountCharge of same size as fChargeCountsMinMax
335  fCountCharge.clear();
336  for (UInt_t icountCharge = 0; icountCharge < FairEvtFilter::kChargeLastElement; ++icountCharge){
337  fCountCharge.push_back(0);
338  }
339 
340  // fCountCharge counts up, event matches if all entries in fCountCharge are between the according entries in fChargeCountsMinMax
341  if (fVerbose >11){
342  std::cout << "fCountCharge after initialization: " << fCountCharge << "\n";
343  }
344 }
std::vector< std::pair< Int_t, Int_t > > fGroupIdCountsMinMax
void FairEvtFilter::PrintAllTParticleInEvent ( )
inherited

Definition at line 52 of file FairEvtFilter.cxx.

References FairEvtFilter::fParticleList, and particle.

Referenced by PndEvtFilterOnInvMassCounts::EventMatches(), and EventMatches().

52  {
53  // sanity checks
54  if (0==fParticleList){ return; }
55  for (Int_t iPart=0; iPart<fParticleList->GetEntries(); ++iPart) {
56  TParticle *particle = (TParticle*)fParticleList->At(iPart);
57  particle->Print();
58  }
59 }
const int particle
TClonesArray * fParticleList
void FairEvtFilterOnSingleParticleCounts::SetDefaultBoundaries ( )
private

Definition at line 105 of file FairEvtFilterOnSingleParticleCounts.cxx.

References fChargeCountsMinMax, fGeomMinMax, fMomMinMax, FairEvtFilter::kChargeLastElement, FairEvtFilter::kGeomLastElement, and FairEvtFilter::kMomLastElement.

Referenced by FairEvtFilterOnSingleParticleCounts().

105  {
106  // initialize a default fChargeCountsMinMax with kChargeLastElement entries
107  for (UInt_t icount = 0; icount < FairEvtFilter::kChargeLastElement; ++icount){
108  fChargeCountsMinMax.push_back(std::pair<Int_t, Int_t> (0,99999));
109  }
110 
111  // initialize a default fMomMinMax with kMomLastElement entries
112  for (UInt_t icount = 0; icount < FairEvtFilter::kMomLastElement; ++icount){
113  fMomMinMax.push_back(std::pair<Double_t, Double_t> (-99999.,99999.));
114  }
115 
116  // initialize a default fGeomMinMax with kGeomLastElement entries
117  for (UInt_t icount = 0; icount < FairEvtFilter::kGeomLastElement; ++icount){
118  fGeomMinMax.push_back(std::pair<Double_t, Double_t> (-99999.,99999.));
119  }
120 }
std::vector< std::pair< Int_t, Int_t > > fChargeCountsMinMax
std::vector< std::pair< Double_t, Double_t > > fGeomMinMax
std::vector< std::pair< Double_t, Double_t > > fMomMinMax
Bool_t FairEvtFilter::SetListOfParticles ( TClonesArray *  ParticleList)
inlineinherited

Definition at line 57 of file FairEvtFilter.h.

References FairEvtFilter::fParticleList.

Referenced by FairFilteredPrimaryGenerator::GenerateEvent().

57 {fParticleList=ParticleList; return kTRUE;};
TClonesArray * fParticleList
void FairEvtFilter::SetTestMode ( Int_t *  arrayPtr,
Int_t  nEntries 
)
inlineinherited

Definition at line 91 of file FairEvtFilter.h.

References FairEvtFilter::fEvtNrsToAccept, and FairEvtFilter::fTestMode.

91  {
92  //turns on the test mode with the declared fEvtNrsToAccept
93  fTestMode=kTRUE;
94  std::set<Int_t> evtNrsToAccept (arrayPtr,arrayPtr+nEntries);
95  fEvtNrsToAccept=evtNrsToAccept;
96  }
std::set< Int_t > fEvtNrsToAccept
void FairEvtFilter::SetVerbose ( Int_t  verbose = 12)
inlineinherited

Definition at line 67 of file FairEvtFilter.h.

References FairEvtFilter::fVerbose, and verbose.

Referenced by sim_filter_inv_mass().

67  {
68  if ( verbose >= 0 ){
69  fVerbose = verbose;
70  std::cout << "FairEvtFilter: fVerbose is now set to " << verbose << "\n";
71  } else {
72  std::cout << "\n\n\n -WARNING from FairEvtFilter: verbose must be a positive number! Check your SetVerbose call!\n\n\n";
73  }
74  }
#define verbose
void FairEvtFilter::ShowAcceptedEventNumbers ( )
inlineinherited

Definition at line 77 of file FairEvtFilter.h.

References FairEvtFilter::fAcceptedEventNumbers.

77  {
78  // for QA
79  // shows fAcceptedEventNumbers that is filled after running a simulation
80  std::cout << "\n fAcceptedEventNumbers" << " = " << fAcceptedEventNumbers;
81  }
std::set< Int_t > fAcceptedEventNumbers
void FairEvtFilter::ShowEvtNrsToAccept ( )
inlineinherited

Definition at line 84 of file FairEvtFilter.h.

References FairEvtFilter::fEvtNrsToAccept.

84  {
85  // for QA
86  //shows fEvtNrsToAccept that has to be set if you want to run your simulation in test mode
87  std::cout << "\n fEvtNrsToAccept" << " = " << fEvtNrsToAccept;
88  }
std::set< Int_t > fEvtNrsToAccept
Bool_t FairEvtFilter::TestPassed ( )
inlineinherited

Definition at line 99 of file FairEvtFilter.h.

References FairEvtFilter::fAcceptedEventNumbers, FairEvtFilter::fEvtNrsToAccept, and FairEvtFilter::fTestMode.

99  {
100  if(kFALSE==fTestMode){
101  std::cout << "\n\n\n WARNING from FairEvtFilter: Test mode not set.\n\n\n";
102  return kFALSE;
103  }
105  //std::cout << "\n\n\n FairEvtFilter: Test passed.\n\n\n";
106  return kTRUE;
107  }else{
108  //std::cout << "\n\n\n FairEvtFilter: Test failed. Check your SetTestMode calls. \n\n\n";
109  return kFALSE;
110  }
111  }
std::set< Int_t > fAcceptedEventNumbers
std::set< Int_t > fEvtNrsToAccept

Member Data Documentation

std::set<Int_t> FairEvtFilter::fAcceptedEventNumbers
protectedinherited
std::vector< std::pair<Int_t,Int_t> > FairEvtFilterOnSingleParticleCounts::fChargeCountsMinMax
private
std::vector<Int_t> FairEvtFilterOnSingleParticleCounts::fCountCharge
private
std::vector<Int_t> FairEvtFilterOnSingleParticleCounts::fCountGroupId
private
TDatabasePDG* FairEvtFilter::fdbPdg
protectedinherited
Int_t FairEvtFilter::fEventNr
protectedinherited

Definition at line 136 of file FairEvtFilter.h.

std::set<Int_t> FairEvtFilter::fEvtNrsToAccept
protectedinherited
Bool_t FairEvtFilterOnSingleParticleCounts::fFilterCharge
private

is kTRUE if any filter on electrical charge is set (also kTRUE if filter on the total number of particles is set - charge state kAll)

Definition at line 250 of file FairEvtFilterOnSingleParticleCounts.h.

Referenced by AndMinMaxCharge(), EventMatches(), and FilterActive().

Bool_t FairEvtFilterOnSingleParticleCounts::fFilterGeom
private

is kTRUE if any geometry filter is set

Definition at line 252 of file FairEvtFilterOnSingleParticleCounts.h.

Referenced by AndMinMaxGeom(), and EventMatches().

Bool_t FairEvtFilterOnSingleParticleCounts::fFilterMom
private

is kTRUE if any momentum filter is set

Definition at line 251 of file FairEvtFilterOnSingleParticleCounts.h.

Referenced by AndMinMaxMom(), and EventMatches().

Bool_t FairEvtFilterOnSingleParticleCounts::fFilterPdg
private

is kTRUE if any pdg / group ID filter is set

Definition at line 249 of file FairEvtFilterOnSingleParticleCounts.h.

Referenced by AndMinMaxPdgCodes(), EventMatches(), and FilterActive().

std::vector< std::pair<Double_t,Double_t> > FairEvtFilterOnSingleParticleCounts::fGeomMinMax
private
std::vector< std::pair<Int_t,Int_t> > FairEvtFilterOnSingleParticleCounts::fGroupIdCountsMinMax
private
std::vector< std::pair<Double_t,Double_t> > FairEvtFilterOnSingleParticleCounts::fMomMinMax
private
TClonesArray* FairEvtFilter::fParticleList
protectedinherited
PdgGroupId FairEvtFilterOnSingleParticleCounts::fPdgGroupId
private
Bool_t FairEvtFilter::fTestMode
protectedinherited

Definition at line 133 of file FairEvtFilter.h.

Referenced by FairEvtFilter::SetTestMode(), and FairEvtFilter::TestPassed().

Int_t FairEvtFilter::fVerbose
protectedinherited
const Int_t FairEvtFilterOnSingleParticleCounts::kInvalidPdgCode = 0
staticprivate

constant holding an integer number which is not used as a pdg code this serves as a place holder when constructing a vector from integers given by the user

Definition at line 257 of file FairEvtFilterOnSingleParticleCounts.h.

Referenced by AndMinMaxPdgCodes().

const Double_t FairEvtFilter::kNoChargeSpecified = -999.9
staticprotectedinherited

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