evgb::GENIEHelper Class Reference

#include "GENIEHelper.h"

Inheritance diagram for evgb::GENIEHelper:

Public Member Functions
	GENIEHelper (fhicl::ParameterSet const &pset, TGeoManager *rootGeom, std::string const &rootFile, double const &detectorMass)
	~GENIEHelper ()
void	Initialize ()
bool	Stop ()
bool	Sample (simb::MCTruth &truth, simb::MCFlux &flux, simb::GTruth &gtruth)
double	TotalHistFlux ()
double	TotalExposure () const
double	SpillExposure () const
std::string	FluxType () const
std::string	DetectorLocation () const
std::vector< TH1D * >	FluxHistograms () const
double	TotalMass () const
genie::EventRecord *	GetGenieEventRecord ()
TRandom3 *	GetHelperRandom ()
genie::GFluxI *	GetFluxDriver (bool base=true)
std::string	GetTuneName () const
std::string	GetEventGeneratorList () const

Private Member Functions
void	RegularizeFluxType ()
void	SqueezeFilePatterns ()
void	AtmoFluxCheck ()
void	HistogramFluxCheck ()
void	InitializeGeometry ()
void	InitializeFiducialSelection ()
void	InitializeRockBoxSelection ()
void	InitializeFluxDriver ()
void	ConfigGeomScan ()
void	SetMaxPathOutInfo ()
void	BuildFluxRotation ()
void	ExpandFluxPaths ()
void	ExpandFluxFilePatternsDirect ()
void	ExpandFluxFilePatternsIFDH ()
bool	StringToBool (std::string v)
void	SetGXMLPATH ()
void	SetGMSGLAYOUT ()
void	StartGENIEMessenger (std::string prodmode)
void	ReadXSecTable ()

Private Attributes
TGeoManager *	fGeoManager
	pointer to ROOT TGeoManager More...
std::string	fGeoFile
	name of file containing the Geometry description More...
genie::EventRecord *	fGenieEventRecord
	last generated event More...
genie::GeomAnalyzerI *	fGeomD
genie::GFluxI *	fFluxD
	real flux driver More...
genie::GFluxI *	fFluxD2GMCJD
	flux driver passed to genie GMCJDriver, might be GFluxBlender More...
genie::GMCJDriver *	fDriver
ifdh_ns::ifdh *	fIFDH
	(optional) flux file handling More...
TRandom3 *	fHelperRandom
	random # generator for GENIEHelper More...
bool	fUseHelperRndGen4GENIE
	use fHelperRandom for gRandom during Sample() More...
evgb::EvtTimeShiftI *	fTimeShifter
	generator for time offset within a spill More...
std::string	fFluxType
std::string	fFluxSearchPaths
	colon separated set of path stems More...
std::vector< std::string >	fFluxFilePatterns
	wildcard patterns files containing histograms or ntuples, or txt More...
std::vector< std::string >	fSelectedFluxFiles
	flux files selected after wildcard expansion and subset selection More...
int	fMaxFluxFileMB
	maximum size of flux files (MB) More...
int	fMaxFluxFileNumber
	maximum # of flux files More...
std::string	fFluxCopyMethod
	"DIRECT" = old direct access method, otherwise = ifdh approach schema ("" okay) More...
std::string	fFluxCleanup
	"ALWAYS", "/var/tmp", "NEVER" More...
std::string	fBeamName
	name of the beam we are simulating More...
std::string	fFluxRotCfg
	how to interpret fFluxRotValues More...
std::vector< double >	fFluxRotValues
	parameters for rotation More...
TRotation *	fFluxRotation
	rotation for atmos / astro flux coord systems More...
std::string	fTopVolume
	top volume in the ROOT geometry in which to generate events More...
std::string	fWorldVolume
	name of the world volume in the ROOT geometry More...
std::string	fDetLocation
	name of flux window location More...
std::vector< TH1D * >	fFluxHistograms
	histograms for each nu species More...
double	fFluxUpstreamZ
	z where flux starts from (if non-default, simple/ntuple only) More...
double	fEventsPerSpill
double	fPOTPerSpill
	number of pot per spill More...
double	fHistEventsPerSpill
	number of events per spill for histogram fluxes - changes each spill More...
int	fSpillEvents
	total events for this spill More...
double	fSpillExposure
	total exposure (i.e. pot) for this spill More...
double	fTotalExposure
	pot used from flux ntuple More...
double	fMonoEnergy
	energy of monoenergetic neutrinos More...
std::string	fFunctionalFlux
int	fFunctionalBinning
double	fEmin
double	fEmax
double	fXSecMassPOT
	product of cross section, mass and POT/spill for histogram fluxes More...
double	fTotalHistFlux
	total flux of neutrinos from flux histograms for used flavors More...
TVector3	fBeamDirection
	direction of the beam for histogram fluxes More...
TVector3	fBeamCenter
	center of beam for histogram fluxes - must be in meters More...
double	fBeamRadius
	radius of cylindar for histogram fluxes - must be in meters More...
double	fDetectorMass
	mass of the detector in kg More...
double	fSurroundingMass
double	fGlobalTimeOffset
	overall time shift (ns) added to every particle time More...
double	fRandomTimeOffset
	additional random time shift (ns) added to every particle time More...
std::string	fSpillTimeConfig
	alternative to flat spill distribution More...
bool	fAddGenieVtxTime
	incorporate time from flux window to interaction point and (possibily) proton-on-target to flux window More...
bool	fForceApplyFlxWgt
	apply GFluxI::Weight() before returning event More...
std::vector< int >	fGenFlavors
	pdg codes for flavors to generate More...
double	fAtmoEmin
	atmo: Minimum energy of neutrinos in GeV More...
double	fAtmoEmax
	atmo: Maximum energy of neutrinos in GeV More...
double	fAtmoRl
	atmo: radius of the sphere on where the neutrinos are generated More...
double	fAtmoRt
double	fAtmoSpectralIndex
	atmo: Spectral index for power spectrum generation More...
std::vector< std::string >	fEnvironment
	environmental variables and settings used by genie More...
std::string	fXSecTable
	cross section file (was $GSPLOAD) More...
std::string	fTuneName
	GENIE R-3 Tune name (defines model configuration) More...
std::string	fEventGeneratorList
	control over event topologies, was $GEVGL [Default] More...
std::string	fGXMLPATH
	locations for GENIE XML files More...
std::string	fGMSGLAYOUT
	format for GENIE log message [BASIC]|SIMPLE (SIMPLE=no timestamps) More...
std::string	fGENIEMsgThresholds
	additional XML file setting Messager level thresholds (":" separated) More...
int	fGHepPrintLevel
	GHepRecord::SetPrintLevel(), -1=no-print. More...
std::string	fMixerConfig
	configuration string for genie GFlavorMixerI More...
double	fMixerBaseline
	baseline distance if genie flux can't calculate it More...
bool	fUseBlenderDist
	get neutrino's travel distance from blender (default: true) More...
std::string	fFiducialCut
	configuration for geometry selector More...
std::string	fGeomScan
	configuration for geometry scan to determine max pathlengths More...
std::string	fMaxPathOutInfo
	output info if writing PathLengthList from GeomScan More...
unsigned int	fDebugFlags
	set bits to enable debug info More...

Detailed Description

Definition at line 54 of file GENIEHelper.h.

Constructor & Destructor Documentation

evgb::GENIEHelper::GENIEHelper ( fhicl::ParameterSet const &  pset, TGeoManager *  rootGeom, std::string const &  rootFile, double const &  detectorMass  )

explicit

Determine which flux files to use Do this after random number seed initialization for stability

For atmos_ / astro_ fluxes we might need to set a coordinate system rotation

Set the GENIE environment if using entries in the fEnvironment vector

Definition at line 214 of file GENIEHelper.cxx.

References AtmoFluxCheck(), BuildFluxRotation(), ExpandFluxFilePatternsDirect(), ExpandFluxFilePatternsIFDH(), ExpandFluxPaths(), fAddGenieVtxTime, fAtmoEmax, fAtmoEmin, fAtmoRl, fAtmoRt, fAtmoSpectralIndex, fBeamCenter, fBeamDirection, fBeamRadius, fDebugFlags, fDetectorMass, fDetLocation, fEmax, fEmin, fEnvironment, fEventGeneratorList, fEventsPerSpill, fFiducialCut, fFluxCopyMethod, fFluxRotation, fFluxType, fFluxUpstreamZ, fForceApplyFlxWgt, fFunctionalBinning, fFunctionalFlux, fGenFlavors, fGENIEMsgThresholds, fGeomScan, fGHepPrintLevel, fGlobalTimeOffset, fGMSGLAYOUT, fGXMLPATH, fHelperRandom, fHistEventsPerSpill, fMixerBaseline, fMixerConfig, fMonoEnergy, fPOTPerSpill, fRandomTimeOffset, fSelectedFluxFiles, fSpillEvents, fSpillExposure, fSpillTimeConfig, fSurroundingMass, fTimeShifter, fTopVolume, fTotalExposure, fTuneName, fUseBlenderDist, fWorldVolume, fXSecTable, evgb::EvtTimeShiftFactory::GetEvtTimeShift(), evgb::EvtTimeShiftI::GetRandomGenerator(), evgb::GetRandomNumberSeed(), HistogramFluxCheck(), evgb::EvtTimeShiftFactory::Instance(), evgb::EvtTimeShiftI::IsRandomGeneratorSeeded(), art::left(), evgb::EvtTimeShiftFactory::Print(), evgb::EvtTimeShiftI::PrintConfig(), RegularizeFluxType(), SetGMSGLAYOUT(), SetGXMLPATH(), SqueezeFilePatterns(), and StartGENIEMessenger().

    : fGeoManager        (geoManager)
    , fGeoFile           (rootFile)
    , fGenieEventRecord  (0)
    , fGeomD             (0)
    , fFluxD             (0)
    , fFluxD2GMCJD       (0)
    , fDriver            (0)
    , fIFDH              (0)
    , fHelperRandom      (0)
    , fUseHelperRndGen4GENIE(pset.get< bool                  >("UseHelperRndGen4GENIE",true))
    , fTimeShifter       (0)
    , fFluxType          (pset.get< std::string              >("FluxType")               )
    , fFluxSearchPaths   (pset.get< std::string              >("FluxSearchPaths","")     )
    , fFluxFilePatterns  (pset.get< std::vector<std::string> >("FluxFiles")              )
    , fMaxFluxFileMB     (pset.get< int                      >("MaxFluxFileMB",    2000) ) // 2GB max default
    , fMaxFluxFileNumber (pset.get< int                      >("MaxFluxFileNumber",9999) ) // at most 9999 files
    , fFluxCopyMethod    (pset.get< std::string              >("FluxCopyMethod","DIRECT")) // "DIRECT" = old direct access method
    , fFluxCleanup       (pset.get< std::string              >("FluxCleanup","/var/tmp") ) // "ALWAYS", "NEVER", "/var/tmp"
    , fBeamName          (pset.get< std::string              >("BeamName")               )
    , fFluxRotCfg        (pset.get< std::string              >("FluxRotCfg","none")      )
    , fFluxRotValues     (pset.get< std::vector<double>      >("FluxRotValues", {} )     ) // default empty vector
    , fFluxRotation      (0)
    , fTopVolume         (pset.get< std::string              >("TopVolume")              )
    , fWorldVolume       ("volWorld")
    , fDetLocation       (pset.get< std::string              >("DetectorLocation")       )
    , fFluxUpstreamZ     (pset.get< double                   >("FluxUpstreamZ",  -2.e30) )
    , fEventsPerSpill    (pset.get< double                   >("EventsPerSpill",      0) )
    , fPOTPerSpill       (pset.get< double                   >("POTPerSpill",     0.0) )
    , fHistEventsPerSpill(0.)
    , fSpillEvents       (0)
    , fSpillExposure     (0.)
    , fTotalExposure     (0.)
    , fMonoEnergy        (pset.get< double                   >("MonoEnergy",        2.0) )
    , fFunctionalFlux    (pset.get< std::string              >("FunctionalFlux", "x") )
    , fFunctionalBinning (pset.get< int                      >("FunctionalBinning", 10000) )
    , fEmin              (pset.get< double                   >("FluxEmin", 0) )
    , fEmax              (pset.get< double                   >("FluxEmax", 10) )
    , fBeamRadius        (pset.get< double                   >("BeamRadius",        3.0) )
    , fDetectorMass      (detectorMass)
    , fSurroundingMass   (pset.get< double                   >("SurroundingMass",    0.) )
    , fGlobalTimeOffset  (pset.get< double                   >("GlobalTimeOffset", 1.e4) )
    , fRandomTimeOffset  (pset.get< double                   >("RandomTimeOffset", 1.e4) )
    , fSpillTimeConfig   (pset.get< std::string              >("SpillTimeConfig",    "") )
    , fAddGenieVtxTime   (pset.get< bool                     >("AddGenieVtxTime", false) )
    , fForceApplyFlxWgt  (pset.get< bool                     >("ForceApplyFlxWgt", true) )
    , fGenFlavors        (pset.get< std::vector<int>         >("GenFlavors")             )
    , fAtmoEmin          (pset.get< double                   >("AtmoEmin",          0.1) )
    , fAtmoEmax          (pset.get< double                   >("AtmoEmax",         10.0) )
    , fAtmoRl            (pset.get< double                   >("Rl",               20.0) )
    , fAtmoRt            (pset.get< double                   >("Rt",               20.0) )
    , fAtmoSpectralIndex (pset.get< double                   >("SpectralIndex",     2.0) )
    , fEnvironment       (pset.get< std::vector<std::string> >("Environment")            )
    , fXSecTable         (pset.get< std::string              >("XSecTable",          "") ) //e.g. "gxspl-FNALsmall.xml"
    , fTuneName          (pset.get< std::string              >("TuneName","${GENIE_XSEC_TUNE}") )
    , fEventGeneratorList(pset.get< std::string              >("EventGeneratorList", "Default") )
    , fGXMLPATH          (pset.get< std::string              >("GXMLPATH",           "") )
    , fGMSGLAYOUT        (pset.get< std::string              >("GMSGLAYOUT",         "") ) // [BASIC] or SIMPLE
    , fGENIEMsgThresholds(pset.get< std::string              >("GENIEMsgThresholds", "") ) // : separate list of files
    , fGHepPrintLevel    (pset.get< int                      >("GHepPrintLevel",     -1) ) // see GHepRecord::SetPrintLevel() -1=no-print
    , fMixerConfig       (pset.get< std::string              >("MixerConfig",    "none") )
    , fMixerBaseline     (pset.get< double                   >("MixerBaseline",      0.) )
    , fUseBlenderDist    (pset.get< bool                     >("UseBlenderDist",     true) )
    , fFiducialCut       (pset.get< std::string              >("FiducialCut",    "none") )
    , fGeomScan          (pset.get< std::string              >("GeomScan",    "default") )
    , fDebugFlags        (pset.get< unsigned int             >("DebugFlags",          0) )
  {

    // fEnvironment is (generally) deprecated ... print out any settings
    // but for some make it fatal because it's no longer supported and
    // can lead to unexpected behaviour
    if ( fEnvironment.size() > 0 ) {
      bool fatal = false;
      std::ostringstream fenvout, fenvfatal;
      fenvout << " fEnviroment.size() = " << fEnvironment.size();
      for (size_t i = 0; i < fEnvironment.size(); i += 2) {
        fenvout << std::endl << " [" << std::setw(20) << std::left
                << fEnvironment[i] << "] ==> "
                << fEnvironment[i+1];
        // fatal:     GEVGL, GSPLOAD
        // harmless:  GXMLPATH, GMSGLAYOUT, GMSGCONF, GPRODMODE
        if ( fEnvironment[i] == "GEVGL"   ||
             fEnvironment[i] == "GSPLOAD"    ) {
          fatal = true;
          std::string altparam = "";
          if ( fEnvironment[i] == "GEVGL"   ) altparam = "EventGeneratorList";
          if ( fEnvironment[i] == "GSPLOAD" ) altparam = "XSecTable";

          fenvfatal << std::endl << "use of \"" << fEnvironment[i] << "\""
                    << " is no longer supported as a Environment fcl parameter pair key,"
                    << "\n"
                    << "   please remove it from from the list and use"
                    << " the appropriate direct fcl parameter \""
                    << altparam << "\"" << std::endl;
        }
      }
      mf::LogInfo("GENIEHelper")
        << " Use of fEnvironment parameters is deprecated:\n"
        << fenvout.str();
      if ( fatal ) {
        mf::LogError("GENIEHelper") << fenvfatal.str();
        throw cet::exception("GENIEHelper")
          << " bad use of Environment fcl parameter";
      }
    }

    // for histogram, mono-energetic, functional form fluxes ...
    std::vector<double> beamCenter   (pset.get< std::vector<double> >("BeamCenter")   );
    std::vector<double> beamDirection(pset.get< std::vector<double> >("BeamDirection"));
    fBeamCenter.SetXYZ(beamCenter[0], beamCenter[1], beamCenter[2]);
    fBeamDirection.SetXYZ(beamDirection[0], beamDirection[1], beamDirection[2]);

    // special processing of GSEED (GENIE's random seed)... priority:
    //    if set in .fcl file RandomSeed variable, use that
    //    else if already set in environment use that
    //    else use evgb::GetRandomNumberSeed()
    int dfltseed;
    const char* gseedstr = std::getenv("GSEED");
    if ( gseedstr ) {
      dfltseed = strtol(gseedstr,NULL,0);
    } else {
      dfltseed = evgb::GetRandomNumberSeed();
    }
    int seedval = pset.get< int >("RandomSeed", dfltseed);
    // initialize random # generator for use within GENIEHelper
    mf::LogInfo("GENIEHelper") << "Init HelperRandom with seed " << seedval;
    fHelperRandom = new TRandom3(seedval);

    // clean up user input
    // also classifies flux type to simplify tests
    //  e.g.  atmo_  tree_
    RegularizeFluxType();


    // for tree-based fluxes
    // (e.g. "tree_numi" (nee "ntuple"), "tree_simple" and "tree_dk2nu")
    // we don't care about order and don't want duplicates
    // for others order might matter
    if ( fFluxType.find("tree_") == 0 ) SqueezeFilePatterns();

    ExpandFluxPaths();
    if ( fFluxCopyMethod == "DIRECT" ) ExpandFluxFilePatternsDirect();
    else                               ExpandFluxFilePatternsIFDH();

    if ( fFluxType.find("atmo_")  == 0 ||
         fFluxType.find("astro_") == 0    ) {
      BuildFluxRotation();
    }

    //    they should come in pairs of variable name key, then value

    // Process GXMLPATH extensions first, so they are available
    // when GENIE starts to get initialized; these might be
    // alternative locations for configurations (including
    // the GENIE Messenger system).
    SetGXMLPATH();

    // Also set GENIE log4cpp Messenger layout format before
    // initializing GENIE (can't be changed after singleton is created)
    SetGMSGLAYOUT();

    // Now initialize GENIE Messenger service
    StartGENIEMessenger(pset.get<std::string>("ProductionMode","false"));

    // In case we're printing the event record, how verbose should it be
    genie::GHepRecord::SetPrintLevel(fGHepPrintLevel);

    // Set GENIE's random # seed
    mf::LogInfo("GENIEHelper")
      << "Init genie::utils::app_init::RandGen() with seed " << seedval;
    genie::utils::app_init::RandGen(seedval);

    // special things for atmos fluxes
    if ( fFluxType.find("atmo_") == 0 ) AtmoFluxCheck();

    // make the histogram associations
    if ( fFluxType.find("histogram") == 0 ) HistogramFluxCheck();

    std::string flvlist;
    for ( std::vector<int>::iterator itr = fGenFlavors.begin(); itr != fGenFlavors.end(); itr++ )
      flvlist += Form(" %d",*itr);

    if ( fFluxType.find("mono") == 0 ) {
      fEventsPerSpill = 1;
      mf::LogInfo("GENIEHelper")
        << "Generating monoenergetic (" << fMonoEnergy
        << " GeV) neutrinos with the following flavors: "
        << flvlist;
    } else if ( fFluxType.find("function") == 0    ) {
      fEventsPerSpill = 1;
      mf::LogInfo("GENIEHelper")
        << "Generating neutrinos using the functional form: "
        << fFunctionalFlux << " E [" << fEmin << ":" << fEmax
        << "] GeV with " << fFunctionalBinning << " bins "
        << "with the following flavors: " << flvlist;
    } else if (fFluxType.find("PowerSpectrum") != std::string::npos) {
      mf::LogInfo("GENIEHelper")
        << "Generating neutrinos using the a power spectrum with Spectral index = "
        << fAtmoSpectralIndex
        << ", with the following flavors: " << flvlist
        << "\nThe energy range is between:  " << fAtmoEmin << " GeV and "
        << fAtmoEmax << " GeV."
        << '\n'
        << "  Generation surface of: (" << fAtmoRl << ","
        << fAtmoRt << ")";
    }  else {

      // flux methods other than "mono" and "function" require files
      std::string fileliststr;
      if ( fSelectedFluxFiles.empty() ) {
        fileliststr = "NO FLUX FILES FOUND!";
        mf::LogWarning("GENIEHelper")  << fileliststr;
      }
      else {
        std::vector<std::string>::iterator sitr = fSelectedFluxFiles.begin();
        for ( ; sitr != fSelectedFluxFiles.end(); sitr++) {
          fileliststr +=  "\n\t";
          fileliststr += *sitr;
        }
      }
      mf::LogInfo("GENIEHelper")
        << "Generating flux with the following flavors: " << flvlist
        << "\nand these file patterns: " << fileliststr;

    }

    // disallow conflicting settings here
    if ( ( fEventsPerSpill != 0 && fPOTPerSpill != 0 ) ||
         ( fEventsPerSpill == 0 && fPOTPerSpill == 0 )    ) {
      throw cet::exception("GENIEHelper")
        << "one or the other of EventsPerSpill ("
        << fEventsPerSpill << ") or "
        << "POTPerSpill ("
        << fPOTPerSpill << ") needs to be zero (but not both)";
    }

    if ( fEventsPerSpill != 0 ) {
      mf::LogInfo("GENIEHelper") << "Generating " << fEventsPerSpill
                                 << " events for each spill";

    } else {
      mf::LogInfo("GENIEHelper") << "Using " << fPOTPerSpill
                                 << " pot for each spill";
    }

    // how to distribute events in time
    if ( fSpillTimeConfig != "" ) {
      evgb::EvtTimeShiftFactory& timeShiftFactory = evgb::EvtTimeShiftFactory::Instance();
      fTimeShifter = timeShiftFactory.GetEvtTimeShift(fSpillTimeConfig);
      if ( fTimeShifter ) {
        if ( ! fTimeShifter->IsRandomGeneratorSeeded() ) {
          mf::LogInfo("GENIEHelper")
            << "Initialize spill time random generator with seed " << seedval;
          fTimeShifter->GetRandomGenerator()->SetSeed(seedval);
        }
        fTimeShifter->PrintConfig();
      } else {
        timeShiftFactory.Print();
      }
    }

    return;
  }

evgb::GENIEHelper::~GENIEHelper

(

)

Definition at line 488 of file GENIEHelper.cxx.

References e, fDriver, ff, fFluxCleanup, fFluxD, fGenieEventRecord, fGeomD, fHelperRandom, fIFDH, fMaxPathOutInfo, fSelectedFluxFiles, and fTotalExposure.

  {
    // user request writing out the scan of the geometry
    if ( fGeomD && fMaxPathOutInfo != "" ) {
      genie::geometry::ROOTGeomAnalyzer* rgeom =
        dynamic_cast<genie::geometry::ROOTGeomAnalyzer*>(fGeomD);

      string filename = "maxpathlength.xml";
      mf::LogInfo("GENIEHelper")
        << "Saving MaxPathLengths as: \"" << filename << "\"";

      const genie::PathLengthList& maxpath = rgeom->GetMaxPathLengths();

      maxpath.SaveAsXml(filename);
      // append extra info to file
      std::ofstream mpfile(filename.c_str(), std::ios_base::app);
      mpfile
        << std::endl
        << "<!-- this file is only relevant for a setup compatible with:"
        << std::endl
        << fMaxPathOutInfo
        << std::endl
        << "-->"
        << std::endl;
      mpfile.close();
    }  // finished writing max path length XML file (if requested)

    // protect against lack of driver due to not getting to Initialize()
    // (called from module's beginRun() method)
    if ( ! fDriver || ! fFluxD ) {
      mf::LogInfo("GENIEHelper")
        << "~GENIEHelper called, but previously failed to construct "
        << ( (fDriver) ? "":" genie::GMCJDriver" )
        << ( (fFluxD)  ? "":" genie::GFluxI" );
    } else {

      double probscale = fDriver->GlobProbScale();
      double rawpots   = 0;

      // rather than ask individual flux drivers for info
      // use the unified interface

      genie::flux::GFluxExposureI* fexposure =
        dynamic_cast<genie::flux::GFluxExposureI*>(fFluxD);
      if ( fexposure ) {
        rawpots = fexposure->GetTotalExposure();
      }
      genie::flux::GFluxFileConfigI* ffileconfig =
        dynamic_cast<genie::flux::GFluxFileConfigI*>(fFluxD);
      if ( ffileconfig ) {
        ffileconfig->PrintConfig();
      }

      mf::LogInfo("GENIEHelper")
        << " Total Exposure " << fTotalExposure
        << " GMCJDriver GlobProbScale " << probscale
        << " FluxDriver base pots " << rawpots
        << " corrected POTS " << rawpots/TMath::Max(probscale,1.0e-100);
    }

    // clean up owned genie object (other genie obj are ref ptrs)
    delete fGenieEventRecord;
    delete fDriver;
    delete fHelperRandom;

#ifndef NO_IFDH_LIB
  #ifdef USE_IFDH_SERVICE
    art::ServiceHandle<IFDH> ifdhp;
    if ( fFluxCleanup.find("ALWAYS")   == 0 ) {
      ifdhp->cleanup();
    } else if ( fFluxCleanup.find("/var/tmp") == 0 ) {
      auto ffitr = fSelectedFluxFiles.begin();
      for ( ; ffitr != fSelectedFluxFiles.end(); ++ffitr ) {
        std::string ff = *ffitr;
        if ( ff.find("/var/tmp") == 0 ) {
          mf::LogDebug("GENIEHelper") << "delete " << ff;
          ifdhp->rm(ff);
        }
      }
    }
  #else
    if ( fIFDH ) {
      if        ( fFluxCleanup.find("ALWAYS")   == 0 ) {
        fIFDH->cleanup();
      } else if ( fFluxCleanup.find("/var/tmp") == 0 ) {
        auto ffitr = fSelectedFluxFiles.begin();
        for ( ; ffitr != fSelectedFluxFiles.end(); ++ffitr ) {
          std::string ff = *ffitr;
          if ( ff.find("/var/tmp") == 0 ) {
            mf::LogDebug("GENIEHelper") << "delete " << ff;
            fIFDH->rm(ff);
          }
        }
      }
      delete fIFDH;
      fIFDH = 0;
    }
  #endif
#endif

  }

Member Function Documentation

void evgb::GENIEHelper::AtmoFluxCheck ( )

private

Speical pre-checks for atmo_ fluxes

Definition at line 779 of file GENIEHelper.cxx.

References fAtmoEmax, fAtmoEmin, fAtmoRl, fAtmoRt, fEventsPerSpill, fFluxType, fGenFlavors, and fSelectedFluxFiles.

Referenced by GENIEHelper().

  {

    if ( fGenFlavors.size() != fSelectedFluxFiles.size() ) {
      mf::LogInfo("GENIEHelper")
        <<  "ERROR: The number of generated neutrino flavors ("
        << fGenFlavors.size()
        << ") doesn't correspond to the number of files ("
        << fSelectedFluxFiles.size() << ")!!!";
      throw cet::exception("GENIEHelper")
        << "ERROR: atmo_ flavors != files";
    } else {
      std::ostringstream atmofluxmatch;
      for (size_t indx=0; indx < fGenFlavors.size(); ++indx ) {
        atmofluxmatch << "   " << std::setw(3) << fGenFlavors[indx]
                      << " " << fSelectedFluxFiles[indx] << "\n";
      }
      mf::LogInfo("GENIEHelper")
        <<  "atmo flux assignment : \n" << atmofluxmatch.str();
    }

    if ( fEventsPerSpill != 1 ) {
      mf::LogInfo("GENIEHelper")
        <<  "ERROR: For Atmospheric Neutrino generation,"
        << " EventPerSpill need to be 1!!";
      throw cet::exception("GENIEHelper")
        << "ERROR: " << fFluxType << " EventsPerSpill wasn't 1 ("
        << fEventsPerSpill << ")";
    }

    std::ostringstream atmofluxinfo;

    if (fFluxType.find("FLUKA") != std::string::npos ){
      atmofluxinfo << "  The fluxes are from FLUKA";
    }
    else if (fFluxType.find("BARTOL") != std::string::npos ||
             fFluxType.find("BGLRS")  != std::string::npos    ){
      atmofluxinfo << "  The fluxes are from BARTOL/BGLRS";
    }
    else if (fFluxType.find("HONDA") != std::string::npos ||
             fFluxType.find("HAKKM") != std::string::npos    ){
      atmofluxinfo << "  The fluxes are from HONDA/HAKKM";
    }
    else if (fFluxType.find("PowerSpectrum") != std::string::npos){
      atmofluxinfo << "  The fluxes are interpolated HONDA";
    }
    else {
      mf::LogInfo("GENIEHelper")
        << "Unknown atmo_ flux simulation: " << fFluxType;
      throw cet::exception("GENIEHelper")
        << "ERROR: bad atmo_ flux type " << fFluxType;
    }

    atmofluxinfo
      << '\n'
      << "  The energy range is between:  " << fAtmoEmin << " GeV and "
      << fAtmoEmax << " GeV."
      << '\n'
      << "  Generation surface of: (" << fAtmoRl << ","
      << fAtmoRt << ")";

    mf::LogInfo("GENIEHelper") << atmofluxinfo.str();

  }

void evgb::GENIEHelper::BuildFluxRotation ( )

private

Definition at line 1818 of file GENIEHelper.cxx.

References fFluxRotation, fFluxRotCfg, fFluxRotValues, and scale.

Referenced by GENIEHelper().

  {
    // construct fFluxRotation matrix
    //     from fFluxRotCfg + fFluxRotValues
    if ( fFluxRotCfg == "" ||
         ( fFluxRotCfg.find("none") != std::string::npos ) ) return;

    size_t nval = fFluxRotValues.size();

    bool verbose = ( fFluxRotCfg.find("verbose") != std::string::npos );
    if (verbose) {
      std::ostringstream indata;
      indata << "BuildFluxRotation: Cfg \"" << fFluxRotCfg << "\"\n"
             << " " << nval << " values\n";
      for (size_t i=0; i<nval; ++i) {
        indata << "   [" << std::setw(2) << i << "] " << fFluxRotValues[i] << "\n";
      }
      mf::LogInfo("GENIEHelper") << indata.str();
    }

    // interpret as a full 3x3 array
    if ( fFluxRotCfg.find("newxyz") != std::string::npos ||
         fFluxRotCfg.find("3x3")    != std::string::npos    ) {
      if ( nval == 9 ) {
        TRotation fTempRot;
        TVector3 newX = TVector3(fFluxRotValues[0],
                                 fFluxRotValues[1],
                                 fFluxRotValues[2]);
        TVector3 newY = TVector3(fFluxRotValues[3],
                                 fFluxRotValues[4],
                                 fFluxRotValues[5]);
        TVector3 newZ = TVector3(fFluxRotValues[6],
                                 fFluxRotValues[7],
                                 fFluxRotValues[8]);
        fTempRot.RotateAxes(newX,newY,newZ);
        // weirdly necessary; frame vs. obj rotation
        fFluxRotation = new TRotation(fTempRot.Inverse());
        return;
      } else {
        throw cet::exception("BadFluxRotation")
          << "specified: " << fFluxRotCfg << "\n"
          << " but nval=" << nval << ", need 9";
      }
    }

    // another possibility  ... series of rotations around particular axes
    if ( fFluxRotCfg.find("series") != std::string::npos ) {
      TRotation fTempRot;
      // if series then cfg should also have series of rot{X|Y|Z}{deg|rad}
      std::vector<std::string> strs = genie::utils::str::Split(fFluxRotCfg," ,;(){}[]");
      size_t nrot = -1;
      for (size_t j=0; j<strs.size(); ++j) {
        std::string what = strs[j];
        if ( what == ""        ) continue;
        // lower case for convenience
        std::transform(what.begin(),what.end(),what.begin(),::tolower);
        if ( what == "series"  ) continue;
        if ( what == "verbose" ) continue;
        if ( what.find("rot") != 0 ) {
          mf::LogWarning("GENIEHelper")
            << "processing series rotation saw keyword \"" << what << "\" -- ignoring";
          continue;
        }
        char axis = what[3];
        // check that axis is sensibly x, y or z
        if ( axis != 'x' && axis != 'y' && axis != 'z' ) {
          throw cet::exception("BadFluxRotation")
            << "specified: " << fFluxRotCfg << "\n"
            << " keyword '" << what << "': bad axis '" << axis << "'";
        }
        std::string units = what.substr(4);
        // don't worry if written fully as "radians" or "degrees"
        if (units.size() > 3) units.erase(3);
        if ( units != "" && units != "rad" && units != "deg" ) {
          throw cet::exception("BadFluxRotation")
            << "specified: " << fFluxRotCfg << "\n"
            << " keyword '" << what << "': bad units '" << units << "'";
        }
        // no units?  assume degrees
        double scale = (( units == "rad" ) ? 1.0 : TMath::DegToRad() );

        ++nrot;
        if ( nrot >= nval ) {
          // not enough values
          throw cet::exception("BadFluxRotation")
            << "specified: " << fFluxRotCfg << "\n"
            << " asking for rotation [" << nrot << "] "
            << what << " but nval=" << nval;
        }
        double rot = scale * fFluxRotValues[nrot];
        if ( axis == 'x' ) fTempRot.RotateX(rot);
        if ( axis == 'y' ) fTempRot.RotateY(rot);
        if ( axis == 'z' ) fTempRot.RotateZ(rot);

      } // loop over tokens in cfg string

      // weirdly necessary; frame vs. obj rotation
      fFluxRotation = new TRotation(fTempRot.Inverse());

      if ( nrot+1 != nval ) {
        // call out possible user mistake
        mf::LogWarning("GENIEHelper")
          << "BuildFluxRotation only used " << nrot+1 << " of "
          << nval << " FluxRotValues";
      }
      return;
    }

    // could put other interpretations here ...

    throw cet::exception("BadFluxRotation")
      << "specified: " << fFluxRotCfg << "\n"
      << " nval=" << nval << ", but don't know how to interpret that";

  }

void evgb::GENIEHelper::ConfigGeomScan ( )

private

Definition at line 1474 of file GENIEHelper.cxx.

References util::abs(), fDriver, fFluxD, fGeomD, fGeomScan, and SetMaxPathOutInfo().

Referenced by Initialize().

  {

    // trim any leading whitespace
    if( fGeomScan.find_first_not_of(" \t\n") != 0)
      fGeomScan.erase( 0, fGeomScan.find_first_not_of(" \t\n")  );

    if ( fGeomScan.find("default") == 0 ) return;

    genie::geometry::ROOTGeomAnalyzer* rgeom =
      dynamic_cast<genie::geometry::ROOTGeomAnalyzer*>(fGeomD);

    if ( !rgeom ) {
      throw cet::exception("GENIEHelper") << "fGeomD wasn't of type "
                                          << "genie::geometry::ROOTGeomAnalyzer*";
    }

    // convert string to lowercase

    // parse out string
    vector<string> strtok = genie::utils::str::Split(fGeomScan," ");
    // first value is a string, others should be numbers unless "file:"
    string scanmethod = strtok[0];

    // convert key string to lowercase (but not potential file name)
    std::transform(scanmethod.begin(),scanmethod.end(),scanmethod.begin(),::tolower);

    if ( scanmethod.find("file") == 0 ) {
      // xml expand path before passing in
      string filename = strtok[1];
      string fullname = genie::utils::xml::GetXMLFilePath(filename);
      mf::LogInfo("GENIEHelper")
        << "ConfigGeomScan getting MaxPathLengths from \"" << fullname << "\"";
      fDriver->UseMaxPathLengths(fullname);
      return;
    }

    vector<double> vals;
    for ( size_t indx=1; indx < strtok.size(); ++indx ) {
      const string& valstr1 = strtok[indx];
      if ( valstr1 != "" ) vals.push_back(atof(valstr1.c_str()));
    }
    size_t nvals = vals.size();
    // pad it out to at least 4 entries to avoid index issues
    for ( size_t nadd = 0; nadd < 4-nvals; ++nadd ) vals.push_back(0);

    double safetyfactor = 0;
    int    writeout = 0;
    if (        scanmethod.find("box") == 0 ) {
      // use box method
      int np = (int)vals[0];
      int nr = (int)vals[1];
      if ( nvals >= 3 ) safetyfactor = vals[2];
      if ( nvals >= 4 ) writeout     = vals[3];
      // protect against too small values
      if ( np <= 10 ) np = rgeom->ScannerNPoints();
      if ( nr <= 10 ) nr = rgeom->ScannerNRays();
      mf::LogInfo("GENIEHelper")
        << "ConfigGeomScan scan using box " << np << " points, "
        << nr << " rays";
      rgeom->SetScannerNPoints(np);
      rgeom->SetScannerNRays(nr);
    } else if ( scanmethod.find("flux") == 0 ) {
      // use flux method
      int np = (int)vals[0];
      if ( nvals >= 2 ) safetyfactor = vals[1];
      if ( nvals >= 3 ) writeout     = vals[2];
      // sanity check, need some number of rays to explore the geometry
      // negative now means adjust rays to enu_max (GENIE svn 3676)
      if ( abs(np) <= 100 ) {
        int npnew = rgeom->ScannerNParticles();
        if ( np < 0 ) npnew = -abs(npnew);
        mf::LogWarning("GENIEHelper")
          << "Too few rays requested for geometry scan: " << np
          << ", use: " << npnew << "instead";
        np = npnew;
      }
      mf::LogInfo("GENIEHelper")
        << "ConfigGeomScan scan using " << np << " flux particles"
        << ( (np>0) ? "" : " with ray energy pushed to flux driver maximum" );
      rgeom->SetScannerFlux(fFluxD);
      rgeom->SetScannerNParticles(np);
    }
    else{
      // unknown
      throw cet::exception("GENIEHelper") << "fGeomScan unknown method: \""
                                          << fGeomScan << "\"";
    }
    if ( safetyfactor > 0 ) {
      mf::LogInfo("GENIEHelper")
        << "ConfigGeomScan setting safety factor to " << safetyfactor;
      rgeom->SetMaxPlSafetyFactor(safetyfactor);
    }
    if ( writeout != 0 ) SetMaxPathOutInfo();
  }

std::string evgb::GENIEHelper::DetectorLocation ( ) const

inline

Definition at line 77 of file GENIEHelper.h.

{ return fDetLocation;    }

void evgb::GENIEHelper::ExpandFluxFilePatternsDirect ( )

private

Definition at line 1953 of file GENIEHelper.cxx.

References fFluxFilePatterns, fFluxSearchPaths, fFluxType, fHelperRandom, fMaxFluxFileMB, fMaxFluxFileNumber, fSelectedFluxFiles, and geo::vect::indices().

Referenced by GENIEHelper().

  {
    // Using the the fFluxSearchPaths list of directories, apply the
    // user supplied pattern as a suffix to find the flux files.
    // The userpattern might include simple wildcard globs (in contrast
    // to proper regexp patterns).

    // After expanding the list to individual files, randomize them
    // and start selecting until a size limit is about to be
    // exceeded (though a minimum there needs to be one file, not
    // matter the limit).

    // older GENIE (pre r3669) can't handle vectors/sets of file names
    // but an only accept a pattern that will resolve to files.  This
    // collection could exceed the desired size threshold, but for
    // pick the collection that expands to the largest list

    bool randomizeFiles = false;
    if ( fFluxType.find("tree_") == 0 ) randomizeFiles = true;

    std::vector<std::string> dirs;
    cet::split_path(fFluxSearchPaths,dirs);
    if ( dirs.empty() ) dirs.push_back(std::string()); // at least null string

    glob_t g;
    int flags = GLOB_TILDE;   // expand ~ home directories

    std::ostringstream patterntext;  // for info/error messages
    std::ostringstream dirstext;     // for info/error messages

    std::vector<std::string>::const_iterator uitr = fFluxFilePatterns.begin();
    int ipatt = 0;

    for ( ; uitr != fFluxFilePatterns.end(); ++uitr, ++ipatt ) {
      std::string userpattern = *uitr;
      patterntext << "\n\t" << userpattern;

      std::vector<std::string>::const_iterator ditr = dirs.begin();
      for ( ; ditr != dirs.end(); ++ditr ) {
        std::string dalt = *ditr;
        // if non-null, does it end with a "/"?  if not add one
        size_t len = dalt.size();
        if ( len > 0 && dalt.rfind('/') != len-1 ) dalt.append("/");
        if ( uitr == fFluxFilePatterns.begin() ) dirstext << "\n\t" << dalt;

        std::string filepatt = dalt + userpattern;

        glob(filepatt.c_str(),flags,NULL,&g);
        if ( g.gl_pathc > 0 ) flags |= GLOB_APPEND; // next glob() will append to list

      }  // loop over FluxSearchPaths dirs
    }  // loop over user patterns

    std::ostringstream paretext;
    std::ostringstream flisttext;

    int nfiles = g.gl_pathc;

    if ( nfiles == 0 ) {
      paretext << "\n  expansion resulted in a null list for flux files";

    } else if ( ! randomizeFiles ) {
      // some sets of files should be left in order
      // and no size limitations imposed ... just copy the list

      paretext << "\n  list of files will be processed in order";
      for (int i=0; i<nfiles; ++i) {
        std::string afile(g.gl_pathv[i]);
        fSelectedFluxFiles.push_back(afile);

        flisttext << "[" << setw(3) << i << "] "
                  << afile << "\n";
      }
    } else {

      // now pull from the list randomly
      // do this by assigning a random number to each;
      // ordering that list; and pulling in that order

      paretext << "list of " << nfiles
               << " will be randomized and pared down to "
               << fMaxFluxFileMB << " MB or "
               << fMaxFluxFileNumber << " files";

      double* order = new double[nfiles];
      int* indices  = new int[nfiles];
      fHelperRandom->RndmArray(nfiles,order);
      // assign random # for their relative order

      TMath::Sort((int)nfiles,order,indices,false);

      long long int sumBytes = 0; // accumulated size in bytes
      long long int maxBytes = (long long int)fMaxFluxFileMB * 1024ll * 1024ll;

      FileStat_t fstat;
      for (int i=0; i < TMath::Min(nfiles,fMaxFluxFileNumber); ++i) {
        int indx = indices[i];
        std::string afile(g.gl_pathv[indx]);
        bool keep = true;

        gSystem->GetPathInfo(afile.c_str(),fstat);
        sumBytes += fstat.fSize;
        // skip those that would push sum above total
        // but always accept at least one (the first)
        if ( sumBytes > maxBytes && i != 0 ) keep = false;

        flisttext << "[" << setw(3) << i << "] "
                  << "=> g[" << setw(3) << indx << "] "
                  << ((keep)?"keep":"skip") << " "
                  << setw(6) << (sumBytes/(1024ll*1024ll)) << " "
                  << afile << "\n";

        if ( keep ) fSelectedFluxFiles.push_back(afile);
        else break;  // <voice name=Scotty> Captain, she can't take any more</voice>

      }
      delete [] order;
      delete [] indices;

    }

    mf::LogInfo("GENIEHelper")
      << "ExpandFluxFilePatternsDirect initially found " << nfiles
      << " files for user patterns:"
      << patterntext.str() << "\n  using FluxSearchPaths of: "
      << dirstext.str() <<  "\n" << paretext.str();
      //<< "\"" << cet::getenv("FW_SEARCH_PATH") << "\"";

    mf::LogDebug("GENIEHelper") << "\n" << flisttext.str();

    // done with glob list
    globfree(&g);

    // no null path allowed for at least these
    if ( fFluxType.find("tree_") == 0 ) {
      size_t nfiles = fSelectedFluxFiles.size();
      if ( nfiles == 0 ) {
        mf::LogError("GENIEHelper")
          << "For \"" << fFluxType <<"\" "
          << "(e.g. \"dk2nu\', \"ntuple\" (\"numi\") or \"simple\") "
          << " specification must resolve to at least one file"
          << "\n  none were found. DIRECT user pattern(s): "
          << patterntext.str()
          << "\n  using FluxSearchPaths of: "
          << dirstext.str();

        throw cet::exception("NoFluxFiles")
          << "no flux files found for: " << patterntext.str() << "\n"
          << " in: " << dirstext.str();

      }
    }

  } // ExpandFluxFilePatternsDirect

void evgb::GENIEHelper::ExpandFluxFilePatternsIFDH ( )

private

Definition at line 2109 of file GENIEHelper.cxx.

References fFluxCopyMethod, fFluxFilePatterns, fFluxSearchPaths, fFluxType, fHelperRandom, fIFDH, fMaxFluxFileMB, fMaxFluxFileNumber, fSelectedFluxFiles, and geo::vect::indices().

Referenced by GENIEHelper().

  {
    // Using the the FluxSearchPaths list of directories, apply the
    // user supplied pattern as a suffix to find the flux files.
    // The userpattern might include simple wildcard globs (in contrast
    // to proper regexp patterns).

    // After expanding the list to individual files, randomize them
    // and start selecting until a size limit is about to be
    // exceeded (though a minimum there needs to be one file, not
    // matter the limit).

    // Use the IFDH interface to get the list of files and sizes;
    // after sorting/selecting use IFDH to make a local copy

#ifdef NO_IFDH_LIB
    std::ostringstream fmesg;
    std::string marker = "%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%\n";
    fmesg << marker
          << __FILE__ << ":" << __LINE__
          << "\nno IFDH implemented on this platform\n"
          << marker;
    // make sure the message goes everywhere
    std::cout << fmesg.str() << std::flush;
    std::cerr << fmesg.str();
    throw cet::exception("Attempt to use ifdh class") << fmesg.str();
    assert(0);
#else
    // if "method" just an identifier and not a scheme then clear it
    if ( fFluxCopyMethod.find("IFDH") == 0 ) fFluxCopyMethod = "";

    bool randomizeFiles = false;
    if ( fFluxType.find("tree_") == 0 ) randomizeFiles = true;

  #ifdef USE_IFDH_SERVICE
    art::ServiceHandle<IFDH> ifdhp;
  #else
    if ( ! fIFDH ) fIFDH = new ifdh_ns::ifdh;
  #endif

    std::string spaths = fFluxSearchPaths;

    const char* ifdh_debug_env = std::getenv("IFDH_DEBUG_LEVEL");
    if ( ifdh_debug_env ) {
      mf::LogInfo("GENIEHelper") << "IFDH_DEBUG_LEVEL: " << ifdh_debug_env;
      fIFDH->set_debug(ifdh_debug_env);
    }

    // filenames + size
    std::vector<std::pair<std::string,long>>
      partiallist, fulllist, selectedlist, locallist;

    std::ostringstream patterntext;  // stringification of vector of patterns
    std::ostringstream fulltext;     // for info on full list of matches
    std::ostringstream selectedtext; // for info on selected files
    std::ostringstream localtext;    // for info on local files
    fulltext << "search paths: " << spaths;

    //std::vector<std::string>::const_iterator uitr = fFluxFilePatterns.begin();

    // loop over possible patterns
    // IFDH handles various stems but done a list of globs
    size_t ipatt=0;
    auto uitr = fFluxFilePatterns.begin();
    for ( ; uitr != fFluxFilePatterns.end(); ++uitr, ++ipatt ) {
      std::string userpattern = *uitr;
      patterntext << "\npattern [" << std::setw(3) << ipatt << "] " << userpattern;
      fulltext    << "\npattern [" << std::setw(3) << ipatt << "] " << userpattern;

  #ifdef USE_IFDH_SERVICE
      partiallist = ifdhp->findMatchingFiles(spaths,userpattern);
  #else
      partiallist = fIFDH->findMatchingFiles(spaths,userpattern);
  #endif
      fulllist.insert(fulllist.end(),partiallist.begin(),partiallist.end());

      // make a complete list ...
      fulltext << " found " << partiallist.size() << " files";
      for (auto pitr = partiallist.begin(); pitr != partiallist.end(); ++pitr) {
        fulltext << "\n  " << std::setw(10) << pitr->second << " " << pitr->first;
      }

      partiallist.clear();
    }  // loop over user patterns

    size_t nfiles = fulllist.size();

    mf::LogInfo("GENIEHelper")
      << "ExpandFluxFilePatternsIFDH initially found " << nfiles << " files";
    mf::LogDebug("GENIEHelper")
      << fulltext.str();

    if ( nfiles == 0 ) {
      selectedtext << "\n  expansion resulted in a null list for flux files";
    } else if ( ! randomizeFiles ) {
      // some sets of files should be left in order
      // and no size limitations imposed ... just copy the list

      selectedtext << "\n  list of files will be processed in order";
      selectedlist.insert(selectedlist.end(),fulllist.begin(),fulllist.end());

    } else {

      // for list needing size based trimming and randomization ...
      // pull from the list randomly until a cummulative limit is reached
      // do this by assigning a random number to each;
      // ordering that list; and pulling in that order

      selectedtext << "list of " << nfiles
                   << " will be randomized and pared down to "
                   << fMaxFluxFileMB << " MB or "
                   << fMaxFluxFileNumber << " files";

      double* order = new double[nfiles];
      int* indices  = new int[nfiles];
      fHelperRandom->RndmArray(nfiles,order);
      // assign random # for their relative order

      TMath::Sort((int)nfiles,order,indices,false);

      long long int sumBytes = 0; // accumulated size in bytes
      long long int maxBytes = (long long int)fMaxFluxFileMB * 1024ll * 1024ll;

      for (size_t i=0; i < TMath::Min(nfiles,size_t(fMaxFluxFileNumber)); ++i) {
        int indx = indices[i];
        bool keep = true;

        auto p = fulllist[indx]; // this pair <name,size>
        sumBytes += p.second;
        // skip those that would push sum above total
        // but always accept at least one (the first)
        if ( sumBytes > maxBytes && i != 0 ) keep = false;

        selectedtext << "\n[" << setw(3) << i << "] "
                     << "=> [" << setw(3) << indx << "] "
                     << ((keep)?"keep":"SKIP") << " "
                     << std::setw(6) << (sumBytes/(1024ll*1024ll)) << " MB "
                     << p.first;

        if ( keep ) selectedlist.push_back(p);
        else break;  // <voice name=Scotty> Captain, she can't take any more</voice>

      }
      delete [] order;
      delete [] indices;
    }

    mf::LogInfo("GENIEHelper")
      << selectedtext.str();

    // have a selected list of remote files
    // get paths to local copies
  #ifdef USE_IFDH_SERVICE
    locallist = ifdhp->fetchSharedFiles(selectedlist,fFluxCopyMethod);
  #else
    locallist = fIFDH->fetchSharedFiles(selectedlist,fFluxCopyMethod);
  #endif

    localtext << "final list of files:";
    size_t i=0;
    for (auto litr = locallist.begin(); litr != locallist.end(); ++litr, ++i) {
        fSelectedFluxFiles.push_back(litr->first);
        localtext << "\n\t[" << std::setw(3) << i << "]\t" << litr->first;
      }

    mf::LogInfo("GENIEHelper")
      << localtext.str();

    // no null path allowed for at least these
    if ( fFluxType.find("tree_") == 0 ) {
      size_t nfiles = fSelectedFluxFiles.size();
      if ( nfiles == 0 ) {
        mf::LogError("GENIEHelper")
          << "For \"" << fFluxType <<"\" "
          << "(e.g. \"dk2nu\', \"ntuple\" (\"numi\") or \"simple\") "
          << "specification must resolve to at least one file"
          << "\n  none were found. IFDH user pattern(s): "
          << patterntext.str()
          << "\n  using FluxSearchPaths of: "
          << spaths;

        throw cet::exception("NoFluxFiles")
          << "no flux files found for: " << patterntext.str() << "\n"
          << " in " << spaths;

      }
    }
#endif  // 'else' code only if NO_IFDH_LIB not defined
  } // ExpandFluxFilePatternsIFDH

void evgb::GENIEHelper::ExpandFluxPaths ( )

private

Definition at line 1934 of file GENIEHelper.cxx.

References fFluxCopyMethod, and fFluxSearchPaths.

Referenced by GENIEHelper().

  {
    // expand any wildcards in the paths variable
    // if unset and using the old DIRECT method allow it to fall back
    // to using FW_SEARCH_PATH ... but not for the new ifdhc approach

    std::string initial = fFluxSearchPaths;

    if ( fFluxCopyMethod == "DIRECT" && fFluxSearchPaths == "" ) {
      fFluxSearchPaths = cet::getenv("FW_SEARCH_PATH");
    }
    fFluxSearchPaths = gSystem->ExpandPathName(fFluxSearchPaths.c_str());

    mf::LogInfo("GENIEHelper")
      << "ExpandFluxPaths initially: \"" << initial << "\"\n"
      << "             final result: \"" << fFluxSearchPaths << "\"\n"
      << "                    using: \"" << fFluxCopyMethod << "\" method";
  }

std::vector<TH1D*> evgb::GENIEHelper::FluxHistograms ( ) const

inline

Definition at line 81 of file GENIEHelper.h.

{ return fFluxHistograms; }

std::string evgb::GENIEHelper::FluxType ( ) const

inline

Definition at line 76 of file GENIEHelper.h.

{ return fFluxType;       }

std::string evgb::GENIEHelper::GetEventGeneratorList ( ) const

inline

Definition at line 96 of file GENIEHelper.h.

{ return fEventGeneratorList; }

genie::GFluxI* evgb::GENIEHelper::GetFluxDriver ( bool  base = true )

inline

Definition at line 91 of file GENIEHelper.h.

Referenced by evgen::TestGENIEHelper::produce(), and evgen::GENIEGen::produce().

      { return ( (base) ? fFluxD : fFluxD2GMCJD ); }

genie::EventRecord* evgb::GENIEHelper::GetGenieEventRecord ( )

inline

Definition at line 84 of file GENIEHelper.h.

Referenced by evgen::TestGENIEHelper::produce().

{ return fGenieEventRecord; }

TRandom3* evgb::GENIEHelper::GetHelperRandom ( )

inline

Definition at line 87 of file GENIEHelper.h.

{ return fHelperRandom; }

std::string evgb::GENIEHelper::GetTuneName ( ) const

inline

Definition at line 95 of file GENIEHelper.h.

{ return fTuneName; }

void evgb::GENIEHelper::HistogramFluxCheck ( )

private

Definition at line 846 of file GENIEHelper.cxx.

References fBeamCenter, fBeamDirection, fBeamRadius, fFluxHistograms, fGenFlavors, fSelectedFluxFiles, and fTotalHistFlux.

Referenced by GENIEHelper().

  {

    mf::LogInfo("GENIEHelper")
      << "setting beam direction and center at "
      << fBeamDirection.X() << " " << fBeamDirection.Y()
      << " " << fBeamDirection.Z()
      << " (" << fBeamCenter.X() << "," << fBeamCenter.Y()
      << "," << fBeamCenter.Z()
      << ") with radius " << fBeamRadius;

    TDirectory *savedir = gDirectory;

    fFluxHistograms.clear();

    TFile tf((*fSelectedFluxFiles.begin()).c_str());
    tf.ls();

    for ( std::vector<int>::iterator flvitr = fGenFlavors.begin(); flvitr != fGenFlavors.end(); flvitr++){
      const char* histname = "none";
      switch ( *flvitr ) {
        case  12: histname = "nue";      break;
        case -12: histname = "nuebar";   break;
        case  14: histname = "numu";     break;
        case -14: histname = "numubar";  break;
        case  16: histname = "nutau";    break;
        case -16: histname = "nutaubar"; break;
        default: {
          throw cet::exception("GENIEHelper")
            <<  "ERROR: no support for histogram flux with flavor PDG="
            << *flvitr;
        }
      }
      fFluxHistograms.push_back(dynamic_cast<TH1D *>(tf.Get(histname)));
    }

    for ( unsigned int i = 0; i < fFluxHistograms.size(); ++i ) {
      fFluxHistograms[i]->SetDirectory(savedir);
      fTotalHistFlux += fFluxHistograms[i]->Integral();
    }

    mf::LogInfo("GENIEHelper")
      << "total histogram flux over desired flavors = "
      << fTotalHistFlux;

  }

void evgb::GENIEHelper::Initialize

(

)

Definition at line 606 of file GENIEHelper.cxx.

References ConfigGeomScan(), e, fDetectorMass, fDriver, fEventGeneratorList, fEventsPerSpill, fFluxD, fFluxD2GMCJD, fFluxType, fGeomD, fHelperRandom, fHistEventsPerSpill, fPOTPerSpill, fSpillEvents, fSpillExposure, fSurroundingMass, fTotalExposure, fTotalHistFlux, fTuneName, fXSecMassPOT, InitializeFluxDriver(), InitializeGeometry(), ReadXSecTable(), and evgb::SetEventGeneratorListAndTune().

Referenced by evgen::GENIEGen::beginJob(), and evgen::TestGENIEHelper::beginRun().

  {
#ifdef GENIE_PRE_R3
  // trigger early initialization of PDG database & GENIE message service
  // just to get it out of the way and not intermixed with other output
  genie::PDGLibrary::Instance();
#else
  // get the GENIE banner out of the way
  // no longer can use genie::PDGLibrary::Instance() to do this
  // because that must happen, in some cases in v3_02_xx, after the tune
  // is determined
  // banner is triggered by first use of GENIE Messenger
  // avoid using GENIE macros (possible conflict with mf macros)
  // LOG("GENIE",pInfo) << "Trigger GENIE banner";
  (*genie::Messenger::Instance())("GENIE") << log4cpp::Priority::INFO
                                           << "Trigger GENIE banner";
#endif


    mf::LogInfo("GENIEHelper") << "GENIE EventGeneratorList using \""
                               << fEventGeneratorList << "\"";

    evgb::SetEventGeneratorListAndTune(fEventGeneratorList,fTuneName);

#ifdef GENIE_PRE_R3
    // no tunes ... so can report current setting if different from fcl
#else
    genie::RunOpt* grunopt = genie::RunOpt::Instance();

    // RunOpt's SetEventGeneratorList() wasn't introduced until R-3
    // so above call could not have modified it in RunOpt
    std::string currentEvtGenListName = grunopt->EventGeneratorList();
    if ( currentEvtGenListName != fEventGeneratorList ) {
      mf::LogInfo("GENIEHelper") << " EventGeneratorList name changed from \""
                                 << fEventGeneratorList << "\" to \""
                                 << currentEvtGenListName << "\""
                                 << " by evgb::SetEventGeneratorListAndTune";
      fEventGeneratorList = currentEvtGenListName;
    }

    std::string currentTuneName = grunopt->Tune()->Name();
    if ( currentTuneName != fTuneName ) {
      mf::LogInfo("GENIEHelper") << " TuneName name changed from \""
                                 << fTuneName << "\" to \""
                                 << currentTuneName << "\""
                                 << " by evgb::SetEventGeneratorListAndTune";
      fTuneName = currentTuneName;
    }
#endif


    fDriver = new genie::GMCJDriver(); // this needs to be before ConfigGeomScan
    // let the driver know which list
    // (for R-3 this is in _addition_ to setting it in RunOpt)
    // pre-R-3 one could not set it in RunOpt directly (ie. without parsing
    //    command line args), thus don't try to fetch it from there.
    // post-R-3 one *-could-* use
    //   fDriver->SetEventGeneratorList(RunOpt::Instance()->EventGeneratorList());
    fDriver->SetEventGeneratorList(fEventGeneratorList);

    // Figure out which cross section file to use
    // post R-2_8_0 this actually triggers reading the file
    ReadXSecTable();

    // initialize the Geometry and Flux drivers
    InitializeGeometry();
    InitializeFluxDriver();

    fDriver->UseFluxDriver(fFluxD2GMCJD);
    fDriver->UseGeomAnalyzer(fGeomD);

    // must come after creation of Geom, Flux and GMCJDriver
    ConfigGeomScan();      // could trigger fDriver->UseMaxPathLengths(*xmlfile*)

    fDriver->Configure();  // could trigger GeomDriver::ComputeMaxPathLengths()
    fDriver->UseSplines();
    fDriver->ForceSingleProbScale();

    if ( fFluxType.find("histogram") == 0 && fEventsPerSpill < 0.01 ) {
      // fluxes are assumed to be given in units of neutrinos/cm^2/1e20POT/energy
      // integral over all fluxes removes energy dependence
      // histograms should have bin width that reflects the value of the /energy bit
      // ie if /energy = /50MeV then the bin width should be 50 MeV

      // determine product of pot/spill, mass, and cross section
      // events = flux * pot * 10^-38 cm^2 (xsec) * (mass detector (in kg) / nucleon mass (in kg))
      fXSecMassPOT  = 1.e-38*1.e-20;
      fXSecMassPOT *= fPOTPerSpill*(fDetectorMass+fSurroundingMass)/(1.67262158e-27);

      mf::LogInfo("GENIEHelper") << "Number of events per spill will be based on poisson mean of "
                                 << fXSecMassPOT*fTotalHistFlux;

      fHistEventsPerSpill = fHelperRandom->Poisson(fXSecMassPOT*fTotalHistFlux);
    }

    // set the pot/event counters to zero
    fSpillEvents   = 0;
    fSpillExposure = 0.;
    fTotalExposure = 0.;

    // If the flux driver knows how to keep track of exposure (time,pots)
    // reset it now as some might have been used in determining
    // the geometry maxpathlength or internally scanning for weights.
    // Should be happen for all cases since R-2_8_0 .. but no harm doing it ourselves

    fFluxD->Clear("CycleHistory");

    return;
  }

void evgb::GENIEHelper::InitializeFiducialSelection ( )

private

User defined fiducial volume cut [0][M]<SHAPE>:val1,val2,... "0" means reverse the cut (i.e. exclude the volume) "M" means the coordinates are given in the ROOT geometry "master" system and need to be transformed to "top vol" system <SHAPE> can be any of "zcyl" "box" "zpoly" "sphere" [each takes different # of args] This must be followed by a ":" and a list of values separated by punctuation (allowed separators: commas , parentheses () braces {} or brackets [] ) Value mapping: zcly:x0,y0,radius,zmin,zmax - cylinder along z at (x0,y0) capped at z's box:xmin,ymin,zmin,xmax,ymax,zmax - box w/ upper & lower extremes zpoly:nfaces,x0,y0,r_in,phi,zmin,zmax - nfaces sided polygon in x-y plane

Examples: 1) 0mbox:0,0,0.25,1,1,8.75 exclude (i.e. reverse) a box in master coordinates w/ corners (0,0,0.25) (1,1,8.75) 2) mzpoly:6,(2,-1),1.75,0,{0.25,8.75} six sided polygon in x-y plane, centered at x,y=(2,-1) w/ inscribed radius 1.75 no rotation (so first face is in y-z plane +r from center, i.e. hex sits on point) limited to the z range of {0.25,8.75} in the master ROOT geom coordinates 3) zcly:(3,4),5.5,-2,10 a cylinder oriented parallel to the z axis in the "top vol" coordinates at x,y=(3,4) with radius 5.5 and z range of {-2,10}

Definition at line 930 of file GENIEHelper.cxx.

References fFiducialCut, fGeomD, and InitializeRockBoxSelection().

Referenced by InitializeGeometry().

  {
    genie::GeomAnalyzerI* geom_driver = fGeomD; // GENIEHelper name -> gNuMIExptEvGen name
    std::string fidcut = fFiducialCut;   // ditto

    if( fidcut.find_first_not_of(" \t\n") != 0) // trim any leading whitespace
      fidcut.erase( 0, fidcut.find_first_not_of(" \t\n")  );

    // convert string to lowercase
    std::transform(fidcut.begin(),fidcut.end(),fidcut.begin(),::tolower);

    if ( "" == fidcut || "none" == fidcut ) return;

    if ( fidcut.find("rock") != string::npos ) {
      // deal with RockBox separately than basic shapes
      InitializeRockBoxSelection();
      return;
    }

    // below is as it is in $GENIE/src/support/numi/EvGen/gNuMIExptEvGen
    // except the change in message logger from log4cpp (GENIE) to cet's MessageLogger used by art

    //       sphere:x0,y0,z0,radius                - sphere of fixed radius at (x0,y0,z0)
    genie::geometry::ROOTGeomAnalyzer * rgeom =
      dynamic_cast<genie::geometry::ROOTGeomAnalyzer *>(geom_driver);
    if ( ! rgeom ) {
      mf::LogWarning("GENIEHelpler")
        << "Can not create GeomVolSelectorFiduction,"
        << " geometry driver is not ROOTGeomAnalyzer";
      return;
    }

    mf::LogInfo("GENIEHelper") << "fiducial cut: " << fidcut;

    // for now, only fiducial no "rock box"
    genie::geometry::GeomVolSelectorFiducial* fidsel =
      new genie::geometry::GeomVolSelectorFiducial();

    fidsel->SetRemoveEntries(true);  // drop segments that won't be considered

    vector<string> strtok = genie::utils::str::Split(fidcut,":");
    if ( strtok.size() != 2 ) {
      mf::LogWarning("GENIEHelper")
        << "Can not create GeomVolSelectorFiduction,"
        << " no \":\" separating type from values.  nsplit=" << strtok.size();
      for ( unsigned int i=0; i < strtok.size(); ++i )
        mf::LogWarning("GENIEHelper")
          << "strtok[" << i << "] = \"" << strtok[i] << "\"";
      return;
    }

    // parse out optional "x" and "m"
    string stype = strtok[0];
    bool reverse = ( stype.find("0") != string::npos );
    bool master  = ( stype.find("m") != string::npos );  // action after values are set

    // parse out values
    vector<double> vals;
    vector<string> valstrs = genie::utils::str::Split(strtok[1]," ,;(){}[]");
    vector<string>::const_iterator iter = valstrs.begin();
    for ( ; iter != valstrs.end(); ++iter ) {
      const string& valstr1 = *iter;
      if ( valstr1 != "" ) vals.push_back(atof(valstr1.c_str()));
    }
    size_t nvals = vals.size();
    // pad it out to at least 7 entries to avoid index issues if used
    for ( size_t nadd = 0; nadd < 7-nvals; ++nadd ) vals.push_back(0);

    //std::cout << "ivals = [";
    //for (unsigned int i=0; i < nvals; ++i) {
    //  if (i>0) cout << ",";
    //  std::cout << vals[i];
    //}
    //std::cout << "]" << std::endl;

    // std::vector elements are required to be adjacent so we can treat address as ptr

    if        ( stype.find("zcyl")   != string::npos ) {
      // cylinder along z direction at (x0,y0) radius zmin zmax
      if ( nvals < 5 )
        mf::LogError("GENIEHelper") << "MakeZCylinder needs 5 values, not " << nvals
                                    << " fidcut=\"" << fidcut << "\"";
      fidsel->MakeZCylinder(vals[0],vals[1],vals[2],vals[3],vals[4]);

    } else if ( stype.find("box")    != string::npos ) {
      // box (xmin,ymin,zmin) (xmax,ymax,zmax)
      if ( nvals < 6 )
        mf::LogError("GENIEHelper") << "MakeBox needs 6 values, not " << nvals
                                    << " fidcut=\"" << fidcut << "\"";
      double xyzmin[3] = { vals[0], vals[1], vals[2] };
      double xyzmax[3] = { vals[3], vals[4], vals[5] };
      fidsel->MakeBox(xyzmin,xyzmax);

    } else if ( stype.find("zpoly")  != string::npos ) {
      // polygon along z direction nfaces at (x0,y0) radius phi zmin zmax
      if ( nvals < 7 )
        mf::LogError("GENIEHelper") << "MakeZPolygon needs 7 values, not " << nvals
                                    << " fidcut=\"" << fidcut << "\"";
      int nfaces = (int)vals[0];
      if ( nfaces < 3 )
        mf::LogError("GENIEHelper") << "MakeZPolygon needs nfaces>=3, not " << nfaces
                                    << " fidcut=\"" << fidcut << "\"";
      fidsel->MakeZPolygon(nfaces,vals[1],vals[2],vals[3],vals[4],vals[5],vals[6]);

    } else if ( stype.find("sphere") != string::npos ) {
      // sphere at (x0,y0,z0) radius
      if ( nvals < 4 )
        mf::LogError("GENIEHelper") << "MakeZSphere needs 4 values, not " << nvals
                                    << " fidcut=\"" << fidcut << "\"";
      fidsel->MakeSphere(vals[0],vals[1],vals[2],vals[3]);

    } else {
      mf::LogError("GENIEHelper")
        << "Can not create GeomVolSelectorFiduction for shape \"" << stype << "\"";
    }

    if ( master  ) {
      fidsel->ConvertShapeMaster2Top(rgeom);
      mf::LogInfo("GENIEHelper") << "Convert fiducial volume from master to topvol coords";
    }
    if ( reverse ) {
      fidsel->SetReverseFiducial(true);
      mf::LogInfo("GENIEHelper") << "Reverse sense of fiducial volume cut";
    }

    rgeom->AdoptGeomVolSelector(fidsel);

  }

void evgb::GENIEHelper::InitializeFluxDriver ( )

private

Definition at line 1180 of file GENIEHelper.cxx.

References genie::flux::GPowerSpectrumAtmoFlux::AddFluxFile(), fAtmoEmax, fAtmoEmin, fAtmoRl, fAtmoRt, fAtmoSpectralIndex, fBeamCenter, fBeamDirection, fBeamRadius, fDebugFlags, fDetLocation, fEmax, fEmin, fFluxD, fFluxD2GMCJD, fFluxHistograms, fFluxRotation, fFluxType, fFluxUpstreamZ, fFunctionalBinning, fFunctionalFlux, fGenFlavors, fMixerBaseline, fMixerConfig, fMonoEnergy, fSelectedFluxFiles, genie::flux::GPowerSpectrumAtmoFlux::LoadFluxData(), genie::flux::GPowerSpectrumAtmoFlux::SetFlavors(), genie::flux::GPowerSpectrumAtmoFlux::SetMaxEnergy(), genie::flux::GPowerSpectrumAtmoFlux::SetMinEnergy(), genie::flux::GPowerSpectrumAtmoFlux::SetRadii(), genie::flux::GPowerSpectrumAtmoFlux::SetSpectralIndex(), genie::flux::GPowerSpectrumAtmoFlux::SetUserCoordSystem(), w, and weight.

Referenced by Initialize().

  {

    // simplify a lot of things ...
    // but for now this part only handles the 3 ntuple styles
    // that support the GFluxFileConfig mix-in
    // not the atmos, histo or mono versions

    std::string fluxName = "";

    // what looks like the start of a fully qualified class name
    // or one of our tree_ classes "numi" "simple" "dk2nu"
    // but only know how to configure those that dervie from:
    //    genie::flux::GFluxFileConfigI*
    if ( fFluxType.find("genie::flux::")   != std::string::npos )
      fluxName = fFluxType;
    else if ( fFluxType.find("tree_numi")    == 0 )
      fluxName = "genie::flux::GNuMIFlux";
    else if ( fFluxType.find("tree_simple")  == 0 )
      fluxName = "genie::flux::GSimpleNtpFlux";
    else if ( fFluxType.find("tree_dk2nu")   == 0 )
      fluxName = "genie::flux::GDk2NuFlux";

    if ( fluxName != "" ) {
      // any fall through to hopefully be handled below ...
      genie::flux::GFluxDriverFactory& fluxDFactory =
        genie::flux::GFluxDriverFactory::Instance();
      fFluxD = fluxDFactory.GetFluxDriver(fluxName);
      if ( ! fFluxD ) {
        mf::LogInfo("GENIEHelper")
          << "genie::flux::GFluxDriverFactory had not result for '"
          << fluxName << "' (fFluxType was '" << fFluxType << "'";
        // fall through in hopes someone later picks it up
      } else {
        // got something
        // for the ones that support GFluxFileConfigI (numi,simple,dk2nu)
        // initialize them
        genie::flux::GFluxFileConfigI* ffileconfig =
          dynamic_cast<genie::flux::GFluxFileConfigI*>(fFluxD);
        if ( ffileconfig ) {
          ffileconfig->LoadBeamSimData(fSelectedFluxFiles,fDetLocation);
          ffileconfig->PrintConfig();
          // initialize to only use neutrino flavors requested by user
          genie::PDGCodeList probes;
          for ( std::vector<int>::iterator flvitr = fGenFlavors.begin(); flvitr != fGenFlavors.end(); flvitr++ )
            probes.push_back(*flvitr);
          ffileconfig->SetFluxParticles(probes);
          if ( TMath::Abs(fFluxUpstreamZ) < 1.0e30 ) ffileconfig->SetUpstreamZ(fFluxUpstreamZ);
        }
      }
    } // is genie::flux:: or tree_{numi|simple|dk2nu}

    if ( fFluxType.find("histogram") == 0 ) {

      genie::flux::GCylindTH1Flux* histFlux = new genie::flux::GCylindTH1Flux();

      // now add the different fluxes - fluxes were added to the vector in the same
      // order that the flavors appear in fGenFlavors
      int ctr = 0;
      for ( std::vector<int>::iterator i = fGenFlavors.begin(); i != fGenFlavors.end(); i++ ) {
        histFlux->AddEnergySpectrum(*i, fFluxHistograms[ctr]);
        ++ctr;
      } //end loop to add flux histograms to driver

      histFlux->SetNuDirection(fBeamDirection);
      histFlux->SetBeamSpot(fBeamCenter);
      histFlux->SetTransverseRadius(fBeamRadius);

      fFluxD = histFlux; // dynamic_cast<genie::GFluxI *>(histFlux);

    } // end if using a histogram
    else if ( fFluxType.find("mono") == 0 ) {

      // weight each species equally in the generation
      double weight = 1./(1.*fGenFlavors.size());
      //make a map of pdg to weight codes
      std::map<int, double> pdgwmap;
      for ( std::vector<int>::iterator i = fGenFlavors.begin(); i != fGenFlavors.end(); i++ )
        pdgwmap[*i] = weight;

      genie::flux::GMonoEnergeticFlux *monoflux = new genie::flux::GMonoEnergeticFlux(fMonoEnergy, pdgwmap);
      monoflux->SetDirectionCos(fBeamDirection.X(), fBeamDirection.Y(), fBeamDirection.Z());
      monoflux->SetRayOrigin(fBeamCenter.X(), fBeamCenter.Y(), fBeamCenter.Z());
      fFluxD = monoflux; // dynamic_cast<genie::GFluxI *>(monoflux);

    } //end if using monoenergetic beam
    else if ( fFluxType.find("function") == 0 ) {

      genie::flux::GCylindTH1Flux* histFlux = new genie::flux::GCylindTH1Flux();
      TF1* input_func = new TF1("input_func", fFunctionalFlux.c_str(), fEmin, fEmax);
      TH1D* spectrum = new TH1D("spectrum", "neutrino flux", fFunctionalBinning, fEmin, fEmax);
      spectrum->Add(input_func);

      for ( std::vector<int>::iterator i = fGenFlavors.begin(); i != fGenFlavors.end(); i++ ) {
        histFlux->AddEnergySpectrum(*i, spectrum);
      } //end loop to add flux histograms to driver

      histFlux->SetNuDirection(fBeamDirection);
      histFlux->SetBeamSpot(fBeamCenter);
      histFlux->SetTransverseRadius(fBeamRadius);

      fFluxD = histFlux; // dynamic_cast<genie::GFluxI *>(histFlux);
      delete input_func;
    } //end if using function beam

    else if(fFluxType.find("PowerSpectrum") != std::string::npos){
      genie::flux::GPowerSpectrumAtmoFlux *power_flux = new genie::flux::GPowerSpectrumAtmoFlux;
      power_flux->SetSpectralIndex(fAtmoSpectralIndex);
      power_flux->SetFlavors(fGenFlavors);
      power_flux->SetMinEnergy(fAtmoEmin);
      power_flux->SetMaxEnergy(fAtmoEmax);

      mf::LogInfo("GENIEHelper") << "Setting Emin=" << fEmin << " ; Emax=" << fEmax;

      for ( size_t j = 0; j < fGenFlavors.size(); ++j ) {
        int         flavor  = fGenFlavors[j];
        std::string flxfile = fSelectedFluxFiles[j];
        power_flux->AddFluxFile(flavor,flxfile);
      }


      if ( fFluxRotation ){
        power_flux->SetUserCoordSystem(*fFluxRotation);
        std::ostringstream atmoCfgText;
        const int w=13, p=6;
        auto old_p = atmoCfgText.precision(p);
        atmoCfgText << "\n UserCoordSystem rotation:\n"
                    << "  [ "
                    << std::setw(w) << fFluxRotation->XX() << " "
                    << std::setw(w) << fFluxRotation->XY() << " "
                    << std::setw(w) << fFluxRotation->XZ() << " ]\n"
                    << "  [ "
                    << std::setw(w) << fFluxRotation->YX() << " "
                    << std::setw(w) << fFluxRotation->YY() << " "
                    << std::setw(w) << fFluxRotation->YZ() << " ]\n"
                    << "  [ "
                    << std::setw(w) << fFluxRotation->ZX() << " "
                    << std::setw(w) << fFluxRotation->ZY() << " "
                    << std::setw(w) << fFluxRotation->ZZ() << " ]\n";
        atmoCfgText.precision(old_p);
        mf::LogInfo("GENIEHelper") << atmoCfgText.str();
      }

      power_flux->LoadFluxData();

      // configure flux generation surface:
      power_flux->SetRadii(fAtmoRl, fAtmoRt);

      fFluxD = power_flux;
    }

    // Using the atmospheric fluxes
    else if ( fFluxType.find("atmo_") == 0 ) {

      // Instantiate appropriate concrete flux driver
      genie::flux::GAtmoFlux *atmo_flux_driver = 0;

      if ( fFluxType.find("FLUKA") != std::string::npos ) {
        genie::flux::GFLUKAAtmoFlux * fluka_flux =
          new genie::flux::GFLUKAAtmoFlux;
        atmo_flux_driver = dynamic_cast<genie::flux::GAtmoFlux *>(fluka_flux);
      }
      if ( fFluxType.find("BARTOL") != std::string::npos ||
           fFluxType.find("BGLRS")  != std::string::npos    ) {
        genie::flux::GBGLRSAtmoFlux * bartol_flux =
          new genie::flux::GBGLRSAtmoFlux;
        atmo_flux_driver = dynamic_cast<genie::flux::GAtmoFlux *>(bartol_flux);
      }
#if __GENIE_RELEASE_CODE__ >= GRELCODE(2,12,2)
      if (fFluxType.find("atmo_HONDA") != std::string::npos ||
          fFluxType.find("atmo_HAKKM") != std::string::npos    ) {
        genie::flux::GHAKKMAtmoFlux * honda_flux =
          new genie::flux::GHAKKMAtmoFlux;
        atmo_flux_driver = dynamic_cast<genie::flux::GAtmoFlux *>(honda_flux);
      }
#endif

      atmo_flux_driver->ForceMinEnergy(fAtmoEmin);
      atmo_flux_driver->ForceMaxEnergy(fAtmoEmax);
      if ( fFluxRotation ) atmo_flux_driver->SetUserCoordSystem(*fFluxRotation);

      std::ostringstream atmoCfgText;
      atmoCfgText << "Configuration for " << fFluxType
                  << ", Rl " << fAtmoRl << " Rt " << fAtmoRt;
      for ( size_t j = 0; j < fGenFlavors.size(); ++j ) {
        int         flavor  = fGenFlavors[j];
        std::string flxfile = fSelectedFluxFiles[j];
        atmo_flux_driver->AddFluxFile(flavor,flxfile); // pre-R-2_11_0 was SetFluxFile()
        atmoCfgText << "\n  FLAVOR: " << std::setw(3) << flavor
                    << "  FLUX FILE: " <<  flxfile;
      }
      if ( fFluxRotation ) {
        const int w=13, p=6;
        auto old_p = atmoCfgText.precision(p);
        atmoCfgText << "\n UserCoordSystem rotation:\n"
                    << "  [ "
                    << std::setw(w) << fFluxRotation->XX() << " "
                    << std::setw(w) << fFluxRotation->XY() << " "
                    << std::setw(w) << fFluxRotation->XZ() << " ]\n"
                    << "  [ "
                    << std::setw(w) << fFluxRotation->YX() << " "
                    << std::setw(w) << fFluxRotation->YY() << " "
                    << std::setw(w) << fFluxRotation->YZ() << " ]\n"
                    << "  [ "
                    << std::setw(w) << fFluxRotation->ZX() << " "
                    << std::setw(w) << fFluxRotation->ZY() << " "
                    << std::setw(w) << fFluxRotation->ZZ() << " ]\n";
          atmoCfgText.precision(old_p);
      }
      mf::LogInfo("GENIEHelper") << atmoCfgText.str();

      atmo_flux_driver->LoadFluxData();

      // configure flux generation surface:
      atmo_flux_driver->SetRadii(fAtmoRl, fAtmoRt);

      fFluxD = atmo_flux_driver;//dynamic_cast<genie::GFluxI *>(atmo_flux_driver);

    } //end if using atmospheric fluxes

    if ( ! fFluxD ) {
      mf::LogError("GENIEHelper")
        << "Failed to contruct base flux driver for FluxType '"
        << fFluxType << "'";
      throw cet::exception("GENIEHelper")
        << "Failed to contruct base flux driver for FluxType '"
        << fFluxType << "'\n"
        << __FILE__ << ":" << __LINE__ << "\n";
    }

    //
    // Is the user asking to do flavor mixing?
    //
    fFluxD2GMCJD = fFluxD;  // default: genie's GMCJDriver uses the bare flux generator
    if( fMixerConfig.find_first_not_of(" \t\n") != 0) // trim any leading whitespace
      fMixerConfig.erase( 0, fMixerConfig.find_first_not_of(" \t\n")  );
    std::string keyword = fMixerConfig.substr(0,fMixerConfig.find_first_of(" \t\n"));
    if ( keyword != "none" && keyword != "" ) {
      // Wrap the true flux driver up in the adapter to allow flavor mixing
      genie::flux::GFlavorMixerI* mixer = 0;
      // here is where we map MixerConfig string keyword to actual class
      // first is a special case that is part of GENIE proper
      if ( keyword == "map" || keyword == "swap" || keyword == "fixedfrac" )
        mixer = new genie::flux::GFlavorMap();
      // if it wasn't one of the predefined known mixers then
      // see if the factory knows about it and can create one
      // assuming the keyword (first token) is the class name
      if ( ! mixer ) {
        genie::flux::GFlavorMixerFactory& mixerFactory =
          genie::flux::GFlavorMixerFactory::Instance();
        mixer = mixerFactory.GetFlavorMixer(keyword);
        if ( mixer ) {
          // remove class name from config string
          fMixerConfig.erase(0,keyword.size());
          // trim any leading whitespace
          if ( fMixerConfig.find_first_not_of(" \t\n") != 0 )
            fMixerConfig.erase( 0, fMixerConfig.find_first_not_of(" \t\n")  );
        } else {
          const std::vector<std::string>& knownMixers =
            mixerFactory.AvailableFlavorMixers();
          std::ostringstream mixers;
          for (unsigned int j=0; j < knownMixers.size(); ++j ) {
            mixers << "\n   [" << std::setw(2) << j << "]  " << knownMixers[j];
          }
          mf::LogWarning("GENIEHelper")
            << " GFlavorMixerFactory known mixers: " << mixers.str();
        }
      }
      // configure the mixer
      if ( mixer ) mixer->Config(fMixerConfig);
      else {
        mf::LogWarning("GENIEHelper")
          << "GENIEHelper MixerConfig keyword was \"" << keyword
          << "\" but that did not map to a class; " << std::endl
          << "GFluxBlender in use, but no mixer";
      }

      genie::GFluxI* realFluxD = fFluxD;
      genie::flux::GFluxBlender* blender = new genie::flux::GFluxBlender();
      blender->SetBaselineDist(fMixerBaseline);
      blender->AdoptFluxGenerator(realFluxD);
      blender->AdoptFlavorMixer(mixer);
      fFluxD2GMCJD = blender;
      if ( fDebugFlags & 0x01 ) {
        if ( mixer ) mixer->PrintConfig();
        blender->PrintConfig();
        std::cout << std::flush;
      }
    }

    return;
  }

void evgb::GENIEHelper::InitializeGeometry ( )

private

Definition at line 894 of file GENIEHelper.cxx.

References fDebugFlags, fGeoManager, fGeomD, fTopVolume, fWorldVolume, and InitializeFiducialSelection().

Referenced by Initialize().

  {
    genie::geometry::ROOTGeomAnalyzer *rgeom =
      new genie::geometry::ROOTGeomAnalyzer(fGeoManager);

    // pass some of the debug flag bits on to the geometry manager
    int geomFlags = ( fDebugFlags >> 16 ) & 0xFF ;
    if ( geomFlags ) {
      int keep = ( geomFlags >> 7 );
      mf::LogInfo("GENIEHelper")
        << "InitializeGeometry set debug 0x"
        << std::hex << geomFlags << std::dec
        << " keepSegPath " << keep;
      rgeom->SetDebugFlags(geomFlags);
      if ( keep ) rgeom->SetKeepSegPath(true);
    }

    // get the world volume name from the geometry
    fWorldVolume = fGeoManager->GetTopVolume()->GetName();

    // the detector geometry uses cgs units.
    rgeom->SetLengthUnits(genie::units::centimeter);
    rgeom->SetDensityUnits(genie::units::gram_centimeter3);
    rgeom->SetTopVolName(fTopVolume.c_str());
    rgeom->SetMixtureWeightsSum(1.);
    mf::LogInfo("GENIEHelper")
      << "GENIEHelper::InitializeGeometry using TopVolume '"
      << fTopVolume << "'";
    //  casting to the GENIE geometry driver interface
    fGeomD        = rgeom; // dynamic_cast<genie::GeomAnalyzerI *>(rgeom);
    InitializeFiducialSelection();

    return;
  }

void evgb::GENIEHelper::InitializeRockBoxSelection ( )

private

Definition at line 1085 of file GENIEHelper.cxx.

References e, fFiducialCut, fGeomD, fTopVolume, and fWorldVolume.

Referenced by InitializeFiducialSelection().

  {
    genie::GeomAnalyzerI* geom_driver = fGeomD; // GENIEHelper name -> gNuMIExptEvGen name
    std::string fidcut = fFiducialCut;   // ditto

    if( fidcut.find_first_not_of(" \t\n") != 0) // trim any leading whitespace
      fidcut.erase( 0, fidcut.find_first_not_of(" \t\n")  );

    // convert string to lowercase
    std::transform(fidcut.begin(),fidcut.end(),fidcut.begin(),::tolower);

    genie::geometry::ROOTGeomAnalyzer * rgeom =
      dynamic_cast<genie::geometry::ROOTGeomAnalyzer *>(geom_driver);
    if ( ! rgeom ) {
      mf::LogWarning("GENIEHelpler")
        << "Can not create GeomVolSelectorRockBox,"
        << " geometry driver is not ROOTGeomAnalyzer";
      return;
    }

    mf::LogInfo("GENIEHelper") << "fiducial (rock) cut: " << fidcut;

    // for now, only fiducial no "rock box"
    genie::geometry::GeomVolSelectorRockBox* rocksel =
      new genie::geometry::GeomVolSelectorRockBox();

    vector<string> strtok = genie::utils::str::Split(fidcut,":");
    if ( strtok.size() != 2 ) {
      mf::LogWarning("GENIEHelper")
        << "Can not create GeomVolSelectorRockBox,"
        << " no \":\" separating type from values.  nsplit=" << strtok.size();
      for ( unsigned int i=0; i < strtok.size(); ++i )
        mf::LogWarning("GENIEHelper")
          << "strtok[" << i << "] = \"" << strtok[i] << "\"";
      return;
    }

    string stype = strtok[0];

    // parse out values
    vector<double> vals;
    vector<string> valstrs = genie::utils::str::Split(strtok[1]," ,;(){}[]\t\n\r");
    vector<string>::const_iterator iter = valstrs.begin();
    for ( ; iter != valstrs.end(); ++iter ) {
      const string& valstr1 = *iter;
      if ( valstr1 != "" ) {
        double aval = atof(valstr1.c_str());
        mf::LogDebug("GENIEHelper") << "rock value [" << vals.size() << "] "
                                    << aval;
        vals.push_back(aval);
      }
    }
    size_t nvals = vals.size();

    rocksel->SetRemoveEntries(true);  // drop segments that won't be considered

    // assume coordinates are in the *master* (not "top volume") system
    // need to set fTopVolume to fWorldVolume as Sample() will keep setting it
    fTopVolume = fWorldVolume;
    rgeom->SetTopVolName(fTopVolume.c_str());

    if ( nvals < 6 ) {
      throw cet::exception("GENIEHelper") << "rockbox needs at "
                                          << "least 6 values, found "
                                          << nvals << "in \""
                                          << strtok[1] << "\"";

    }
    double xyzmin[3] = { vals[0], vals[1], vals[2] };
    double xyzmax[3] = { vals[3], vals[4], vals[5] };

    bool   rockonly  = true;
    double wallmin   = 800.;   // geometry in cm, ( 8 meter buffer)
    double dedx      = 2.5 * 1.7e-3; // GeV/cm, rho=2.5, 1.7e-3 ~ rock like loss
    double fudge     = 1.05;

    if ( nvals >=  7 ) rockonly = vals[6];
    if ( nvals >=  8 ) wallmin  = vals[7];
    if ( nvals >=  9 ) dedx     = vals[8];
    if ( nvals >= 10 ) fudge    = vals[9];

    rocksel->SetRockBoxMinimal(xyzmin,xyzmax);
    rocksel->SetMinimumWall(wallmin);
    rocksel->SetDeDx(dedx/fudge);

    // if not rock-only then make a tiny exclusion bubble
    // call to MakeBox shouldn't be necessary
    //  should be done by SetRockBoxMinimal but for some GENIE versions isn't
    if ( ! rockonly ) rocksel->MakeSphere(0,0,0,1.0e-10);
    else              rocksel->MakeBox(xyzmin,xyzmax);

    rgeom->AdoptGeomVolSelector(rocksel);
  }

void evgb::GENIEHelper::ReadXSecTable ( )

private

determine which cross section table to use fully expand the path

Definition at line 2449 of file GENIEHelper.cxx.

References fEnvironment, fGXMLPATH, and fXSecTable.

Referenced by Initialize().

  {

    // priority order:
    //    fcl fEnvironment GSPLOAD
    //    fcl XSecTable
    //    $GSPLOAD in environment
    //    default 'gxspl-FNALsmall.xml'

    if ( fXSecTable == "" ) {
      // stand-alone value is not set
      const char* gspload_alt = std::getenv("GSPLOAD");
      if ( ! gspload_alt ) {
        const char* gspload_dflt = "gxspl-FNALsmall.xml";  // fall back
        gspload_alt = gspload_dflt;
      } else {
        throw cet::exception("$GSPLOAD")
          << "using env variable $GSPLOAD: " << gspload_alt
          << ", use fcl parameter 'XSecTable' instead.";
      }
      fXSecTable = std::string(gspload_alt);
    }

    // find GSPLOAD in the vector, if it exists
    int indxGSPLOAD   = -1;
    for (size_t i = 0; i < fEnvironment.size(); i += 2) {
      if ( fEnvironment[i].find("GSPLOAD") == 0 ) {
        indxGSPLOAD = i;
        throw cet::exception("UsingGSPLOAD")
          << "using Environment fcl parameter GSPLOAD: "
          << fEnvironment[indxGSPLOAD+1]
          << ", use fcl parameter 'XSecTable' instead. "
          << __FILE__ << ":" << __LINE__ << "\n";
        continue;
      }
    }

    if ( indxGSPLOAD < 0 ) {
      // nothing in fcl parameters
      indxGSPLOAD=fEnvironment.size();
      fEnvironment.push_back("GSPLOAD");
      fEnvironment.push_back(fXSecTable);
    } else {
      fXSecTable = fEnvironment[indxGSPLOAD+1];  // get two in sync
    }

    // currently GENIE doesn't internally use GXMLPATH when looking for
    // spline files, but instead wants a fully expanded path.
    // Do the expansion here using the extended GXMLPATH list
    // of locations (which included $FW_SEARCH_PATH)
    mf::LogDebug("GENIEHelper") << "GSPLOAD as originally: " << fXSecTable;

    // RWH cet::search_path returns "" if the input string is actually
    // the full path to the file .. this is not really what one wants,
    // one just wan't the full path to the file; seems to work if
    // "/" is made to be another possible PATH
    cet::search_path spGXML( "/:" + fGXMLPATH );
    std::string fullpath;
    spGXML.find_file(fXSecTable, fullpath);

    if ( fullpath == "" ) {
      mf::LogError("GENIEHelper")
        << "could not resolve full path for spline file XSecTable/GSPLOAD "
        << "\"" << fXSecTable << "\" using: "
        << fGXMLPATH;
      throw cet::exception("UnresolvedGSPLOAD")
        << "can't find XSecTable/GSPLOAD file: " << fXSecTable;
    }
    fXSecTable = fullpath;
    fEnvironment[indxGSPLOAD+1] = fXSecTable;  // get two in sync

    mf::LogInfo("GENIEHelper")
      << "XSecTable/GSPLOAD full path \"" << fXSecTable << "\"";

    TStopwatch xtime;
    xtime.Start();

    // can't use gSystem->Unsetenv() as it is really gSystem->Setenv(name,"")
    unsetenv("GSPLOAD");  // MUST!!! ensure that it isn't set externally
    genie::utils::app_init::XSecTable(fXSecTable,true);

    xtime.Stop();
    mf::LogInfo("GENIEHelper")
      << "Time to read GENIE XSecTable: "
      << " Real " << xtime.RealTime() << " s,"
      << " CPU " << xtime.CpuTime() << " s"
      << " from " << fXSecTable;

  }

void evgb::GENIEHelper::RegularizeFluxType ( )

private

Definition at line 717 of file GENIEHelper.cxx.

References fFluxType.

Referenced by GENIEHelper().

  {
    // Regularize fFluxType string to sensible setting

    std::string tmpFluxType = fFluxType;

    // remove lead/trailing whitespace
    size_t ftlead  = tmpFluxType.find_first_not_of(" \t\n");
    if ( ftlead )    tmpFluxType.erase( 0, ftlead );
    size_t ftlen   = tmpFluxType.length();
    size_t fttrail = tmpFluxType.find_last_not_of(" \t\n");
    if ( fttrail != ftlen ) tmpFluxType.erase( fttrail+1, ftlen );

    // strip off leading catagories ... we'll put them back later
    // so we don't accidently allow arbitrary strings
    if ( tmpFluxType.find("atmos_") == 0 ) tmpFluxType.erase(0,6);
    if ( tmpFluxType.find("atmo_")  == 0 ) tmpFluxType.erase(0,5);
    if ( tmpFluxType.find("tree_")  == 0 ) tmpFluxType.erase(0,5);

    // make reasonable inferences of what the user intended

    // simple fluxes
    if ( tmpFluxType.find("hist")   != std::string::npos ) tmpFluxType = "histogram";
    if ( tmpFluxType.find("func")   != std::string::npos ) tmpFluxType = "function";
    if ( tmpFluxType.find("mono")   != std::string::npos ) tmpFluxType = "mono";
    // Atmospheric fluxes
    // prior to R-2_11_0 BGLRS was "BARTOL" and HAKKM was "HONDA"
    if ( tmpFluxType.find("FLUKA")  != std::string::npos ) tmpFluxType = "atmo_FLUKA";
    if ( tmpFluxType.find("BARTOL") != std::string::npos ) tmpFluxType = "atmo_BGLRS";
    if ( tmpFluxType.find("BGLRS")  != std::string::npos ) tmpFluxType = "atmo_BGLRS";
    if ( tmpFluxType.find("HONDA")  != std::string::npos ) tmpFluxType = "atmo_HAKKM";
    if ( tmpFluxType.find("HAKKM")  != std::string::npos ) tmpFluxType = "atmo_HAKKM";
    if ( tmpFluxType.find("POWER")  != std::string::npos ) tmpFluxType = "atmo_PowerSpectrum";
    // TTree-based fluxes (old "ntuple" is really "numi")
    //    we're allowed to randomize the order here, and squeeze out duplicates
    if ( tmpFluxType.find("simple") != std::string::npos ) tmpFluxType = "tree_simple";
    if ( tmpFluxType.find("ntuple") != std::string::npos ) tmpFluxType = "tree_numi";
    if ( tmpFluxType.find("numi")   != std::string::npos ) tmpFluxType = "tree_numi";
    if ( tmpFluxType.find("dk2nu")  != std::string::npos ) tmpFluxType = "tree_dk2nu";

    fFluxType = tmpFluxType;
  }

bool evgb::GENIEHelper::Sample	(	simb::MCTruth &	truth,
		simb::MCFlux &	flux,
		simb::GTruth &	gtruth
	)

Definition at line 1662 of file GENIEHelper.cxx.

References bin, simb::MCParticle::E(), fAddGenieVtxTime, fAtmoSpectralIndex, fDebugFlags, simb::MCFlux::fdk2gen, fDriver, fEventsPerSpill, fFluxD, fFluxD2GMCJD, fFluxHistograms, simb::MCFlux::fFluxType, fFluxType, fForceApplyFlxWgt, simb::MCFlux::fgen2vtx, fGenFlavors, fGenieEventRecord, simb::MCFlux::fgenx, simb::MCFlux::fgeny, simb::MCFlux::fgenz, fGeoManager, fGHepPrintLevel, fGlobalTimeOffset, fHelperRandom, evgb::FillGTruth(), evgb::FillMCFlux(), evgb::FillMCTruth(), fRandomTimeOffset, fSpillEvents, fSpillExposure, fTimeShifter, fTopVolume, fTotalExposure, fTuneName, fUseBlenderDist, fUseHelperRndGen4GENIE, fWorldVolume, simb::MCTruth::GetNeutrino(), simb::kHistPlusFocus, evgb::kNue, evgb::kNueBar, evgb::kNuMu, evgb::kNuMuBar, evgb::kNuTau, evgb::kNuTauBar, simb::MCNeutrino::Nu(), simb::MCFlux::SetFluxGen(), evgb::EvtTimeShiftI::TimeOffset(), and util::to_string().

Referenced by evgen::TestGENIEHelper::produce(), and evgen::GENIEGen::produce().

  {
    // set the top volume for the geometry
    fGeoManager->SetTopVolume(fGeoManager->FindVolumeFast(fTopVolume.c_str()));

    if ( fGenieEventRecord ) delete fGenieEventRecord;

    // ART Framework plays games with gRandom, undo that if requested
    TRandom* old_gRandom = gRandom;
    if (fUseHelperRndGen4GENIE) gRandom = fHelperRandom;

    fGenieEventRecord = fDriver->GenerateEvent();

    if (fForceApplyFlxWgt) {
      //const issue:  double flxweight = fDriver->FluxDriver().Weight();
      //genie::GFluxI& flx = const_cast<genie::GFluxI&>(fDriver->FluxDriver());
      //double flxweight = flx.Weight();
      // use the flux driver handle we already have, rather than fetch const one
      //   from the GMCJDriver fDriver
      double flxweight = fFluxD->Weight();
      double curweight = fGenieEventRecord->Weight();
      fGenieEventRecord->SetWeight(flxweight*curweight);
    }

    if (fUseHelperRndGen4GENIE) gRandom = old_gRandom;

    // now check if we produced a viable event record
    bool viableInteraction = true;
    if ( ! fGenieEventRecord ) viableInteraction = false;

    // update the spill total information, then check to see
    // if we got an event record that was valid

    genie::flux::GFluxExposureI* fexposure =
      dynamic_cast<genie::flux::GFluxExposureI*>(fFluxD);
    if ( fexposure ) {
      fSpillExposure =
        (fexposure->GetTotalExposure()/fDriver->GlobProbScale()) - fTotalExposure;
    }
    // use GENIE2ART code to fill MCFlux
    evgb::FillMCFlux(fFluxD,flux);

    // if no interaction generated return false
    if ( ! viableInteraction ) return false;

    // fill the MCTruth & GTruth information as we have a good interaction
    // these two objects are enough to reconstruct the GENIE record
    // use the new external functions in GENIE2ART

    // choose a time within the spill (ns) to shift the vertex times by:
    double timeoffset = 0;
    if ( ! fTimeShifter ) {
      timeoffset = fHelperRandom->Uniform()*fRandomTimeOffset;
    } else {
      timeoffset = fTimeShifter->TimeOffset();
    }
    // mf::LogInfo("GENIEHelper") << "TimeShifter adding " << timeoffset;
    double spilltime  = fGlobalTimeOffset + timeoffset;

    std::unordered_map<std::string, std::string> genConfig;

    if(fFluxType.find("PowerSpectrum") != std::string::npos){
      genConfig.emplace("SpectralIndex", std::to_string(fAtmoSpectralIndex));
    }

    evgb::FillMCTruth(fGenieEventRecord, spilltime, truth,
                      __GENIE_RELEASE__, fTuneName, fAddGenieVtxTime );
    evgb::FillGTruth(fGenieEventRecord, gtruth);

    // check to see if we are using flux ntuples but want to
    // make n events per spill
    if ( ( fEventsPerSpill > 0 ) &&
         ( fFluxType.find("tree_") == 0 ) ) ++fSpillEvents;

    // now check if using either histogram or mono fluxes, using
    // either n events per spill or basing events on POT per spill for the
    // histogram case
    if ( fFluxType.find("histogram") == 0 ) {
      // set the flag in the parent object that says the
      // fluxes came from histograms and fill related values
      flux.fFluxType = simb::kHistPlusFocus;

      // save the fluxes - fluxes were added to the vector in the same
      // order that the flavors appear in fGenFlavors
      int ctr = 0;
      int bin = fFluxHistograms[0]->FindBin(truth.GetNeutrino().Nu().E());
      std::vector<double> fluxes(6, 0.);
      for(std::vector<int>::iterator i = fGenFlavors.begin(); i != fGenFlavors.end(); i++){
        if(*i ==  12) fluxes[kNue]      = fFluxHistograms[ctr]->GetBinContent(bin);
        if(*i == -12) fluxes[kNueBar]   = fFluxHistograms[ctr]->GetBinContent(bin);
        if(*i ==  14) fluxes[kNuMu]     = fFluxHistograms[ctr]->GetBinContent(bin);
        if(*i == -14) fluxes[kNuMuBar]  = fFluxHistograms[ctr]->GetBinContent(bin);
        if(*i ==  16) fluxes[kNuTau]    = fFluxHistograms[ctr]->GetBinContent(bin);
        if(*i == -16) fluxes[kNuTauBar] = fFluxHistograms[ctr]->GetBinContent(bin);
        ++ctr;
      }

      // get the flux for each neutrino flavor of this energy
      flux.SetFluxGen(fluxes[kNue],   fluxes[kNueBar],
                      fluxes[kNuMu],  fluxes[kNuMuBar],
                      fluxes[kNuTau], fluxes[kNuTauBar]);

      ++fSpillEvents;
    }
    else if ( fFluxType.find("mono")     == 0 ||
              fFluxType.find("function") == 0    ){
      ++fSpillEvents;
    }
    else if ( fFluxType.find("atmo_FLUKA")  == 0 ||
              fFluxType.find("atmo_BARTOL") == 0 ||
              fFluxType.find("atmo_PowerSpectrum") == 0 ||
              fFluxType.find("atmo_BGLRS")  == 0 ||
              fFluxType.find("atmo_HAKKM")  == 0 ||
              fFluxType.find("atmo_HONDA")  == 0    ) {
      if ( fEventsPerSpill > 0 ) ++fSpillEvents;
      flux.fFluxType = simb::kHistPlusFocus;
    }


    // fill these after the Pack[NuMI|Simple]Flux because those
    // will Reset() the values at the start
    TLorentzVector *vertex = fGenieEventRecord->Vertex();
    TLorentzVector nuray_pos = fFluxD->Position();
    TVector3 ray2vtx = nuray_pos.Vect() - vertex->Vect();
    flux.fgenx    = nuray_pos.X();
    flux.fgeny    = nuray_pos.Y();
    flux.fgenz    = nuray_pos.Z();
    flux.fgen2vtx = ray2vtx.Mag();

    genie::flux::GFluxBlender* blender =
      dynamic_cast<genie::flux::GFluxBlender*>(fFluxD2GMCJD);
    if ( blender ) {
      if ( fUseBlenderDist ) flux.fdk2gen = blender->TravelDist();
      // / if mixing flavors print the state of the blender
      if ( fDebugFlags & 0x02 ) blender->PrintState();
    }

    if ( fDebugFlags & 0x04 ) {
      mf::LogInfo("GENIEHelper") << "vertex loc " << vertex->X() << ","
                                 << vertex->Y() << "," << vertex->Z() << std::endl
                                 << " flux ray start " << nuray_pos.X() << ","
                                 << nuray_pos.Y() << "," << nuray_pos.Z() << std::endl
                                 << " ray2vtx = " << flux.fgen2vtx
                                 << " dk2ray = " << flux.fdk2gen;
    }
    if ( fGHepPrintLevel >= 0 ) {
      std::cout << *fGenieEventRecord;
    }

    // set the top volume of the geometry back to the world volume
    fGeoManager->SetTopVolume(fGeoManager->FindVolumeFast(fWorldVolume.c_str()));

    return true;
  }

void evgb::GENIEHelper::SetGMSGLAYOUT ( )

private

GMSGLAYOUT ([BASIC}|SIMPLE) control GENIE's layout of log4cpp message SIMPLE lacks the timestamp; this must be set in the environment at the time the log4cpp Messenger singleton is created

Definition at line 2359 of file GENIEHelper.cxx.

References fEnvironment, and fGMSGLAYOUT.

Referenced by GENIEHelper().

  {

    // start with GMSGLAYOUT set from pset "GMSGLAYOUT" value

    // find it in the vector, if it exists
    // this will override the top level fcl value
    for (size_t i = 0; i < fEnvironment.size(); i += 2) {
      if ( fEnvironment[i].find("GMSGLAYOUT") == 0 ) {
        fGMSGLAYOUT = fEnvironment[i+1];
        break;
      }
    }

    // now set it externally, if we have a value
    if ( fGMSGLAYOUT != "" ) {
      gSystem->Setenv("GMSGLAYOUT",fGMSGLAYOUT.c_str());
    }
  }

void evgb::GENIEHelper::SetGXMLPATH ( )

private

GXMLPATH is where GENIE will look for alternative XML configurations (including message service threshold files)

Definition at line 2300 of file GENIEHelper.cxx.

References fEnvironment, and fGXMLPATH.

Referenced by GENIEHelper().

  {

    // priority order
    //   (fcl file paths):(existing user environment):(FW_SEARCH_PATH)
    //
    //  fcl parameters might be the explicit GXMLPATH value
    //     or a pair in the fEnvironment vector
    //  but the final "official result" should the fGXMLPATH

    // start with fGXMLPATH set from pset "GXMLPATH" value

    // find it in the vector, if it exists
    int indxGXMLPATH = -1;
    for (size_t i = 0; i < fEnvironment.size(); i += 2) {
      if ( fEnvironment[i].find("GXMLPATH") == 0 ) {
        if ( fGXMLPATH != "" ) fGXMLPATH += ":";
        fGXMLPATH += fEnvironment[i+1];
        indxGXMLPATH = i;
        /*
        throw cet::exception("UsingGXMLPATH")
          << "using Environment fcl parameter GXMLPATH: "
          << fEnvironment[indxGXMLPATH+1]
          << ", use fcl parameter GXMLPATH instead.";
        */
        break;
      }
    }

    const char* gxmlpath_env = std::getenv("GXMLPATH");
    if ( gxmlpath_env ) {
      if ( fGXMLPATH != "" ) fGXMLPATH += ":";
      fGXMLPATH += gxmlpath_env;
    }
    const char* fwpath_env = std::getenv("FW_SEARCH_PATH");
    if ( fwpath_env ) {
      if ( fGXMLPATH != "" ) fGXMLPATH += ":";
      fGXMLPATH += fwpath_env;
    }

    // refresh fEnvironment (in case fEnvironment is still being used)
    if ( indxGXMLPATH < 0 ) {
       // nothing in fcl parameters
      indxGXMLPATH=fEnvironment.size();
      fEnvironment.push_back("GXMLPATH");
      fEnvironment.push_back(fGXMLPATH);
    } else {
      // update value
      fEnvironment[indxGXMLPATH+1] = fGXMLPATH;
    }

    // now set it externally for use by GENIE
    gSystem->Setenv("GXMLPATH",fGXMLPATH.c_str());

  }

void evgb::GENIEHelper::SetMaxPathOutInfo ( )

private

Definition at line 1571 of file GENIEHelper.cxx.

References fBeamName, fDetLocation, fFiducialCut, fFluxType, fGeoFile, fGeomScan, fMaxPathOutInfo, fSelectedFluxFiles, fTopVolume, and fWorldVolume.

Referenced by ConfigGeomScan().

  {
    // create an info string based on:
    // ROOT geometry, TopVolume, FiducialCut, GeomScan, Flux

    mf::LogInfo("GENIEHelper") << "about to create MaxPathOutInfo";

    fMaxPathOutInfo = "\n";
    fMaxPathOutInfo += "   FluxType:     " + fFluxType + "\n";
    fMaxPathOutInfo += "   BeamName:     " + fBeamName + "\n";
    fMaxPathOutInfo += "   FluxFiles:    ";
    std::vector<string>::iterator ffitr = fSelectedFluxFiles.begin();
    for ( ; ffitr != fSelectedFluxFiles.end() ; ++ffitr )
      fMaxPathOutInfo += "\n         " + *ffitr;
    fMaxPathOutInfo += "\n";
    fMaxPathOutInfo += "   DetLocation:  " + fDetLocation + "\n";
    fMaxPathOutInfo += "   ROOTFile:     " + fGeoFile     + "\n";
    fMaxPathOutInfo += "   WorldVolume:  " + fWorldVolume + "\n";
    fMaxPathOutInfo += "   TopVolume:    " + fTopVolume   + "\n";
    fMaxPathOutInfo += "   FiducialCut:  " + fFiducialCut + "\n";
    fMaxPathOutInfo += "   GeomScan:     " + fGeomScan    + "\n";

    mf::LogInfo("GENIEHelper") << "MaxPathOutInfo: \""
                               << fMaxPathOutInfo << "\"";

  }

double evgb::GENIEHelper::SpillExposure ( ) const

inline

Definition at line 75 of file GENIEHelper.h.

{ return fSpillExposure;  }

void evgb::GENIEHelper::SqueezeFilePatterns ( )

private

Definition at line 761 of file GENIEHelper.cxx.

References fFluxFilePatterns.

Referenced by GENIEHelper().

  {
    // for "numi" (nee "ntuple"), "simple" and "dk2nu" squeeze the patterns
    // so there are no duplicates; for the others we want to preserve order

    // convert vector<> to a set<> and back to vector<>
    // to avoid having duplicate patterns in the list
    std::set<std::string> fluxpattset(fFluxFilePatterns.begin(),
                                      fFluxFilePatterns.end());
    //std::copy(fFluxFilePatterns.begin(),fFluxFilePatterns.end(),
    //          std::inserter(fluxpattset,fluxpattset.begin()));
    fFluxFilePatterns.clear(); // clear vector, copy unique set back
    std::copy(fluxpattset.begin(),fluxpattset.end(),
              std::back_inserter(fFluxFilePatterns));
  }

void evgb::GENIEHelper::StartGENIEMessenger ( std::string  prodmode )

private

start with fGENIEMsgThresholds from pset "GENIEMsgThresholds" value allow fEnvironment $GMSGCONF and $GPRODMODE to expand it function arg might also trigger addition of Messenger_production.xml (pre-R-2_9_0) or Messenger_whisper.xml

$GPRODMODE used to trigger Messenger_production.xml with R-2_8_0 one must add it explicitly to $GMSGCONF

Definition at line 2383 of file GENIEHelper.cxx.

References fEnvironment, fGENIEMsgThresholds, and StringToBool().

Referenced by GENIEHelper().

  {


    int indxGPRODMODE = -1;
    int indxGMSGCONF  = -1;

    for (size_t i = 0; i < fEnvironment.size(); i += 2) {
      if ( fEnvironment[i].find("GPRODMODE") == 0 ) {
        indxGPRODMODE = i;
        continue;
      }
      if ( fEnvironment[i].find("GMSGCONF") == 0 ) {
        indxGMSGCONF = i;
        continue;
      }
    }
    // GMSGCONF set in fEnvironment tack it on to current setting
    if ( indxGMSGCONF >= 0 ) {
      if ( fGENIEMsgThresholds != "" ) fGENIEMsgThresholds += ":";
      fGENIEMsgThresholds += fEnvironment[indxGMSGCONF+1];
    } else {
      // nothing in fcl parameters, make a placeholder
      indxGMSGCONF = fEnvironment.size();
      fEnvironment.push_back("GMSGCONF");
      fEnvironment.push_back("");
    }

    bool prodmode = StringToBool(prodmodestr);
    if ( indxGPRODMODE >= 0 ) {
      // user tried to set GPRODMODE via Environment fcl values
      prodmode |= StringToBool(fEnvironment[indxGPRODMODE+1]);
    }

    if ( prodmode ) {
      // okay they really meant it
      // PREpend "Messenger_whisper.xml" to existing value
      // don't forget the colon if necessary
      std::string newval = "Messenger_whisper.xml";
      if ( fGENIEMsgThresholds != "" ) {
        newval += ":";
        newval += fGENIEMsgThresholds;
      }
      fGENIEMsgThresholds = newval;
    }

    // update fEnvironment value
    if ( indxGMSGCONF >= 0 ) {
      fEnvironment[indxGMSGCONF+1] = fGENIEMsgThresholds;
    }

    mf::LogInfo("GENIEHelper")
      << "StartGENIEMessenger ProdMode=" << ((prodmode)?"yes":"no")
      << " read from: " << fGENIEMsgThresholds;

    genie::utils::app_init::MesgThresholds(fGENIEMsgThresholds);

  }

bool evgb::GENIEHelper::Stop

(

)

Definition at line 1599 of file GENIEHelper.cxx.

References fAtmoRt, fDriver, fEventsPerSpill, fFluxD, fFluxType, fGenFlavors, fHelperRandom, fHistEventsPerSpill, fPOTPerSpill, fSpillEvents, fSpillExposure, fTotalExposure, fTotalHistFlux, and fXSecMassPOT.

Referenced by evgen::TestGENIEHelper::produce(), and evgen::GENIEGen::produce().

  {
    //   std::cout << "in GENIEHelper::Stop(), fEventsPerSpill = " << fEventsPerSpill
    //      << " fPOTPerSpill = " << fPOTPerSpill << " fSpillExposure = " << fSpillExposure
    //      << " fSpillEvents = " << fSpillEvents
    //      << " fHistEventsPerSpill = " << fHistEventsPerSpill << std::endl;

    // determine if we should keep throwing neutrinos for
    // this spill or move on

    if ( fEventsPerSpill > 0 ) {
      if ( fSpillEvents < fEventsPerSpill ) return false;
    } else {
      // exposure (POT) based
      if ( fFluxType.find("tree_") == 0 ) {
        if ( fSpillExposure < fPOTPerSpill ) {
          return false;
        }
      }
      else if ( fFluxType.find("histogram") == 0 ) {
        if ( fSpillEvents < fHistEventsPerSpill ) return false;
        else fSpillExposure = fPOTPerSpill;
      }
    }

    if ( fFluxType.find("atmo_")  == 0 ) {
      // the exposure for atmo is in SECONDS. In order to get seconds,
      // it needs to be normalized by 1e4 to take into account the units
      // discrepency between AtmoFluxDriver(/m2) and Generate(/cm2)
      // and it need to be normalized by the generation surface area since
      // it's not taken into accoutn in the flux driver

      long int nNeutrinos;

      if(fFluxType.find("PowerSpectrum") != std::string::npos){
        nNeutrinos = dynamic_cast<genie::flux::GPowerSpectrumAtmoFlux *>(fFluxD)->NFluxNeutrinos();
        fTotalExposure = nNeutrinos/fGenFlavors.size()/fDriver->GlobProbScale();
        mf::LogInfo("GENIEHelper")
        << "===> Atmo Pscale*Ngen/Nflavours = " << fTotalExposure;
      }
      else{
        nNeutrinos = dynamic_cast<genie::flux::GAtmoFlux *>(fFluxD)->NFluxNeutrinos();
        fTotalExposure = nNeutrinos * 1.0e4 / (TMath::Pi() * fAtmoRt*fAtmoRt);
        mf::LogInfo("GENIEHelper")
        << "===> Atmo EXPOSURE = " << fTotalExposure << " seconds.";
      }

      

      

    } else {
      fTotalExposure += fSpillExposure;
    }

    // made it to here, means need to reset the counters
    fSpillEvents   = 0;
    fSpillExposure = 0.;
    fHistEventsPerSpill = fHelperRandom->Poisson(fXSecMassPOT*fTotalHistFlux);
    return true;
  }

bool evgb::GENIEHelper::StringToBool ( std::string  v )

private

Definition at line 2542 of file GENIEHelper.cxx.

Referenced by StartGENIEMessenger().

  {
    if (v == "true")   return true;  // C++ style
    if (v == "false")  return false;
    if (v == "kTRUE")  return true;  // ROOT style
    if (v == "kFALSE") return false;
    if (v == "TRUE")   return true;  // Some other reasonable variations...
    if (v == "FALSE")  return false;
    if (v == "True")   return true;
    if (v == "False")  return false;
    if (v == "on")     return true;
    if (v == "off")    return false;
    if (v == "On")     return true;
    if (v == "Off")    return false;
    if (v == "ON")     return true;
    if (v == "OFF")    return false;
    if (v == "YES")    return true;
    if (v == "NO")     return false;
    if (v == "Yes")    return true;
    if (v == "No")     return false;
    if (v == "yes")    return true;
    if (v == "no")     return false;
    if (v == "1")      return true;
    if (v == "0")      return false;

    return false;  // by default invalid strings are false
  }

double evgb::GENIEHelper::TotalExposure ( ) const

inline

Definition at line 71 of file GENIEHelper.h.

Referenced by evgen::GENIEGen::beginSubRun(), evgen::TestGENIEHelper::endSubRun(), and evgen::GENIEGen::endSubRun().

{ return fTotalExposure;  }

double evgb::GENIEHelper::TotalHistFlux

(

)

Definition at line 591 of file GENIEHelper.cxx.

References fFluxType, and fTotalHistFlux.

  {
    if ( fFluxType.find("mono")     == 0 ||
         fFluxType.find("function") == 0 ||
         fFluxType.find("tree_")    == 0 ||
         fFluxType.find("atmo_")    == 0    ) {
      // shouldn't be asking for this for any of the above
      // perhaps not-not "function"?
      return -999.;
    }

    return fTotalHistFlux;
  }

double evgb::GENIEHelper::TotalMass ( ) const

inline

Definition at line 82 of file GENIEHelper.h.

{ return fDetectorMass+fSurroundingMass; }

Member Data Documentation

bool evgb::GENIEHelper::fAddGenieVtxTime

private

incorporate time from flux window to interaction point and (possibily) proton-on-target to flux window

Definition at line 183 of file GENIEHelper.h.

Referenced by GENIEHelper(), and Sample().

double evgb::GENIEHelper::fAtmoEmax

private

atmo: Maximum energy of neutrinos in GeV

Definition at line 188 of file GENIEHelper.h.

Referenced by AtmoFluxCheck(), GENIEHelper(), and InitializeFluxDriver().

double evgb::GENIEHelper::fAtmoEmin

private

atmo: Minimum energy of neutrinos in GeV

Definition at line 187 of file GENIEHelper.h.

Referenced by AtmoFluxCheck(), GENIEHelper(), and InitializeFluxDriver().

double evgb::GENIEHelper::fAtmoRl

private

atmo: radius of the sphere on where the neutrinos are generated

Definition at line 189 of file GENIEHelper.h.

Referenced by AtmoFluxCheck(), GENIEHelper(), and InitializeFluxDriver().

double evgb::GENIEHelper::fAtmoRt

private

atmo: radius of the transvere (perpendicular) area on the sphere where the neutrinos are generated

Definition at line 190 of file GENIEHelper.h.

Referenced by AtmoFluxCheck(), GENIEHelper(), InitializeFluxDriver(), and Stop().

double evgb::GENIEHelper::fAtmoSpectralIndex

private

atmo: Spectral index for power spectrum generation

Definition at line 192 of file GENIEHelper.h.

Referenced by GENIEHelper(), InitializeFluxDriver(), and Sample().

TVector3 evgb::GENIEHelper::fBeamCenter

private

center of beam for histogram fluxes - must be in meters

Definition at line 175 of file GENIEHelper.h.

Referenced by GENIEHelper(), HistogramFluxCheck(), and InitializeFluxDriver().

TVector3 evgb::GENIEHelper::fBeamDirection

private

direction of the beam for histogram fluxes

Definition at line 174 of file GENIEHelper.h.

Referenced by GENIEHelper(), HistogramFluxCheck(), and InitializeFluxDriver().

std::string evgb::GENIEHelper::fBeamName

private

name of the beam we are simulating

Definition at line 149 of file GENIEHelper.h.

Referenced by SetMaxPathOutInfo().

double evgb::GENIEHelper::fBeamRadius

private

radius of cylindar for histogram fluxes - must be in meters

Definition at line 176 of file GENIEHelper.h.

Referenced by GENIEHelper(), HistogramFluxCheck(), and InitializeFluxDriver().

unsigned int evgb::GENIEHelper::fDebugFlags

private

set bits to enable debug info

Definition at line 208 of file GENIEHelper.h.

Referenced by GENIEHelper(), InitializeFluxDriver(), InitializeGeometry(), and Sample().

double evgb::GENIEHelper::fDetectorMass

private

mass of the detector in kg

Definition at line 177 of file GENIEHelper.h.

Referenced by GENIEHelper(), and Initialize().

std::string evgb::GENIEHelper::fDetLocation

private

name of flux window location

Definition at line 156 of file GENIEHelper.h.

Referenced by GENIEHelper(), InitializeFluxDriver(), and SetMaxPathOutInfo().

genie::GMCJDriver* evgb::GENIEHelper::fDriver

private

Definition at line 130 of file GENIEHelper.h.

Referenced by ConfigGeomScan(), Initialize(), Sample(), Stop(), and ~GENIEHelper().

double evgb::GENIEHelper::fEmax

private

Definition at line 171 of file GENIEHelper.h.

Referenced by GENIEHelper(), and InitializeFluxDriver().

double evgb::GENIEHelper::fEmin

private

Definition at line 170 of file GENIEHelper.h.

Referenced by GENIEHelper(), and InitializeFluxDriver().

std::vector<std::string> evgb::GENIEHelper::fEnvironment

private

environmental variables and settings used by genie

Definition at line 194 of file GENIEHelper.h.

Referenced by GENIEHelper(), ReadXSecTable(), SetGMSGLAYOUT(), SetGXMLPATH(), and StartGENIEMessenger().

std::string evgb::GENIEHelper::fEventGeneratorList

private

control over event topologies, was $GEVGL [Default]

Definition at line 197 of file GENIEHelper.h.

Referenced by GENIEHelper(), and Initialize().

double evgb::GENIEHelper::fEventsPerSpill

private

number of events to generate in each spill if not using POT/spill. If using Atmo, set to 1

Definition at line 160 of file GENIEHelper.h.

Referenced by AtmoFluxCheck(), GENIEHelper(), Initialize(), Sample(), and Stop().

std::string evgb::GENIEHelper::fFiducialCut

private

configuration for geometry selector

Definition at line 205 of file GENIEHelper.h.

Referenced by GENIEHelper(), InitializeFiducialSelection(), InitializeRockBoxSelection(), and SetMaxPathOutInfo().

std::string evgb::GENIEHelper::fFluxCleanup

private

"ALWAYS", "/var/tmp", "NEVER"

Definition at line 148 of file GENIEHelper.h.

Referenced by ~GENIEHelper().

std::string evgb::GENIEHelper::fFluxCopyMethod

private

"DIRECT" = old direct access method, otherwise = ifdh approach schema ("" okay)

Definition at line 147 of file GENIEHelper.h.

Referenced by ExpandFluxFilePatternsIFDH(), ExpandFluxPaths(), and GENIEHelper().

genie::GFluxI* evgb::GENIEHelper::fFluxD

private

real flux driver

Definition at line 128 of file GENIEHelper.h.

Referenced by ConfigGeomScan(), Initialize(), InitializeFluxDriver(), Sample(), Stop(), and ~GENIEHelper().

genie::GFluxI* evgb::GENIEHelper::fFluxD2GMCJD

private

flux driver passed to genie GMCJDriver, might be GFluxBlender

Definition at line 129 of file GENIEHelper.h.

Referenced by Initialize(), InitializeFluxDriver(), and Sample().

std::vector<std::string> evgb::GENIEHelper::fFluxFilePatterns

private

wildcard patterns files containing histograms or ntuples, or txt

Definition at line 143 of file GENIEHelper.h.

Referenced by ExpandFluxFilePatternsDirect(), ExpandFluxFilePatternsIFDH(), and SqueezeFilePatterns().

std::vector<TH1D *> evgb::GENIEHelper::fFluxHistograms

private

histograms for each nu species

Definition at line 157 of file GENIEHelper.h.

Referenced by HistogramFluxCheck(), InitializeFluxDriver(), and Sample().

TRotation* evgb::GENIEHelper::fFluxRotation

private

rotation for atmos / astro flux coord systems

Definition at line 152 of file GENIEHelper.h.

Referenced by BuildFluxRotation(), GENIEHelper(), and InitializeFluxDriver().

std::string evgb::GENIEHelper::fFluxRotCfg

private

how to interpret fFluxRotValues

Definition at line 150 of file GENIEHelper.h.

Referenced by BuildFluxRotation().

std::vector<double> evgb::GENIEHelper::fFluxRotValues

private

parameters for rotation

Definition at line 151 of file GENIEHelper.h.

Referenced by BuildFluxRotation().

std::string evgb::GENIEHelper::fFluxSearchPaths

private

colon separated set of path stems

Definition at line 142 of file GENIEHelper.h.

Referenced by ExpandFluxFilePatternsDirect(), ExpandFluxFilePatternsIFDH(), and ExpandFluxPaths().

std::string evgb::GENIEHelper::fFluxType

private

histogram, gsimple, dk2nu, ntuple/gnumi, atmos_XXXX atmo_{FLUKA|BARTOL/BGLRS|HONDA/HAKKM}

Definition at line 139 of file GENIEHelper.h.

Referenced by AtmoFluxCheck(), ExpandFluxFilePatternsDirect(), ExpandFluxFilePatternsIFDH(), GENIEHelper(), Initialize(), InitializeFluxDriver(), RegularizeFluxType(), Sample(), SetMaxPathOutInfo(), Stop(), and TotalHistFlux().

double evgb::GENIEHelper::fFluxUpstreamZ

private

z where flux starts from (if non-default, simple/ntuple only)

Definition at line 159 of file GENIEHelper.h.

Referenced by GENIEHelper(), and InitializeFluxDriver().

bool evgb::GENIEHelper::fForceApplyFlxWgt

private

apply GFluxI::Weight() before returning event

Definition at line 184 of file GENIEHelper.h.

Referenced by GENIEHelper(), and Sample().

int evgb::GENIEHelper::fFunctionalBinning

private

Definition at line 169 of file GENIEHelper.h.

Referenced by GENIEHelper(), and InitializeFluxDriver().

std::string evgb::GENIEHelper::fFunctionalFlux

private

Definition at line 168 of file GENIEHelper.h.

Referenced by GENIEHelper(), and InitializeFluxDriver().

std::vector<int> evgb::GENIEHelper::fGenFlavors

private

pdg codes for flavors to generate

Definition at line 186 of file GENIEHelper.h.

Referenced by AtmoFluxCheck(), GENIEHelper(), HistogramFluxCheck(), InitializeFluxDriver(), Sample(), and Stop().

genie::EventRecord* evgb::GENIEHelper::fGenieEventRecord

private

last generated event

Definition at line 126 of file GENIEHelper.h.

Referenced by Sample(), and ~GENIEHelper().

std::string evgb::GENIEHelper::fGENIEMsgThresholds

private

additional XML file setting Messager level thresholds (":" separated)

Definition at line 200 of file GENIEHelper.h.

Referenced by GENIEHelper(), and StartGENIEMessenger().

std::string evgb::GENIEHelper::fGeoFile

private

name of file containing the Geometry description

Definition at line 124 of file GENIEHelper.h.

Referenced by SetMaxPathOutInfo().

TGeoManager* evgb::GENIEHelper::fGeoManager

private

pointer to ROOT TGeoManager

Definition at line 123 of file GENIEHelper.h.

Referenced by InitializeGeometry(), and Sample().

genie::GeomAnalyzerI* evgb::GENIEHelper::fGeomD

private

Definition at line 127 of file GENIEHelper.h.

Referenced by ConfigGeomScan(), Initialize(), InitializeFiducialSelection(), InitializeGeometry(), InitializeRockBoxSelection(), and ~GENIEHelper().

std::string evgb::GENIEHelper::fGeomScan

private

configuration for geometry scan to determine max pathlengths

Definition at line 206 of file GENIEHelper.h.

Referenced by ConfigGeomScan(), GENIEHelper(), and SetMaxPathOutInfo().

int evgb::GENIEHelper::fGHepPrintLevel

private

GHepRecord::SetPrintLevel(), -1=no-print.

Definition at line 201 of file GENIEHelper.h.

Referenced by GENIEHelper(), and Sample().

double evgb::GENIEHelper::fGlobalTimeOffset

private

overall time shift (ns) added to every particle time

Definition at line 180 of file GENIEHelper.h.

Referenced by GENIEHelper(), and Sample().

std::string evgb::GENIEHelper::fGMSGLAYOUT

private

format for GENIE log message [BASIC]|SIMPLE (SIMPLE=no timestamps)

Definition at line 199 of file GENIEHelper.h.

Referenced by GENIEHelper(), and SetGMSGLAYOUT().

std::string evgb::GENIEHelper::fGXMLPATH

private

locations for GENIE XML files

Definition at line 198 of file GENIEHelper.h.

Referenced by GENIEHelper(), ReadXSecTable(), and SetGXMLPATH().

TRandom3* evgb::GENIEHelper::fHelperRandom

private

random # generator for GENIEHelper

Definition at line 135 of file GENIEHelper.h.

Referenced by ExpandFluxFilePatternsDirect(), ExpandFluxFilePatternsIFDH(), GENIEHelper(), Initialize(), Sample(), Stop(), and ~GENIEHelper().

double evgb::GENIEHelper::fHistEventsPerSpill

private

number of events per spill for histogram fluxes - changes each spill

Definition at line 163 of file GENIEHelper.h.

Referenced by GENIEHelper(), Initialize(), and Stop().

ifdh_ns::ifdh* evgb::GENIEHelper::fIFDH

private

(optional) flux file handling

Definition at line 133 of file GENIEHelper.h.

Referenced by ExpandFluxFilePatternsIFDH(), and ~GENIEHelper().

int evgb::GENIEHelper::fMaxFluxFileMB

private

maximum size of flux files (MB)

Definition at line 145 of file GENIEHelper.h.

Referenced by ExpandFluxFilePatternsDirect(), and ExpandFluxFilePatternsIFDH().

int evgb::GENIEHelper::fMaxFluxFileNumber

private

maximum # of flux files

Definition at line 146 of file GENIEHelper.h.

Referenced by ExpandFluxFilePatternsDirect(), and ExpandFluxFilePatternsIFDH().

std::string evgb::GENIEHelper::fMaxPathOutInfo

private

output info if writing PathLengthList from GeomScan

Definition at line 207 of file GENIEHelper.h.

Referenced by SetMaxPathOutInfo(), and ~GENIEHelper().

double evgb::GENIEHelper::fMixerBaseline

private

baseline distance if genie flux can't calculate it

Definition at line 203 of file GENIEHelper.h.

Referenced by GENIEHelper(), and InitializeFluxDriver().

std::string evgb::GENIEHelper::fMixerConfig

private

configuration string for genie GFlavorMixerI

Definition at line 202 of file GENIEHelper.h.

Referenced by GENIEHelper(), and InitializeFluxDriver().

double evgb::GENIEHelper::fMonoEnergy

private

energy of monoenergetic neutrinos

Definition at line 167 of file GENIEHelper.h.

Referenced by GENIEHelper(), and InitializeFluxDriver().

double evgb::GENIEHelper::fPOTPerSpill

private

number of pot per spill

Definition at line 162 of file GENIEHelper.h.

Referenced by GENIEHelper(), Initialize(), and Stop().

double evgb::GENIEHelper::fRandomTimeOffset

private

additional random time shift (ns) added to every particle time

Definition at line 181 of file GENIEHelper.h.

Referenced by GENIEHelper(), and Sample().

std::vector<std::string> evgb::GENIEHelper::fSelectedFluxFiles

private

flux files selected after wildcard expansion and subset selection

Definition at line 144 of file GENIEHelper.h.

Referenced by AtmoFluxCheck(), ExpandFluxFilePatternsDirect(), ExpandFluxFilePatternsIFDH(), GENIEHelper(), HistogramFluxCheck(), InitializeFluxDriver(), SetMaxPathOutInfo(), and ~GENIEHelper().

int evgb::GENIEHelper::fSpillEvents

private

total events for this spill

Definition at line 164 of file GENIEHelper.h.

Referenced by GENIEHelper(), Initialize(), Sample(), and Stop().

double evgb::GENIEHelper::fSpillExposure

private

total exposure (i.e. pot) for this spill

Definition at line 165 of file GENIEHelper.h.

Referenced by GENIEHelper(), Initialize(), Sample(), and Stop().

std::string evgb::GENIEHelper::fSpillTimeConfig

private

alternative to flat spill distribution

Definition at line 182 of file GENIEHelper.h.

Referenced by GENIEHelper().

double evgb::GENIEHelper::fSurroundingMass

private

mass of material surrounding the detector that is intercepted by the cylinder for the histogram flux in kg

Definition at line 178 of file GENIEHelper.h.

Referenced by GENIEHelper(), and Initialize().

evgb::EvtTimeShiftI* evgb::GENIEHelper::fTimeShifter

private

generator for time offset within a spill

Definition at line 137 of file GENIEHelper.h.

Referenced by GENIEHelper(), and Sample().

std::string evgb::GENIEHelper::fTopVolume

private

top volume in the ROOT geometry in which to generate events

Definition at line 154 of file GENIEHelper.h.

Referenced by GENIEHelper(), InitializeGeometry(), InitializeRockBoxSelection(), Sample(), and SetMaxPathOutInfo().

double evgb::GENIEHelper::fTotalExposure

private

pot used from flux ntuple

Definition at line 166 of file GENIEHelper.h.

Referenced by GENIEHelper(), Initialize(), Sample(), Stop(), and ~GENIEHelper().

double evgb::GENIEHelper::fTotalHistFlux

private

total flux of neutrinos from flux histograms for used flavors

Definition at line 173 of file GENIEHelper.h.

Referenced by HistogramFluxCheck(), Initialize(), Stop(), and TotalHistFlux().

std::string evgb::GENIEHelper::fTuneName

private

GENIE R-3 Tune name (defines model configuration)

Definition at line 196 of file GENIEHelper.h.

Referenced by GENIEHelper(), Initialize(), and Sample().

bool evgb::GENIEHelper::fUseBlenderDist

private

get neutrino's travel distance from blender (default: true)

Definition at line 204 of file GENIEHelper.h.

Referenced by GENIEHelper(), and Sample().

bool evgb::GENIEHelper::fUseHelperRndGen4GENIE

private

use fHelperRandom for gRandom during Sample()

Definition at line 136 of file GENIEHelper.h.

Referenced by Sample().

std::string evgb::GENIEHelper::fWorldVolume

private

name of the world volume in the ROOT geometry

Definition at line 155 of file GENIEHelper.h.

Referenced by GENIEHelper(), InitializeGeometry(), InitializeRockBoxSelection(), Sample(), and SetMaxPathOutInfo().

double evgb::GENIEHelper::fXSecMassPOT

private

product of cross section, mass and POT/spill for histogram fluxes

Definition at line 172 of file GENIEHelper.h.

Referenced by Initialize(), and Stop().

std::string evgb::GENIEHelper::fXSecTable

private

cross section file (was $GSPLOAD)

Definition at line 195 of file GENIEHelper.h.

Referenced by GENIEHelper(), and ReadXSecTable().

---

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

• nugen/EventGeneratorBase/GENIE/GENIEHelper.h
• nugen/EventGeneratorBase/GENIE/GENIEHelper.cxx