genie::evtlib::EventLibraryInterface Class Reference

Reads pre-generated events produced by an external generator. On an event-by-event basis, it can accept GENIE input specifying the neutrino and target IDs and neutrino energy and, therefore, it can re-use the upstream GENIE flux and geometry tools. More...

#include <EventLibraryInterface.h>

Inheritance diagram for genie::evtlib::EventLibraryInterface:

Public Member Functions
	EventLibraryInterface ()
	EventLibraryInterface (string config)
	~EventLibraryInterface ()
void	ProcessEventRecord (GHepRecord *event) const override
void	Configure (const Registry &config) override
void	Configure (string config) override
Public Member Functions inherited from genie::EventRecordVisitorI
virtual	~EventRecordVisitorI ()
Public Member Functions inherited from genie::Algorithm
virtual	~Algorithm ()
virtual void	FindConfig (void)
virtual const Registry &	GetConfig (void) const
Registry *	GetOwnedConfig (void)
virtual const AlgId &	Id (void) const
	Get algorithm ID. More...
virtual AlgStatus_t	GetStatus (void) const
	Get algorithm status. More...
virtual bool	AllowReconfig (void) const
virtual AlgCmp_t	Compare (const Algorithm *alg) const
	Compare with input algorithm. More...
virtual void	SetId (const AlgId &id)
	Set algorithm ID. More...
virtual void	SetId (string name, string config)
const Algorithm *	SubAlg (const RgKey &registry_key) const
void	AdoptConfig (void)
void	AdoptSubstructure (void)
virtual void	Print (ostream &stream) const
	Print algorithm info. More...

Protected Member Functions
const EvtLibRecord *	GetRecord (const Interaction *interaction) const
void	LoadRecords ()
void	Cleanup ()
void	FillKinematics (const GHepRecord &, Kinematics &kine, int primary_lep_id) const
Protected Member Functions inherited from genie::EventRecordVisitorI
	EventRecordVisitorI ()
	EventRecordVisitorI (string name)
	EventRecordVisitorI (string name, string config)
Protected Member Functions inherited from genie::Algorithm
	Algorithm ()
	Algorithm (string name)
	Algorithm (string name, string config)
void	Initialize (void)
void	DeleteConfig (void)
void	DeleteSubstructure (void)
Registry *	ExtractLocalConfig (const Registry &in) const
Registry *	ExtractLowerConfig (const Registry &in, const string &alg_key) const
	Split an incoming configuration Registry into a block valid for the sub-algo identified by alg_key. More...
template<class T >
bool	GetParam (const RgKey &name, T &p, bool is_top_call=true) const
template<class T >
bool	GetParamDef (const RgKey &name, T &p, const T &def) const
template<class T >
int	GetParamVect (const std::string &comm_name, std::vector< T > &v, bool is_top_call=true) const
	Handle to load vectors of parameters. More...
int	GetParamVectKeys (const std::string &comm_name, std::vector< RgKey > &k, bool is_top_call=true) const
int	AddTopRegistry (Registry *rp, bool owns=true)
	add registry with top priority, also update ownership More...
int	AddLowRegistry (Registry *rp, bool owns=true)
	add registry with lowest priority, also update ownership More...
int	MergeTopRegistry (const Registry &r)
int	AddTopRegisties (const vector< Registry * > &rs, bool owns=false)
	Add registries with top priority, also udated Ownerships. More...

Private Attributes
std::map< Key, const IEvtLibRecordList * >	fRecords
TFile *	fRecordFile

Additional Inherited Members
Static Public Member Functions inherited from genie::Algorithm
static string	BuildParamVectKey (const std::string &comm_name, unsigned int i)
static string	BuildParamVectSizeKey (const std::string &comm_name)
Protected Attributes inherited from genie::Algorithm
bool	fAllowReconfig
bool	fOwnsSubstruc
	true if it owns its substructure (sub-algs,...) More...
AlgId	fID
	algorithm name and configuration set More...
vector< Registry * >	fConfVect
vector< bool >	fOwnerships
	ownership for every registry in fConfVect More...
AlgStatus_t	fStatus
	algorithm execution status More...
AlgMap *	fOwnedSubAlgMp
	local pool for owned sub-algs (taken out of the factory pool) More...

Detailed Description

Reads pre-generated events produced by an external generator. On an event-by-event basis, it can accept GENIE input specifying the neutrino and target IDs and neutrino energy and, therefore, it can re-use the upstream GENIE flux and geometry tools.

Author

February 28, 2020

Copyright (c) 2003-2020, The GENIE Collaboration For the full text of the license visit http://copyright.genie-mc.org

Definition at line 38 of file EventLibraryInterface.h.

Constructor & Destructor Documentation

EventLibraryInterface::EventLibraryInterface

(

)

Definition at line 29 of file EventLibraryInterface.cxx.

                                             :
  EventRecordVisitorI("genie::evtlib::EventLibraryInterface"),
  fRecordFile(0)
{

}

EventLibraryInterface::EventLibraryInterface

(

string

config

)

Definition at line 37 of file EventLibraryInterface.cxx.

                                                          :
  EventRecordVisitorI("genie::evtlib::EventLibraryInterface", config),
  fRecordFile(0)
{

}

EventLibraryInterface::~EventLibraryInterface

(

)

Definition at line 45 of file EventLibraryInterface.cxx.

{
  Cleanup();
}

Member Function Documentation

void EventLibraryInterface::Cleanup ( )

protected

Definition at line 219 of file EventLibraryInterface.cxx.

{
  for(auto it: fRecords) delete it.second;
  fRecords.clear();
  delete fRecordFile;
  fRecordFile = 0;
}

void EventLibraryInterface::Configure ( const Registry &  config )

overridevirtual

Configure the algorithm with an external registry The registry is merged with the top level registry if it is owned, Otherwise a copy of it is added with the highest priority

Reimplemented from genie::Algorithm.

Definition at line 205 of file EventLibraryInterface.cxx.

{
  Algorithm::Configure(config);
  LoadRecords();
}

void EventLibraryInterface::Configure ( string  config )

overridevirtual

Configure the algorithm from the AlgoConfigPool based on param_set string given in input An algorithm contains a vector of registries coming from different xml configuration files, which are loaded according a very precise prioriy This methods will load a number registries in order of priority: 1) "Tunable" parameter set from CommonParametes. This is loaded with the highest prioriry and it is designed to be used for tuning procedure Usage not expected from the user. 2) For every string defined in "CommonParame" the corresponding parameter set will be loaded from CommonParameter.xml 3) parameter set specified by the config string and defined in the xml file of the algorithm 4) if config is not "Default" also the Default parameter set from the same xml file will be loaded Effectively this avoids the repetion of a parameter when it is not changed in the requested configuration

Reimplemented from genie::Algorithm.

Definition at line 212 of file EventLibraryInterface.cxx.

{
  Algorithm::Configure(config);
  LoadRecords();
}

void EventLibraryInterface::FillKinematics ( const GHepRecord &  event, Kinematics &  kine, int  primary_lep_id  ) const

protected

Definition at line 296 of file EventLibraryInterface.cxx.

                                                                       {
  

  const TLorentzVector & p_probe = * event.Particle( 0 ) -> P4() ; 

  const TLorentzVector & p_lep = * event.Particle( primary_lep_id ) -> P4() ; 
  
  const TLorentzVector q = p_probe - p_lep;

  // Initial hadronic state, semi-arbitrary
  const TLorentzVector & p_tgt = * event.Particle( 1 ) -> P4() ;
  
  kine.SetQ2(-q.Mag2(), true);
  kine.SetW((q + p_tgt).Mag(), true);
  kine.Setx(-q.Mag2() / (2*p_tgt.Dot(q)), true);
  kine.Sety( p_tgt.Dot( q ) / p_tgt.Dot( p_probe ), true ) ;

}

const EvtLibRecord * EventLibraryInterface::GetRecord ( const Interaction *  interaction ) const

protected

Definition at line 157 of file EventLibraryInterface.cxx.

{
  const InitialState& init_state = interaction->InitState();

  const double probe_E = init_state.ProbeE(kRfLab);

  if(!init_state.Tgt().IsNucleus()){
    LOG("ELI", pINFO) << "Skippping non-nuclear target " << init_state;
    return 0;
  }
  
  const int tgt_pdgc = init_state.TgtPdg();

  const ProcessInfo& proc = interaction->ProcInfo();

  if(!proc.IsWeakCC() && !proc.IsWeakNC()){
    LOG("ELI", pINFO) << "Skipping unknown process " << proc;
    return 0;
  }

  int probe_pdgc = init_state.ProbePdg();

  // Use nu_mu for NC as a convention internal to this code to index into the
  // records map.
  if(proc.IsWeakNC()){
    if(probe_pdgc > 0) probe_pdgc = kPdgNuMu; 
    else probe_pdgc = kPdgAntiNuMu;
  }

  const Key key(tgt_pdgc, probe_pdgc, proc.IsWeakCC());

  const auto rec_it = fRecords.find(key);

  if(rec_it == fRecords.end()){
    LOG("ELI", pINFO) << "Skippping " << key << " -- not found in library";
    return 0;
  }

  const EvtLibRecord* rec = rec_it->second->GetRecord(probe_E);

  if(!rec){
    LOG("ELI", pINFO) << "Skippping " << key << " at " << probe_E << " GeV -- not found in library";
    return 0;
  }

  return rec;
}

void EventLibraryInterface::LoadRecords ( )

protected

Definition at line 228 of file EventLibraryInterface.cxx.

{
  Cleanup();

  std::string libPath;
  GetParam("EventLibraryPath", libPath);
  Expand(libPath);

  bool onDemand;
  GetParam("OnDemand", onDemand);

  PDGLibrary* pdglib = PDGLibrary::Instance();

  fRecordFile = new TFile(libPath.c_str());
  if(fRecordFile->IsZombie()) exit(1);

  TIter next(fRecordFile->GetListOfKeys());
  while(TObject* dir = next()){
    const std::string& tgtName = dir->GetName();
    const TParticlePDG* tgtPart = pdglib->DBase()->GetParticle(tgtName.c_str());
    if(!tgtPart){
      LOG("ELI", pWARN) << "Unknown nucleus " << tgtName
                        << " found in " << libPath
                        << " -- skipping";
      continue;
    }

    for(int pdg: {kPdgNuE,   kPdgAntiNuE,
                  kPdgNuMu,  kPdgAntiNuMu,
                  kPdgNuTau, kPdgAntiNuTau}){

      for(bool iscc: {true, false}){
        // NCs should be the same for all flavours. Use nu_mu as a convention
        // internal to this code to index into the records map.
        if(!iscc && abs(pdg) != kPdgNuMu) continue;

        std::string nuName = pdglib->Find(pdg)->GetName();
        if(!iscc) nuName = pdg::IsAntiNeutrino(pdg) ? "nu_bar" : "nu";

        const std::string treeName =
          TString::Format("%s/%s/%s/records",
                          tgtName.c_str(),
                          iscc ? "cc" : "nc",
                          nuName.c_str()).Data();

        const Key key(tgtPart->PdgCode(), pdg, iscc);

        TTree* tr = (TTree*)fRecordFile->Get(treeName.c_str());

        if(!tr){
          LOG("ELI", pINFO) << treeName << " not found in "
                            << libPath << " -- skipping";
          continue;
        }

        if(onDemand)
          fRecords[key] = new OnDemandRecordList(tr, treeName);
        else
          fRecords[key] = new SimpleRecordList(tr, treeName);
      } // end for iscc
    } // end for pdg
  } // end for dir

  // Need to keep the record file open for OnDemand, but not Simple
  if(!onDemand){delete fRecordFile; fRecordFile = 0;}
}

void EventLibraryInterface::ProcessEventRecord ( GHepRecord *  event ) const

overridevirtual

Implements genie::EventRecordVisitorI.

Definition at line 51 of file EventLibraryInterface.cxx.

{
// Get event summary constructed by GENIE
//
  Interaction* interaction = event->Summary();
  const InitialState & init_state = interaction->InitState();

  const EvtLibRecord* rec = GetRecord(interaction);
  if(!rec) return; // Reason has already been printed

  
  // at this point we have the event to add to our event record
  // the energy might be not the same
  // so first we correct the probe:
  // we maintain the same direction but we set the energy selected 
  // by the record.
  // We only change the the particle in the event record, 
  // Not the summary so that we keep the history
  // of the event construction
  
  // Neutrino is a parent to the lepton(s)
  GHepParticle * probe = event->Particle(0) ; 
  
  TLorentzVector probe_p4 ( * probe->P4() ) ;

  LOG("ELI", pINFO) << "Difference between requested neutrino E and used neutrino energy: " << probe_p4.E() - rec -> E ;

  probe_p4 *= rec -> E / probe_p4.E() ; 
  probe -> SetMomentum( probe_p4 ) ;

  // because of the selection procedure there might not be enough 
  // events in the library
  // so we add a roation around the z axis which is random for each event
  double alpha = RandomGen::Instance()->RndEvg().Uniform(0, 2*constants::kPi );

  
  // now we simply have to add the particles from the EvtLibRecord
  // But we need to rotate the particles along the direction of the beam
  // so we evaluate this direction once and for all
  TVector3 unit_nudir = probe->P4()->Vect().Unit();


  int firstLep =  -1 ; 
  const auto & parts = rec -> parts ;

  for ( unsigned int i = 0; i < parts.size() ; ++i ) {
    
    const auto & part = parts[i] ;
    
    int pdg = part.pdg;
    
    if(pdg::IsLepton(part.pdg)) {
      if ( firstLep == -1 ) { 
        firstLep = i ;
        interaction -> KinePtr() -> SetFSLeptonP4( part.px, part.py, part.pz, part.E ) ;
      }
    }
    
    // Fix up the outgoing lepton for NC events (due to lepton universality
    // it could be any value in the library)
    
    if( interaction->ProcInfo().IsWeakNC() ) {
      if ( int(i) == firstLep ) {
        // in the case of NC the process is always stored as numu (or numubar)
        // as the couplings are the same
        // here we restore the proper pdg on what we think it is the leading neutrino
        pdg = init_state.ProbePdg();
      }
    }

    TLorentzVector p4(part.px, part.py, part.pz, part.E) ;

    // The next is the rotation of the event according to the common alpha
    p4.RotateZ( alpha ) ;

    // Now we rotate the events so that it matches the direction of the incoming neutrino
    p4.RotateUz(unit_nudir); 

    event->AddParticle(pdg,
                       kIStStableFinalState,
                       0 , 1, -1, -1, // child of the neutrino and nucleus
                       p4,
                       TLorentzVector(0, 0, 0, 0));
    
  }
  
  // we now have the particle in the right direction inside the record we just need to 
  // decorate the summary kinematics
  
  if ( firstLep < 0 ) {
    firstLep = 2 ;
  }
  else { 
    firstLep += 2 ; 
  }
  // these are events from outside, there are no assumptions on the lepton
  // If we are simulating BSM physics there might not even be a lepton in the final state
  // So if there is a lepton we evaluate Q2 based on that lepton
  // If not we evaluate Q2 using the first particle we have available 
  
  FillKinematics( *event, *interaction->KinePtr(), firstLep );
  
}

Member Data Documentation

TFile* genie::evtlib::EventLibraryInterface::fRecordFile

private

Definition at line 67 of file EventLibraryInterface.h.

std::map<Key, const IEvtLibRecordList*> genie::evtlib::EventLibraryInterface::fRecords

private

Definition at line 66 of file EventLibraryInterface.h.

---

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

• Generator/src/Tools/EvtLib/EventLibraryInterface.h
• Generator/src/Tools/EvtLib/EventLibraryInterface.cxx