Public Types | Public Member Functions | Private Member Functions | Private Attributes | List of all members
genie::flux::GNuMIFlux Class Reference

A GENIE flux driver encapsulating the NuMI neutrino flux. It reads-in the official GNUMI neutrino flux ntuples. Supports both geant3 and geant4 formats. More...

#include <GNuMIFlux.h>

Inheritance diagram for genie::flux::GNuMIFlux:
genie::GFluxI genie::flux::GFluxExposureI genie::flux::GFluxFileConfigI

Public Types

enum  EStdFluxWindow {
  kMinosNearDet, kMinosFarDet, kMinosNearRock, kMinosFarRock,
  kMinosNearCenter, kMinosFarCenter
}
 
typedef enum genie::flux::GNuMIFlux::EStdFluxWindow StdFluxWindow_t
 

Public Member Functions

 GNuMIFlux ()
 
 ~GNuMIFlux ()
 
const PDGCodeListFluxParticles (void)
 declare list of flux neutrinos that can be generated (for init. purposes) More...
 
double MaxEnergy (void)
 declare the max flux neutrino energy that can be generated (for init. purposes) More...
 
bool GenerateNext (void)
 generate the next flux neutrino (return false in err) More...
 
int PdgCode (void)
 returns the flux neutrino pdg code More...
 
double Weight (void)
 returns the flux neutrino weight (if any) More...
 
const TLorentzVector & Momentum (void)
 returns the flux neutrino 4-momentum More...
 
const TLorentzVector & Position (void)
 returns the flux neutrino 4-position (note: expect SI rather than physical units) More...
 
bool End (void)
 true if no more flux nu's can be thrown (eg reaching end of beam sim ntuples) More...
 
long int Index (void)
 returns corresponding index for current flux neutrino (e.g. for a flux ntuple returns the current entry number) More...
 
void Clear (Option_t *opt)
 reset state variables based on opt More...
 
void GenerateWeighted (bool gen_weighted)
 set whether to generate weighted or unweighted neutrinos More...
 
const GNuMIFluxPassThroughInfoPassThroughInfo (void)
 GNuMIFluxPassThroughInfo. More...
 
Long64_t GetEntryNumber ()
 index in chain More...
 
double GetDecayDist () const
 dist (user units) from dk to current pos More...
 
void MoveToZ0 (double z0)
 move ray origin to user coord Z0 More...
 
virtual double GetTotalExposure () const
 
virtual long int NFluxNeutrinos () const
 

of rays generated

More...
 
double POT_curr (void)
 current average POT (RWH?) More...
 
double UsedPOTs (void) const
 

of protons-on-target used

More...
 
double SumWeight (void) const
 integrated weight for flux neutrinos looped so far More...
 
void PrintCurrent (void)
 print current entry from leaves More...
 
void PrintConfig ()
 print the current configuration More...
 
std::vector< std::stringGetFileList ()
 list of files currently part of chain More...
 
virtual void LoadBeamSimData (const std::vector< std::string > &filenames, const std::string &det_loc)
 
virtual void GetBranchInfo (std::vector< std::string > &branchNames, std::vector< std::string > &branchClassNames, std::vector< void ** > &branchObjPointers)
 
virtual TTree * GetMetaDataTree ()
 
bool LoadConfig (string cfg)
 load a named configuration More...
 
void SetMaxEnergy (double Ev)
 specify maximum flx neutrino energy More...
 
void SetGenWeighted (bool genwgt=false)
 toggle whether GenerateNext() returns weight=1 flux (initial default false) More...
 
void SetEntryReuse (long int nuse=1)
 

of times to use entry before moving to next

More...
 
void SetTreeName (string name)
 set input tree name (default: "h10") More...
 
void ScanForMaxWeight (void)
 scan for max flux weight (before generating unweighted flux neutrinos) More...
 
void SetMaxWgtScan (double fudge=1.05, long int nentries=2500000)
 
void SetMaxEFudge (double fudge=1.05)
 
void SetApplyWindowTiltWeight (bool apply=true)
 
void SetLengthUnits (double user_units)
 Set units assumed by user. More...
 
double LengthUnits (void) const
 Return user units. More...
 
void SetBeamRotation (TRotation beamrot)
 < beam (0,0,0) relative to user frame, beam direction in user frame More...
 
void SetBeamCenter (TVector3 beam0)
 
TRotation GetBeamRotation () const
 rotation to apply from beam->user More...
 
TVector3 GetBeamCenter () const
 beam origin in user frame More...
 
bool SetFluxWindow (StdFluxWindow_t stdwindow, double padding=0)
 return false if unhandled More...
 
void SetFluxWindow (TVector3 p1, TVector3 p2, TVector3 p3)
 3 points define a plane (by default in user coordinates) More...
 
void GetFluxWindow (TVector3 &p1, TVector3 &p2, TVector3 &p3) const
 3 points define a plane in beam coordinate More...
 
void UseFluxAtNearDetCenter (void)
 force weights at MINOS detector "center" as found in ntuple More...
 
void UseFluxAtFarDetCenter (void)
 
void Beam2UserPos (const TLorentzVector &beamxyz, TLorentzVector &usrxyz) const
 
void Beam2UserDir (const TLorentzVector &beamdir, TLorentzVector &usrdir) const
 
void Beam2UserP4 (const TLorentzVector &beamp4, TLorentzVector &usrp4) const
 
void User2BeamPos (const TLorentzVector &usrxyz, TLorentzVector &beamxyz) const
 
void User2BeamDir (const TLorentzVector &usrdir, TLorentzVector &beamdir) const
 
void User2BeamP4 (const TLorentzVector &usrp4, TLorentzVector &beamp4) const
 
TVector3 FluxWindowNormal ()
 
- Public Member Functions inherited from genie::GFluxI
virtual ~GFluxI ()
 
- Public Member Functions inherited from genie::flux::GFluxExposureI
 GFluxExposureI (genie::flux::Exposure_t etype)
 
virtual ~GFluxExposureI ()
 
const char * GetExposureUnits () const
 what units are returned by GetTotalExposure? More...
 
genie::flux::Exposure_t GetExposureType () const
 
- Public Member Functions inherited from genie::flux::GFluxFileConfigI
 GFluxFileConfigI ()
 
virtual ~GFluxFileConfigI ()
 
virtual void LoadBeamSimData (const std::set< std::string > &filenames, const std::string &det_loc)
 
virtual void LoadBeamSimData (const std::string &filename, const std::string &det_loc)
 
virtual void SetXMLFileBase (std::string xmlbasename="")
 
virtual std::string GetXMLFileBase () const
 
virtual void SetFluxParticles (const PDGCodeList &particles)
 specify list of flux neutrino species More...
 
virtual void SetUpstreamZ (double z0)
 
virtual void SetNumOfCycles (long int ncycle)
 limit cycling through input files More...
 

Private Member Functions

bool GenerateNext_weighted (void)
 
void Initialize (void)
 
void SetDefaults (void)
 
void CleanUp (void)
 
void ResetCurrent (void)
 
void AddFile (TTree *tree, string fname)
 
void CalcEffPOTsPerNu (void)
 

Private Attributes

double fMaxEv
 maximum energy More...
 
bool fEnd
 end condition reached More...
 
std::vector< stringfNuFluxFilePatterns
 (potentially wildcarded) path(s) More...
 
string fNuFluxTreeName
 Tree name "h10" (g3) or "nudata" (g4) More...
 
TChain * fNuFluxTree
 TTree in g3numi or g4numi // REF ONLY! More...
 
string fNuFluxGen
 "g3numi" "g4numi" or "flugg" More...
 
g3numifG3NuMI
 g3numi ntuple More...
 
g4numifG4NuMI
 g4numi ntuple More...
 
fluggfFlugg
 flugg ntuple More...
 
int fNFiles
 number of files in chain More...
 
Long64_t fNEntries
 number of flux ntuple entries More...
 
Long64_t fIEntry
 current flux ntuple entry More...
 
Long64_t fNuTot
 cummulative # of entries (=fNEntries) More...
 
Long64_t fFilePOTs
 

of protons-on-target represented by all files

More...
 
double fWeight
 current neutrino weight, =1 if generating unweighted entries More...
 
double fMaxWeight
 max flux neutrino weight in input file More...
 
double fMaxWgtFudge
 fudge factor for estimating max wgt More...
 
long int fMaxWgtEntries
 

of entries in estimating max wgt

More...
 
double fMaxEFudge
 fudge factor for estmating max enu (0=> use fixed 120GeV) More...
 
long int fNUse
 how often to use same entry in a row More...
 
long int fIUse
 current # of times an entry has been used More...
 
double fSumWeight
 sum of weights for nus thrown so far More...
 
long int fNNeutrinos
 number of flux neutrinos thrown so far More...
 
double fEffPOTsPerNu
 what a entry is worth ... More...
 
double fAccumPOTs
 POTs used so far. More...
 
bool fGenWeighted
 does GenerateNext() give weights? More...
 
bool fApplyTiltWeight
 wgt due to window normal not || beam More...
 
bool fDetLocIsSet
 is a flux location (near/far) set? More...
 
int fUseFluxAtDetCenter
 use flux at near (-1) or far (+1) det center from ntuple? More...
 
double fLengthUnits
 units for coord in user exchanges More...
 
double fLengthScaleB2U
 scale factor beam (cm) –> user More...
 
double fLengthScaleU2B
 scale factor beam user –> (cm) More...
 
TLorentzVector fBeamZero
 beam origin in user coords More...
 
TLorentzRotation fBeamRot
 rotation applied beam –> user coord More...
 
TLorentzRotation fBeamRotInv
 
TVector3 fFluxWindowPtUser [3]
 user points of flux window More...
 
TLorentzVector fFluxWindowBase
 base point for flux window - beam coord More...
 
TLorentzVector fFluxWindowDir1
 extent for flux window (direction 1) More...
 
TLorentzVector fFluxWindowDir2
 extent for flux window (direction 2) More...
 
double fFluxWindowLen1
 
double fFluxWindowLen2
 
TVector3 fWindowNormal
 normal direction for flux window More...
 
TLorentzVector fgX4dkvtx
 decay 4-position beam coord More...
 
GNuMIFluxPassThroughInfofCurEntry
 copy of current ntuple entry info (owned structure) More...
 

Additional Inherited Members

- Static Public Member Functions inherited from genie::flux::GFluxExposureI
static const char * AsString (genie::flux::Exposure_t etype)
 
static genie::flux::Exposure_t StringToEnum (const char *chars, int maxChar=0)
 
- Protected Member Functions inherited from genie::GFluxI
 GFluxI ()
 
- Protected Attributes inherited from genie::flux::GFluxFileConfigI
PDGCodeListfPdgCList
 list of neutrino pdg-codes to generate More...
 
PDGCodeListfPdgCListRej
 list of nu pdg-codes seen but rejected More...
 
std::string fXMLbasename
 XML file that might hold config param_sets. More...
 
long int fNCycles
 

times to cycle through the ntuple(s)

More...
 
long int fICycle
 
double fZ0
 

Detailed Description

A GENIE flux driver encapsulating the NuMI neutrino flux. It reads-in the official GNUMI neutrino flux ntuples. Supports both geant3 and geant4 formats.

GNuMIFlux:

An implementation of the GFluxI interface that provides NuMI flux

http://www.hep.utexas.edu/~zarko/wwwgnumi/v19/

Author
Robert Hatcher <rhatcher fnal.gov> Fermi National Accelerator Laboratory

Jun 27, 2008

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

Definition at line 217 of file GNuMIFlux.h.

Member Typedef Documentation

Member Enumeration Documentation

Enumerator
kMinosNearDet 
kMinosFarDet 
kMinosNearRock 
kMinosFarRock 
kMinosNearCenter 
kMinosFarCenter 

Definition at line 330 of file GNuMIFlux.h.

330  {
331  kMinosNearDet, // appropriate for Near Detector
332  kMinosFarDet, // appropriate for Far Detector
333  kMinosNearRock, // appropriate for Near rock generation
334  kMinosFarRock, // appropriate for Far rock generation
335  kMinosNearCenter, // point location mimic near value in ntuple
336  kMinosFarCenter // point location mimic far value in ntuple
337  } StdFluxWindow_t;
enum genie::flux::GNuMIFlux::EStdFluxWindow StdFluxWindow_t

Constructor & Destructor Documentation

GNuMIFlux::GNuMIFlux ( )

Definition at line 113 of file GNuMIFlux.cxx.

113  :
115 {
116  this->Initialize();
117 }
GFluxExposureI(genie::flux::Exposure_t etype)
GNuMIFlux::~GNuMIFlux ( )

Definition at line 119 of file GNuMIFlux.cxx.

120 {
121  this->CleanUp();
122 }

Member Function Documentation

void GNuMIFlux::AddFile ( TTree *  tree,
string  fname 
)
private

Definition at line 1133 of file GNuMIFlux.cxx.

1134 {
1135  // Add a file to the chain
1136 
1137  ULong64_t nentries = thetree->GetEntries();
1138 
1139  // first/last "evtno" are the proton # of the first/last proton
1140  // that generated a neutrino ... not necessarily true first/last #
1141  // estimate we're probably not off by more than 100 ...
1142  // !Important change
1143  // Some files (due to some intermediate flugg issues) list evtno==0
1144  // when it isn't really true, we need to scan nearby values in case the
1145  // last entry is one of these (otherwise file contributes 0 POTs).
1146  // Also assume quantization of 500 (instead of 100).
1147 
1148  thetree->SetMakeClass(1); // need full ntuple decomposition for
1149  // the SetBranchAddress to work on g4numi ntuples. Works fine
1150  // without it on gnumi (geant3) and flugg ntuples [each with their
1151  // own slight differences] but shouldn't harm them either.
1152 
1153  Int_t evtno = 0;
1154  TBranch* br_evtno = 0;
1155  thetree->SetBranchAddress("evtno",&evtno, &br_evtno);
1156  Int_t evt_1 = 0x7FFFFFFF;
1157  Int_t evt_N = 1;
1158 #define OLDGUESS
1159 #ifdef OLDGUESS
1160  for (int j=0; j<50; ++j) {
1161  thetree->GetEntry(j);
1162  if (evtno != 0) evt_1 = TMath::Min(evtno,evt_1);
1163  thetree->GetEntry(nentries-1 -j );
1164  if (evtno != 0) evt_N = TMath::Max(evtno,evt_N);
1165  }
1166  // looks like low counts are due to "evtno" not including
1167  // protons that miss the actual target (hit baffle, etc)
1168  // this will vary with beam conditions parameters
1169  // so we should round way up, those generating flugg files
1170  // aren't going to quantize less than 1000
1171  // though 500 would probably be okay, 100 would not.
1172  const Int_t nquant = 1000; // 500; // 100
1173  const Double_t rquant = nquant;
1174 #else
1175  for (int j=0; j<50; ++j) {
1176  thetree->GetEntry(j);
1177  if (evtno != 0) evt_1 = TMath::Min(evtno,evt_1);
1178  std::cout << "[" << j << "] evtno=" << evtno << " evt_1=" << evt_1 << std::endl;
1179  }
1180  for (int j=0; j<50; ++j) {
1181  thetree->GetEntry(nentries-1 -j );
1182  if (evtno != 0) evt_N = TMath::Max(evtno,evt_N);
1183  std::cout << "[" << (nentries-1-j) << "] evtno=" << evtno << " evt_N=" << evt_N << std::endl;
1184  }
1185 
1186  Int_t nquant = 1000; // 500;
1187  Double_t rquant = nquant;
1188 #endif
1189 
1190  Int_t est_1 = (TMath::FloorNint(evt_1/rquant))*nquant + 1;
1191  Int_t est_N = (TMath::FloorNint((evt_N-1)/rquant)+1)*nquant;
1192  ULong64_t npots = est_N - est_1 + 1;
1193  LOG("Flux",pNOTICE) //INFO)
1194  << fNuFluxTreeName << "->AddFile() of " << nentries << " entries ["
1195  << evt_1 << ":" << evt_N << "%" << nquant
1196  << "(" << est_1 << ":" << est_N << ")="
1197  << npots <<" POTs] in {" << fNuFluxGen << "} file: " << fname;
1198  fNuTot += nentries;
1199  fFilePOTs += npots;
1200  fNFiles++;
1201 
1202  fNuFluxTree->AddFile(fname.c_str());
1203 }
Long64_t fFilePOTs
of protons-on-target represented by all files
Definition: GNuMIFlux.h:398
TChain * fNuFluxTree
TTree in g3numi or g4numi // REF ONLY!
Definition: GNuMIFlux.h:389
string fNuFluxTreeName
Tree name "h10" (g3) or "nudata" (g4)
Definition: GNuMIFlux.h:388
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
int fNFiles
number of files in chain
Definition: GNuMIFlux.h:394
#define pNOTICE
Definition: Messenger.h:61
string fNuFluxGen
"g3numi" "g4numi" or "flugg"
Definition: GNuMIFlux.h:390
Long64_t fNuTot
cummulative # of entries (=fNEntries)
Definition: GNuMIFlux.h:397
QTextStream & endl(QTextStream &s)
void GNuMIFlux::Beam2UserDir ( const TLorentzVector &  beamdir,
TLorentzVector &  usrdir 
) const

Definition at line 976 of file GNuMIFlux.cxx.

978 {
979  usrdir = fLengthScaleB2U*(fBeamRot*beamdir);
980 }
double fLengthScaleB2U
scale factor beam (cm) –> user
Definition: GNuMIFlux.h:419
TLorentzRotation fBeamRot
rotation applied beam –> user coord
Definition: GNuMIFlux.h:423
void GNuMIFlux::Beam2UserP4 ( const TLorentzVector &  beamp4,
TLorentzVector &  usrp4 
) const

Definition at line 981 of file GNuMIFlux.cxx.

983 {
984  usrp4 = fBeamRot*beamp4;
985 }
TLorentzRotation fBeamRot
rotation applied beam –> user coord
Definition: GNuMIFlux.h:423
void GNuMIFlux::Beam2UserPos ( const TLorentzVector &  beamxyz,
TLorentzVector &  usrxyz 
) const

Definition at line 971 of file GNuMIFlux.cxx.

973 {
974  usrxyz = fLengthScaleB2U*(fBeamRot*beamxyz) + fBeamZero;
975 }
TLorentzVector fBeamZero
beam origin in user coords
Definition: GNuMIFlux.h:422
double fLengthScaleB2U
scale factor beam (cm) –> user
Definition: GNuMIFlux.h:419
TLorentzRotation fBeamRot
rotation applied beam –> user coord
Definition: GNuMIFlux.h:423
void GNuMIFlux::CalcEffPOTsPerNu ( void  )
private

Definition at line 439 of file GNuMIFlux.cxx.

440 {
441  // do this if flux window changes or # of files changes
442 
443  if (!fNuFluxTree) return; // not yet fully configured
444 
445  // effpots = mc_pots * (wgtfunction-area) / window-area / wgt-max-est
446  // wgtfunction-area = pi * radius-det-element^2 = pi * (100.cm)^2
447 
448  // this should match what is used in the CalcEnuWgt()
449  const double kRDET = 100.0; // set to flux per 100 cm radius
450  const double kRDET2 = kRDET * kRDET;
451  double flux_area = fFluxWindowDir1.Vect().Cross(fFluxWindowDir2.Vect()).Mag();
452  LOG("Flux",pNOTICE) << "in CalcEffPOTsPerNu, area = " << flux_area;
453 
454  if ( flux_area < 1.0e-30 ) {
455  LOG("Flux", pWARN)
456  << "CalcEffPOTsPerNu called with flux window area effectively zero";
457  flux_area = 1;
458  }
459  double area_ratio = TMath::Pi() * kRDET2 / flux_area;
460  fEffPOTsPerNu = area_ratio * ( (double)fFilePOTs / (double)fNEntries );
461 }
Long64_t fFilePOTs
of protons-on-target represented by all files
Definition: GNuMIFlux.h:398
TLorentzVector fFluxWindowDir2
extent for flux window (direction 2)
Definition: GNuMIFlux.h:429
Long64_t fNEntries
number of flux ntuple entries
Definition: GNuMIFlux.h:395
TChain * fNuFluxTree
TTree in g3numi or g4numi // REF ONLY!
Definition: GNuMIFlux.h:389
const double e
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
TLorentzVector fFluxWindowDir1
extent for flux window (direction 1)
Definition: GNuMIFlux.h:428
#define pWARN
Definition: Messenger.h:60
double fEffPOTsPerNu
what a entry is worth ...
Definition: GNuMIFlux.h:410
#define pNOTICE
Definition: Messenger.h:61
void GNuMIFlux::CleanUp ( void  )
private

Definition at line 1115 of file GNuMIFlux.cxx.

1116 {
1117  LOG("Flux", pNOTICE) << "Cleaning up...";
1118 
1119  if (fPdgCList) delete fPdgCList;
1120  if (fPdgCListRej) delete fPdgCListRej;
1121  if (fCurEntry) delete fCurEntry;
1122 
1123  if ( fG3NuMI ) delete fG3NuMI;
1124  if ( fG4NuMI ) delete fG4NuMI;
1125  if ( fFlugg ) delete fFlugg;
1126 
1127  LOG("Flux", pNOTICE)
1128  << " flux file cycles: " << fICycle << " of " << fNCycles
1129  << ", entry " << fIEntry << " use: " << fIUse << " of " << fNUse;
1130 }
long int fIUse
current # of times an entry has been used
Definition: GNuMIFlux.h:407
Long64_t fIEntry
current flux ntuple entry
Definition: GNuMIFlux.h:396
g4numi * fG4NuMI
g4numi ntuple
Definition: GNuMIFlux.h:392
PDGCodeList * fPdgCList
list of neutrino pdg-codes to generate
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
long int fNCycles
times to cycle through the ntuple(s)
flugg * fFlugg
flugg ntuple
Definition: GNuMIFlux.h:393
GNuMIFluxPassThroughInfo * fCurEntry
copy of current ntuple entry info (owned structure)
Definition: GNuMIFlux.h:436
PDGCodeList * fPdgCListRej
list of nu pdg-codes seen but rejected
g3numi * fG3NuMI
g3numi ntuple
Definition: GNuMIFlux.h:391
#define pNOTICE
Definition: Messenger.h:61
long int fNUse
how often to use same entry in a row
Definition: GNuMIFlux.h:406
void GNuMIFlux::Clear ( Option_t *  opt)
virtual

reset state variables based on opt

Implements genie::GFluxI.

Definition at line 1009 of file GNuMIFlux.cxx.

1010 {
1011  // Clear the driver state
1012  //
1013  LOG("Flux", pWARN) << "GSimpleNtpFlux::Clear(" << opt << ") called";
1014  // do it in all cases, but EVGDriver/GMCJDriver will pass "CycleHistory"
1015 
1016  fICycle = 0;
1017 
1018  fSumWeight = 0;
1019  fNNeutrinos = 0;
1020  fAccumPOTs = 0;
1021 }
double fSumWeight
sum of weights for nus thrown so far
Definition: GNuMIFlux.h:408
opt
Definition: train.py:196
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
double fAccumPOTs
POTs used so far.
Definition: GNuMIFlux.h:411
#define pWARN
Definition: Messenger.h:60
long int fNNeutrinos
number of flux neutrinos thrown so far
Definition: GNuMIFlux.h:409
bool genie::flux::GNuMIFlux::End ( void  )
inlinevirtual

true if no more flux nu's can be thrown (eg reaching end of beam sim ntuples)

Implements genie::GFluxI.

Definition at line 237 of file GNuMIFlux.h.

237 { return fEnd; }
bool fEnd
end condition reached
Definition: GNuMIFlux.h:385
const PDGCodeList& genie::flux::GNuMIFlux::FluxParticles ( void  )
inlinevirtual

declare list of flux neutrinos that can be generated (for init. purposes)

Implements genie::GFluxI.

Definition at line 230 of file GNuMIFlux.h.

230 { return *fPdgCList; }
PDGCodeList * fPdgCList
list of neutrino pdg-codes to generate
TVector3 genie::flux::GNuMIFlux::FluxWindowNormal ( )
inline

Definition at line 368 of file GNuMIFlux.h.

368 { return fWindowNormal; }
TVector3 fWindowNormal
normal direction for flux window
Definition: GNuMIFlux.h:432
bool GNuMIFlux::GenerateNext ( void  )
virtual

generate the next flux neutrino (return false in err)

Implements genie::GFluxI.

Definition at line 139 of file GNuMIFlux.cxx.

140 {
141 // Get next (unweighted) flux ntuple entry on the specified detector location
142 //
144  while ( true ) {
145  // Check for end of flux ntuple
146  bool end = this->End();
147  if ( end ) {
148  LOG("Flux", pNOTICE) << "GenerateNext signaled End() ";
149  return false;
150  }
151 
152  // Get next weighted flux ntuple entry
153  bool nextok = this->GenerateNext_weighted();
154  if ( fGenWeighted ) return nextok;
155  if ( ! nextok ) continue;
156 
157  /* RWH - debug purposes
158  if ( fNCycles == 0 ) {
159  LOG("Flux", pNOTICE)
160  << "Got flux entry: " << fIEntry
161  << " - Cycle: "<< fICycle << "/ infinite";
162  } else {
163  LOG("Flux", pNOTICE)
164  << "Got flux entry: "<< fIEntry
165  << " - Cycle: "<< fICycle << "/"<< fNCycles;
166  }
167  */
168 
169  // Get fractional weight & decide whether to accept curr flux neutrino
170  double f = this->Weight() / fMaxWeight;
171  //LOG("Flux", pNOTICE)
172  // << "Curr flux neutrino fractional weight = " << f;
173  if (f > 1.) {
174  fMaxWeight = this->Weight() * fMaxWgtFudge; // bump the weight
175  LOG("Flux", pERROR)
176  << "** Fractional weight = " << f
177  << " > 1 !! Bump fMaxWeight estimate to " << fMaxWeight
178  << PassThroughInfo();
179  }
180  double r = (f < 1.) ? rnd->RndFlux().Rndm() : 0;
181  bool accept = ( r < f );
182  if ( accept ) {
183 
184 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
185  LOG("Flux", pNOTICE)
186  << "Generated beam neutrino: "
187  << "\n pdg-code: " << fCurEntry->fgPdgC
188  << "\n p4: " << utils::print::P4AsShortString(&(fCurEntry->fgP4))
189  << "\n x4: " << utils::print::X4AsString(&(fCurEntry->fgX4))
191  << "\n x4: " << utils::print::X4AsString(&(fCurEntry->fgX4User));
192 #endif
193 
194  fWeight = 1.;
195  return true;
196  }
197 
198  //LOG("Flux", pNOTICE)
199  // << "** Rejecting current flux neutrino based on the flux weight only";
200  }
201  return false;
202 }
end
while True: pbar.update(maxval-len(onlies[E][S])) #print iS, "/", len(onlies[E][S]) found = False for...
TLorentzVector fgP4
generated nu 4-momentum beam coord
Definition: GNuMIFlux.h:102
string P4AsShortString(const TLorentzVector *p)
Definition: PrintUtils.cxx:45
#define pERROR
Definition: Messenger.h:59
static RandomGen * Instance()
Access instance.
Definition: RandomGen.cxx:71
bool fGenWeighted
does GenerateNext() give weights?
Definition: GNuMIFlux.h:413
TLorentzVector fgX4User
generated nu 4-position user coord
Definition: GNuMIFlux.h:105
A singleton holding random number generator classes. All random number generation in GENIE should tak...
Definition: RandomGen.h:29
double fMaxWgtFudge
fudge factor for estimating max wgt
Definition: GNuMIFlux.h:402
bool End(void)
true if no more flux nu&#39;s can be thrown (eg reaching end of beam sim ntuples)
Definition: GNuMIFlux.h:237
double Weight(void)
returns the flux neutrino weight (if any)
Definition: GNuMIFlux.h:234
TLorentzVector fgX4
generated nu 4-position beam coord
Definition: GNuMIFlux.h:103
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
const GNuMIFluxPassThroughInfo & PassThroughInfo(void)
GNuMIFluxPassThroughInfo.
Definition: GNuMIFlux.h:252
bool GenerateNext_weighted(void)
Definition: GNuMIFlux.cxx:204
GNuMIFluxPassThroughInfo * fCurEntry
copy of current ntuple entry info (owned structure)
Definition: GNuMIFlux.h:436
double fMaxWeight
max flux neutrino weight in input file
Definition: GNuMIFlux.h:401
int fgPdgC
generated nu pdg-code
Definition: GNuMIFlux.h:98
#define pNOTICE
Definition: Messenger.h:61
TLorentzVector fgP4User
generated nu 4-momentum user coord
Definition: GNuMIFlux.h:104
TRandom3 & RndFlux(void) const
rnd number generator used by flux drivers
Definition: RandomGen.h:71
string X4AsString(const TLorentzVector *x)
Definition: PrintUtils.cxx:57
double fWeight
current neutrino weight, =1 if generating unweighted entries
Definition: GNuMIFlux.h:400
bool GNuMIFlux::GenerateNext_weighted ( void  )
private

user might modify list via SetFluxParticles() in order to reject certain flavors, even if they're found in the file. So don't make a big fuss. Spit out a single message and then stop reporting that flavor as problematic.

Definition at line 204 of file GNuMIFlux.cxx.

205 {
206 // Get next (weighted) flux ntuple entry on the specified detector location
207 //
208 
209  // Check whether a flux ntuple has been loaded
210  if ( ! fG3NuMI && ! fG4NuMI && ! fFlugg ) {
211  LOG("Flux", pFATAL)
212  << "The flux driver has not been properly configured";
213  //return false; // don't do this - creates an infinite loop!
214  exit(1);
215  }
216 
217  // Reuse an entry?
218  //std::cout << " ***** iuse " << fIUse << " nuse " << fNUse
219  // << " ientry " << fIEntry << " nentry " << fNEntries
220  // << " icycle " << fICycle << " ncycle " << fNCycles << std::endl;
221  if ( fIUse < fNUse && fIEntry >= 0 ) {
222  // Reuse this entry
223  fIUse++;
224  } else {
225  // Reset previously generated neutrino code / 4-p / 4-x
226  this->ResetCurrent();
227  // Move on, read next flux ntuple entry
228  fIEntry++;
229  if ( fIEntry >= fNEntries ) {
230  // Ran out of entries @ the current cycle of this flux file
231  // Check whether more (or infinite) number of cycles is requested
232  if ( fICycle < fNCycles || fNCycles == 0 ) {
233  fICycle++;
234  fIEntry=0;
235  } else {
236  LOG("Flux", pWARN)
237  << "No more entries in input flux neutrino ntuple, cycle "
238  << fICycle << " of " << fNCycles;
239  fEnd = true;
240  // assert(0);
241  return false;
242  }
243  }
244 
245  if ( fG3NuMI ) {
248  } else if ( fG4NuMI ) {
251  } else if ( fFlugg ) {
254  } else {
255  LOG("Flux", pERROR) << "No ntuple configured";
256  fEnd = true;
257  //assert(0);
258  return false;
259  }
260 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
261  LOG("Flux",pDEBUG)
262  << "got " << fNNeutrinos << " new fIEntry " << fIEntry
263  << " evtno " << fCurEntry->evtno;
264 #endif
265 
266  fIUse = 1;
267  fCurEntry->pcodes = 0; // fetched entry has geant codes
268  fCurEntry->units = 0; // fetched entry has original units
269 
270  // Convert the current gnumi neutrino flavor mode into a neutrino pdg code
271  // Also convert other particle codes in GNuMIFluxPassThroughInfo to PDG
273  // here we might want to do flavor oscillations or simple mappings
275  }
276 
277  // Check neutrino pdg against declared list of neutrino species declared
278  // by the current instance of the NuMI neutrino flux driver.
279  // No undeclared neutrino species will be accepted at this point as GENIE
280  // has already been configured to handle the specified list.
281  // Make sure that the appropriate list of flux neutrino species was set at
282  // initialization via GNuMIFlux::SetFluxParticles(const PDGCodeList &)
283 
284  // update the # POTs, number of neutrinos
285  // do this HERE (before rejecting flavors that users might be weeding out)
286  // in order to keep the POT accounting correct. This allows one to get
287  // the right normalization for generating only events from the intrinsic
288  // nu_e entries.
290  fNNeutrinos++;
291 
293  /// user might modify list via SetFluxParticles() in order to reject certain
294  /// flavors, even if they're found in the file. So don't make a big fuss.
295  /// Spit out a single message and then stop reporting that flavor as problematic.
296  int badpdg = fCurEntry->fgPdgC;
297  if ( ! fPdgCListRej->ExistsInPDGCodeList(badpdg) ) {
298  fPdgCListRej->push_back(badpdg);
299  LOG("Flux", pWARN)
300  << "Encountered neutrino specie (" << badpdg
301  << " pcodes=" << fCurEntry->pcodes << ")"
302  << " that wasn't in SetFluxParticles() list, "
303  << "\nDeclared list of neutrino species: " << *fPdgCList;
304  }
305  return false;
306  }
307 
308  // Update the curr neutrino weight and energy
309 
310  // Check current neutrino energy against the maximum flux neutrino energy
311  // declared by the current instance of the NuMI neutrino flux driver.
312  // No flux neutrino exceeding that maximum energy will be accepted at this
313  // point as that maximum energy has already been used for normalizing the
314  // interaction probabilities.
315  // Make sure that the appropriate maximum flux neutrino energy was set at
316  // initialization via GNuMIFlux::SetMaxEnergy(double Ev)
317 
319 
320  double Ev = 0;
321  double& wgt_xy = fCurEntry->fgXYWgt;
322  switch ( fUseFluxAtDetCenter ) {
323  case -1: // near detector
324  wgt_xy = fCurEntry->nwtnear;
325  Ev = fCurEntry->nenergyn;
326  break;
327  case +1: // far detector
328  wgt_xy = fCurEntry->nwtfar;
329  Ev = fCurEntry->nenergyf;
330  break;
331  default: // recalculate on x-y window
332  RandomGen * rnd = RandomGen::Instance();
333  fCurEntry->fgX4 += ( rnd->RndFlux().Rndm()*fFluxWindowDir1 +
334  rnd->RndFlux().Rndm()*fFluxWindowDir2 );
335  fCurEntry->CalcEnuWgt(fCurEntry->fgX4,Ev,wgt_xy);
336  break;
337  }
338 
339  if (Ev > fMaxEv) {
340  LOG("Flux", pWARN)
341  << "Flux neutrino energy exceeds declared maximum neutrino energy"
342  << "\nEv = " << Ev << "(> Ev{max} = " << fMaxEv << ")";
343  }
344 
345  // Set the current flux neutrino 4-momentum
346  // this is in *beam* coordinates
347  fgX4dkvtx = TLorentzVector( fCurEntry->vx,
348  fCurEntry->vy,
349  fCurEntry->vz, 0.);
350  // don't use TLorentzVector here for Mag() due to - on metric
351  // normalize direction to 1.0
352  TVector3 dirNu = (fCurEntry->fgX4.Vect() - fgX4dkvtx.Vect()).Unit();
353  fCurEntry->fgP4.SetPxPyPzE( Ev*dirNu.X(),
354  Ev*dirNu.Y(),
355  Ev*dirNu.Z(), Ev);
356 
357  // calculate the weight, potentially includes effect from tilted window
358  // must be done *after* neutrino direction is determined
359  fWeight = fCurEntry->nimpwt * fCurEntry->fgXYWgt; // full weight
360  if ( fApplyTiltWeight ) {
361  double tiltwgt = dirNu.Dot( FluxWindowNormal() );
362  fWeight *= TMath::Abs( tiltwgt );
363  }
364 
365  // update sume of weights
366  fSumWeight += this->Weight();
367 
368  // Set the current flux neutrino 4-position, direction in user coord
371 
372  // if desired, move to user specified user coord z
373  if ( TMath::Abs(fZ0) < 1.0e30 ) this->MoveToZ0(fZ0);
374 
375 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
376  LOG("Flux", pINFO)
377  << "Generated neutrino: " << fIEntry << " " << fCurEntry->evtno
378  << " nenergyn " << fCurEntry->nenergyn
379  << "\n pdg-code: " << fCurEntry->fgPdgC
380  << "\n p4 beam: " << utils::print::P4AsShortString(&fCurEntry->fgP4)
381  << "\n x4 beam: " << utils::print::X4AsString(&fCurEntry->fgX4)
382  << "\n p4 user: " << utils::print::P4AsShortString(&(fCurEntry->fgP4User))
383  << "\n x4 user: " << utils::print::X4AsString(&(fCurEntry->fgX4User));
384 #endif
385  if ( Ev > fMaxEv ) {
386  LOG("Flux", pFATAL)
387  << "Generated neutrino had E_nu = " << Ev << " > " << fMaxEv
388  << " maximum ";
389  assert(0);
390  }
391 
392  return true;
393 }
double fSumWeight
sum of weights for nus thrown so far
Definition: GNuMIFlux.h:408
TLorentzVector fgP4
generated nu 4-momentum beam coord
Definition: GNuMIFlux.h:102
string P4AsShortString(const TLorentzVector *p)
Definition: PrintUtils.cxx:45
long int fIUse
current # of times an entry has been used
Definition: GNuMIFlux.h:407
#define pERROR
Definition: Messenger.h:59
static RandomGen * Instance()
Access instance.
Definition: RandomGen.cxx:71
Long64_t fIEntry
current flux ntuple entry
Definition: GNuMIFlux.h:396
virtual Int_t GetEntry(Long64_t entry)
#define pFATAL
Definition: Messenger.h:56
TLorentzVector fgX4User
generated nu 4-position user coord
Definition: GNuMIFlux.h:105
bool fApplyTiltWeight
wgt due to window normal not || beam
Definition: GNuMIFlux.h:414
TLorentzVector fFluxWindowDir2
extent for flux window (direction 2)
Definition: GNuMIFlux.h:429
bool ExistsInPDGCodeList(int pdg_code) const
Long64_t fNEntries
number of flux ntuple entries
Definition: GNuMIFlux.h:395
void Beam2UserP4(const TLorentzVector &beamp4, TLorentzVector &usrp4) const
Definition: GNuMIFlux.cxx:981
bool fEnd
end condition reached
Definition: GNuMIFlux.h:385
A singleton holding random number generator classes. All random number generation in GENIE should tak...
Definition: RandomGen.h:29
void MakeCopy(const g3numi *)
pull in from g3 ntuple
Definition: GNuMIFlux.cxx:1476
double Weight(void)
returns the flux neutrino weight (if any)
Definition: GNuMIFlux.h:234
TLorentzVector fgX4
generated nu 4-position beam coord
Definition: GNuMIFlux.h:103
void MoveToZ0(double z0)
move ray origin to user coord Z0
Definition: GNuMIFlux.cxx:403
g4numi * fG4NuMI
g4numi ntuple
Definition: GNuMIFlux.h:392
TLorentzVector fgX4dkvtx
decay 4-position beam coord
Definition: GNuMIFlux.h:434
PDGCodeList * fPdgCList
list of neutrino pdg-codes to generate
virtual Int_t GetEntry(Long64_t entry)
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
TLorentzVector fFluxWindowDir1
extent for flux window (direction 1)
Definition: GNuMIFlux.h:428
int CalcEnuWgt(const TLorentzVector &xyz, double &enu, double &wgt_xy) const
Definition: GNuMIFlux.cxx:1729
#define pINFO
Definition: Messenger.h:62
long int fNCycles
times to cycle through the ntuple(s)
flugg * fFlugg
flugg ntuple
Definition: GNuMIFlux.h:393
GNuMIFluxPassThroughInfo * fCurEntry
copy of current ntuple entry info (owned structure)
Definition: GNuMIFlux.h:436
double fAccumPOTs
POTs used so far.
Definition: GNuMIFlux.h:411
#define pWARN
Definition: Messenger.h:60
PDGCodeList * fPdgCListRej
list of nu pdg-codes seen but rejected
void Beam2UserPos(const TLorentzVector &beamxyz, TLorentzVector &usrxyz) const
Definition: GNuMIFlux.cxx:971
double fMaxWeight
max flux neutrino weight in input file
Definition: GNuMIFlux.h:401
virtual Int_t GetEntry(Long64_t entry)
int fgPdgC
generated nu pdg-code
Definition: GNuMIFlux.h:98
g3numi * fG3NuMI
g3numi ntuple
Definition: GNuMIFlux.h:391
int fUseFluxAtDetCenter
use flux at near (-1) or far (+1) det center from ntuple?
Definition: GNuMIFlux.h:416
double fEffPOTsPerNu
what a entry is worth ...
Definition: GNuMIFlux.h:410
TVector3 FluxWindowNormal()
Definition: GNuMIFlux.h:368
TLorentzVector fgP4User
generated nu 4-momentum user coord
Definition: GNuMIFlux.h:104
TRandom3 & RndFlux(void) const
rnd number generator used by flux drivers
Definition: RandomGen.h:71
string X4AsString(const TLorentzVector *x)
Definition: PrintUtils.cxx:57
double fWeight
current neutrino weight, =1 if generating unweighted entries
Definition: GNuMIFlux.h:400
double fMaxEv
maximum energy
Definition: GNuMIFlux.h:384
TLorentzVector fFluxWindowBase
base point for flux window - beam coord
Definition: GNuMIFlux.h:427
long int fNNeutrinos
number of flux neutrinos thrown so far
Definition: GNuMIFlux.h:409
void push_back(int pdg_code)
Definition: PDGCodeList.cxx:58
#define pDEBUG
Definition: Messenger.h:63
void GNuMIFlux::GenerateWeighted ( bool  gen_weighted)
virtual

set whether to generate weighted or unweighted neutrinos

Implements genie::GFluxI.

Definition at line 1023 of file GNuMIFlux.cxx.

1024 {
1025  // Set whether to generate weighted rays
1026  //
1027  fGenWeighted = gen_weighted;
1028 }
bool fGenWeighted
does GenerateNext() give weights?
Definition: GNuMIFlux.h:413
TVector3 GNuMIFlux::GetBeamCenter ( ) const

beam origin in user frame

Definition at line 941 of file GNuMIFlux.cxx.

942 {
943  TVector3 beam0 = fBeamZero.Vect();
944  return beam0;
945 }
TLorentzVector fBeamZero
beam origin in user coords
Definition: GNuMIFlux.h:422
TRotation GNuMIFlux::GetBeamRotation ( ) const

rotation to apply from beam->user

Definition at line 928 of file GNuMIFlux.cxx.

929 {
930  // rotation is really only 3-d vector, but we'll be operating on LorentzV's
931  // give people back the original TRotation ... not pretty
932  // ... it think this is right
933  TRotation rot3;
934  const TLorentzRotation& rot4 = fBeamRot;
935  TVector3 newX(rot4.XX(),rot4.XY(),rot4.XZ());
936  TVector3 newY(rot4.YX(),rot4.YY(),rot4.YZ());
937  TVector3 newZ(rot4.ZX(),rot4.ZY(),rot4.ZZ());
938  rot3.RotateAxes(newX,newY,newZ);
939  return rot3.Inverse();
940 }
TLorentzRotation fBeamRot
rotation applied beam –> user coord
Definition: GNuMIFlux.h:423
void GNuMIFlux::GetBranchInfo ( std::vector< std::string > &  branchNames,
std::vector< std::string > &  branchClassNames,
std::vector< void ** > &  branchObjPointers 
)
virtual

allow caller to copy current status / ntuple entry info in the output file by providing copies of internal info

Assumes that branch object pointers will not change which may require either a copy be made or, if using the class directly for reading the branch, one must force ROOT to not autodelete: myns::MyClassType* fCurrMyClass = new myns::MyClassType; myTree->SetBranchAddress("bname",&fCurMyClass); //? TBranch* b = myTree->GetBranch("bname"); //? b->SetAutoDelete(false);

ensure vectors are sized sufficiently (or use .push_back()) branchNames[i] = "bname" branchClassNames[i] = "myns::MyClassType" branchObjPointers[i] = (void**)

Reimplemented from genie::flux::GFluxFileConfigI.

Definition at line 663 of file GNuMIFlux.cxx.

666 {
667  // allow flux driver to report back current status and/or ntuple entry
668  // info for possible recording in the output file by supplying
669  // the class name, and a pointer to the object that will be filled
670  // as well as a suggested name for the branch.
671 
672  branchNames.push_back("flux");
673  branchClassNames.push_back("genie::flux::GNuMIFluxPassThroughInfo");
674  branchObjPointers.push_back((void**)&fCurEntry);
675 
676 }
GNuMIFluxPassThroughInfo * fCurEntry
copy of current ntuple entry info (owned structure)
Definition: GNuMIFlux.h:436
double GNuMIFlux::GetDecayDist ( ) const

dist (user units) from dk to current pos

Definition at line 395 of file GNuMIFlux.cxx.

396 {
397  // return distance (user units) between dk point and start position
398  // these are in beam units
399  TVector3 x3diff = fCurEntry->fgX4.Vect() - fgX4dkvtx.Vect();
400  return x3diff.Mag() * fLengthScaleB2U;
401 }
TLorentzVector fgX4
generated nu 4-position beam coord
Definition: GNuMIFlux.h:103
TLorentzVector fgX4dkvtx
decay 4-position beam coord
Definition: GNuMIFlux.h:434
GNuMIFluxPassThroughInfo * fCurEntry
copy of current ntuple entry info (owned structure)
Definition: GNuMIFlux.h:436
double fLengthScaleB2U
scale factor beam (cm) –> user
Definition: GNuMIFlux.h:419
Long64_t genie::flux::GNuMIFlux::GetEntryNumber ( )
inline

index in chain

Definition at line 253 of file GNuMIFlux.h.

std::vector< std::string > GNuMIFlux::GetFileList ( )

list of files currently part of chain

Definition at line 2436 of file GNuMIFlux.cxx.

2437 {
2438  std::vector<std::string> flist;
2439  TObjArray *fileElements=fNuFluxTree->GetListOfFiles();
2440  TIter next(fileElements);
2441  TChainElement *chEl=0;
2442  while (( chEl=(TChainElement*)next() )) {
2443  flist.push_back(chEl->GetTitle());
2444  }
2445  return flist;
2446 }
TChain * fNuFluxTree
TTree in g3numi or g4numi // REF ONLY!
Definition: GNuMIFlux.h:389
void GNuMIFlux::GetFluxWindow ( TVector3 &  p1,
TVector3 &  p2,
TVector3 &  p3 
) const

3 points define a plane in beam coordinate

Definition at line 904 of file GNuMIFlux.cxx.

905 {
906  // return flux window points
907  p0 = fFluxWindowPtUser[0];
908  p1 = fFluxWindowPtUser[1];
909  p2 = fFluxWindowPtUser[2];
910 
911 }
TVector3 fFluxWindowPtUser[3]
user points of flux window
Definition: GNuMIFlux.h:426
TTree * GNuMIFlux::GetMetaDataTree ( )
virtual

Reimplemented from genie::flux::GFluxFileConfigI.

Definition at line 677 of file GNuMIFlux.cxx.

677 { return 0; } // there is none
double GNuMIFlux::GetTotalExposure ( ) const
virtual

Implements genie::flux::GFluxExposureI.

Definition at line 125 of file GNuMIFlux.cxx.

126 {
127  // complete the GFluxExposureI interface
128  return UsedPOTs();
129 }
double UsedPOTs(void) const
of protons-on-target used
Definition: GNuMIFlux.cxx:464
long int genie::flux::GNuMIFlux::Index ( void  )
inlinevirtual

returns corresponding index for current flux neutrino (e.g. for a flux ntuple returns the current entry number)

Implements genie::GFluxI.

Definition at line 238 of file GNuMIFlux.h.

238 { return fIEntry; }
Long64_t fIEntry
current flux ntuple entry
Definition: GNuMIFlux.h:396
void GNuMIFlux::Initialize ( void  )
private

Definition at line 1030 of file GNuMIFlux.cxx.

1031 {
1032  LOG("Flux", pNOTICE) << "Initializing GNuMIFlux driver";
1033 
1034  fMaxEv = 0;
1035  fEnd = false;
1036 
1038 
1039  fNuFluxTree = 0;
1040  fG3NuMI = 0;
1041  fG4NuMI = 0;
1042  fFlugg = 0;
1043  fNuFluxTreeName = "";
1044  fNuFluxGen = "";
1045  fNFiles = 0;
1046 
1047  fNEntries = 0;
1048  fIEntry = -1;
1049  fICycle = 0;
1050  fNUse = 1;
1051  fIUse = 999999;
1052 
1053  fNuTot = 0;
1054  fFilePOTs = 0;
1055 
1056  fMaxWeight = -1;
1057  fMaxWgtFudge = 1.05;
1058  fMaxWgtEntries = 2500000;
1059  fMaxEFudge = 0;
1060 
1061  fSumWeight = 0;
1062  fNNeutrinos = 0;
1063  fEffPOTsPerNu = 0;
1064  fAccumPOTs = 0;
1065 
1066  fGenWeighted = false;
1067  fApplyTiltWeight = true;
1068  fUseFluxAtDetCenter = 0;
1069  fDetLocIsSet = false;
1070  // by default assume user length is cm
1072 
1073  this->SetDefaults();
1074  this->ResetCurrent();
1075 }
double fSumWeight
sum of weights for nus thrown so far
Definition: GNuMIFlux.h:408
long int fIUse
current # of times an entry has been used
Definition: GNuMIFlux.h:407
Long64_t fFilePOTs
of protons-on-target represented by all files
Definition: GNuMIFlux.h:398
void SetLengthUnits(double user_units)
Set units assumed by user.
Definition: GNuMIFlux.cxx:1206
Long64_t fIEntry
current flux ntuple entry
Definition: GNuMIFlux.h:396
bool fGenWeighted
does GenerateNext() give weights?
Definition: GNuMIFlux.h:413
double fMaxEFudge
fudge factor for estmating max enu (0=> use fixed 120GeV)
Definition: GNuMIFlux.h:404
bool fApplyTiltWeight
wgt due to window normal not || beam
Definition: GNuMIFlux.h:414
Long64_t fNEntries
number of flux ntuple entries
Definition: GNuMIFlux.h:395
bool fEnd
end condition reached
Definition: GNuMIFlux.h:385
double fMaxWgtFudge
fudge factor for estimating max wgt
Definition: GNuMIFlux.h:402
TChain * fNuFluxTree
TTree in g3numi or g4numi // REF ONLY!
Definition: GNuMIFlux.h:389
string fNuFluxTreeName
Tree name "h10" (g3) or "nudata" (g4)
Definition: GNuMIFlux.h:388
g4numi * fG4NuMI
g4numi ntuple
Definition: GNuMIFlux.h:392
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
double UnitFromString(string u)
Definition: UnitUtils.cxx:18
bool fDetLocIsSet
is a flux location (near/far) set?
Definition: GNuMIFlux.h:415
int fNFiles
number of files in chain
Definition: GNuMIFlux.h:394
flugg * fFlugg
flugg ntuple
Definition: GNuMIFlux.h:393
GNuMIFluxPassThroughInfo * fCurEntry
copy of current ntuple entry info (owned structure)
Definition: GNuMIFlux.h:436
double fAccumPOTs
POTs used so far.
Definition: GNuMIFlux.h:411
double fMaxWeight
max flux neutrino weight in input file
Definition: GNuMIFlux.h:401
long int fMaxWgtEntries
of entries in estimating max wgt
Definition: GNuMIFlux.h:403
g3numi * fG3NuMI
g3numi ntuple
Definition: GNuMIFlux.h:391
int fUseFluxAtDetCenter
use flux at near (-1) or far (+1) det center from ntuple?
Definition: GNuMIFlux.h:416
double fEffPOTsPerNu
what a entry is worth ...
Definition: GNuMIFlux.h:410
#define pNOTICE
Definition: Messenger.h:61
string fNuFluxGen
"g3numi" "g4numi" or "flugg"
Definition: GNuMIFlux.h:390
long int fNUse
how often to use same entry in a row
Definition: GNuMIFlux.h:406
Long64_t fNuTot
cummulative # of entries (=fNEntries)
Definition: GNuMIFlux.h:397
double fMaxEv
maximum energy
Definition: GNuMIFlux.h:384
long int fNNeutrinos
number of flux neutrinos thrown so far
Definition: GNuMIFlux.h:409
double GNuMIFlux::LengthUnits ( void  ) const

Return user units.

Definition at line 1240 of file GNuMIFlux.cxx.

1241 {
1242  // Return the scale factor for lengths the user is getting
1243  double cm = genie::utils::units::UnitFromString("cm");
1244  return fLengthScaleB2U * cm ;
1245 }
static constexpr double cm
Definition: Units.h:68
double UnitFromString(string u)
Definition: UnitUtils.cxx:18
double fLengthScaleB2U
scale factor beam (cm) –> user
Definition: GNuMIFlux.h:419
void GNuMIFlux::LoadBeamSimData ( const std::vector< std::string > &  filenames,
const std::string det_loc 
)
virtual

first is primary method for loading root flux ntuple files and config others are alternatives that can be overloaded but have sensible defaults to fall back to calling the vector version

Implements genie::flux::GFluxFileConfigI.

Definition at line 484 of file GNuMIFlux.cxx.

486 {
487 // Loads in a gnumi beam simulation root file (converted from hbook format)
488 // into the GNuMIFlux driver.
489 
490  bool found_cfg = this->LoadConfig(config);
491  if ( ! found_cfg ) {
492  LOG("Flux", pFATAL)
493  << "LoadBeamSimData could not find XML config \"" << config << "\"\n";
494  exit(1);
495  }
496 
497  fNuFluxFilePatterns = patterns;
498  std::vector<int> nfiles_from_pattern;
499 
500  // create a (sorted) set of file names
501  // this also helps ensure that the same file isn't listed multiple times
502  std::set<std::string> fnames;
503 
504  LOG("Flux", pINFO) << "LoadBeamSimData was passed " << patterns.size()
505  << " patterns";
506 
507  for (size_t ipatt = 0; ipatt < patterns.size(); ++ipatt ) {
508  string pattern = patterns[ipatt];
509  nfiles_from_pattern.push_back(0);
510 
511  LOG("Flux", pNOTICE)
512  << "Loading gnumi flux tree from ROOT file pattern ["
513  << std::setw(3) << ipatt << "] \"" << pattern << "\"";
514 
515  // !WILDCARD only works for file name ... NOT directory
516  string dirname = gSystem->UnixPathName(gSystem->WorkingDirectory());
517  size_t slashpos = pattern.find_last_of("/");
518  size_t fbegin;
519  if ( slashpos != std::string::npos ) {
520  dirname = pattern.substr(0,slashpos);
521  LOG("Flux", pINFO) << "Look for flux using directory " << dirname;
522  fbegin = slashpos + 1;
523  } else { fbegin = 0; }
524 
525  void* dirp = gSystem->OpenDirectory(gSystem->ExpandPathName(dirname.c_str()));
526  if ( dirp ) {
528  pattern.substr(fbegin,pattern.size()-fbegin);
529  TRegexp re(basename.c_str(),kTRUE);
530  const char* onefile;
531  while ( ( onefile = gSystem->GetDirEntry(dirp) ) ) {
532  TString afile = onefile;
533  if ( afile=="." || afile==".." ) continue;
534  if ( basename!=afile && afile.Index(re) == kNPOS ) continue;
535  std::string fullname = dirname + "/" + afile.Data();
536  fnames.insert(fullname);
537  nfiles_from_pattern[ipatt]++;
538  }
539  gSystem->FreeDirectory(dirp);
540  } // legal directory
541  } // loop over patterns
542 
543  size_t indx = 0;
544  std::set<string>::const_iterator sitr = fnames.begin();
545  for ( ; sitr != fnames.end(); ++sitr, ++indx ) {
546  string filename = *sitr;
547  //std::cout << " [" << std::setw(3) << indx << "] \""
548  // << filename << "\"" << std::endl;
549  bool isok = true;
550  // this next test only works for local files, so we can't do that
551  // if we want to stream via xrootd
552  // ! (gSystem->AccessPathName(filename.c_str()));
553  if ( isok ) {
554  // open the file to see what it contains
555  // h10 => g3numi _or_ flugg
556  // nudata => g4numi
557  // for now distinguish between g3numi/flugg using file name
558  TFile* tf = TFile::Open(filename.c_str(),"READ");
559  int isflugg = ( filename.find("flugg") != string::npos ) ? 1 : 0;
560  const std::string tnames[] = { "h10", "nudata" };
561  const std::string gnames[] = { "g3numi","g4numi","flugg","g4flugg"};
562  for (int j = 0; j < 2 ; ++j ) {
563  TTree* atree = (TTree*)tf->Get(tnames[j].c_str());
564  if ( atree ) {
565  const std::string tname_this = tnames[j];
566  const std::string gname_this = gnames[j+2*isflugg];
567  // create chain if none exists
568  if ( ! fNuFluxTree ) {
569  this->SetTreeName(tname_this);
570  fNuFluxTree = new TChain(fNuFluxTreeName.c_str());
571  fNuFluxGen = gname_this;
572  // here we should scan for estimated POTs/file
573  // also maybe the check on max weight
574  }
575  // sanity check for mixing g3/g4/flugg files
576  if ( fNuFluxTreeName != tname_this ||
577  fNuFluxGen != gname_this ) {
578  LOG("Flux", pFATAL)
579  << "Inconsistent flux file types\n"
580  << "The input gnumi flux file \"" << filename
581  << "\"\ncontains a '" << tname_this << "' " << gname_this
582  << "numi ntuple "
583  << "but a '" << fNuFluxTreeName << "' " << fNuFluxGen
584  << " numi ntuple has alread been seen in the chain";
585  exit(1);
586  } // sanity mix/match g3/g4
587  // add the file to the chain
588  this->AddFile(atree,filename);
589  } // found a tree
590  } // loop over either g3 or g4
591  tf->Close();
592  delete tf;
593  } // loop over tree type
594  } // loop over sorted file names
595 
596  if ( fNuFluxTreeName == "" ) {
597  LOG("Flux", pFATAL)
598  << "The input gnumi flux file doesn't exist! Initialization failed!";
599  exit(1);
600  }
601  if ( fNuFluxGen == "g3numi" ) fG3NuMI = new g3numi(fNuFluxTree);
602  if ( fNuFluxGen == "g4numi" ) fG4NuMI = new g4numi(fNuFluxTree);
603  if ( fNuFluxGen == "flugg" ) fFlugg = new flugg(fNuFluxTree);
604 
605  // this will open all files and read header!!
606  fNEntries = fNuFluxTree->GetEntries();
607 
608  if ( fNEntries == 0 ) {
609  LOG("Flux", pERROR)
610  << "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!";
611  LOG("Flux", pERROR)
612  << "Loaded flux tree contains " << fNEntries << " entries";
613  LOG("Flux", pERROR)
614  << "Was passed the file patterns: ";
615  for (size_t ipatt = 0; ipatt < patterns.size(); ++ipatt ) {
616  string pattern = patterns[ipatt];
617  LOG("Flux", pERROR)
618  << " [" << std::setw(3) << ipatt <<"] " << pattern;
619  }
620  LOG("Flux", pERROR)
621  << "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!";
622  } else {
623  LOG("Flux", pNOTICE)
624  << "Loaded flux tree contains " << fNEntries << " entries"
625  << " from " << fnames.size() << " unique files";
626  for (size_t ipatt = 0; ipatt < patterns.size(); ++ipatt ) {
627  string pattern = patterns[ipatt];
628  LOG("Flux", pINFO)
629  << " pattern: " << pattern << " yielded "
630  << nfiles_from_pattern[ipatt] << " files";
631  }
632  }
633 
634  // we have a file we can work with
635  if (!fDetLocIsSet) {
636  LOG("Flux", pERROR)
637  << "LoadBeamSimData left detector location unset";
638  }
639  if (fMaxWeight<=0) {
640  LOG("Flux", pINFO)
641  << "Run ScanForMaxWeight() as part of LoadBeamSimData";
642  this->ScanForMaxWeight();
643  }
644 
645  // current ntuple cycle # (flux ntuples may be recycled)
646  fICycle = 0;
647  // pick a starting entry index [0:fNEntries-1]
648  // pretend we just used up the the previous one
650  fIUse = 9999999;
651  fIEntry = rnd->RndFlux().Integer(fNEntries) - 1;
652 
653  // don't count things we used to estimate max weight
654  fSumWeight = 0;
655  fNNeutrinos = 0;
656  fAccumPOTs = 0;
657 
658  LOG("Flux",pNOTICE) << "about to CalcEffPOTsPerNu";
659  this->CalcEffPOTsPerNu();
660 
661 }
double fSumWeight
sum of weights for nus thrown so far
Definition: GNuMIFlux.h:408
Definition: g3numi.h:15
long int fIUse
current # of times an entry has been used
Definition: GNuMIFlux.h:407
#define pERROR
Definition: Messenger.h:59
static RandomGen * Instance()
Access instance.
Definition: RandomGen.cxx:71
std::string string
Definition: nybbler.cc:12
Long64_t fIEntry
current flux ntuple entry
Definition: GNuMIFlux.h:396
#define pFATAL
Definition: Messenger.h:56
bool LoadConfig(string cfg)
load a named configuration
Definition: GNuMIFlux.cxx:2324
Definition: g4numi.h:18
void SetTreeName(string name)
set input tree name (default: "h10")
Definition: GNuMIFlux.cxx:776
Long64_t fNEntries
number of flux ntuple entries
Definition: GNuMIFlux.h:395
intermediate_table::const_iterator const_iterator
Definition: tf_graph.h:23
string filename
Definition: train.py:213
A singleton holding random number generator classes. All random number generation in GENIE should tak...
Definition: RandomGen.h:29
TChain * fNuFluxTree
TTree in g3numi or g4numi // REF ONLY!
Definition: GNuMIFlux.h:389
string fNuFluxTreeName
Tree name "h10" (g3) or "nudata" (g4)
Definition: GNuMIFlux.h:388
g4numi * fG4NuMI
g4numi ntuple
Definition: GNuMIFlux.h:392
std::vector< string > fNuFluxFilePatterns
(potentially wildcarded) path(s)
Definition: GNuMIFlux.h:387
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
static Config * config
Definition: config.cpp:1054
Definition: flugg.h:15
#define pINFO
Definition: Messenger.h:62
bool fDetLocIsSet
is a flux location (near/far) set?
Definition: GNuMIFlux.h:415
flugg * fFlugg
flugg ntuple
Definition: GNuMIFlux.h:393
double fAccumPOTs
POTs used so far.
Definition: GNuMIFlux.h:411
Q_EXPORT QTSManip setw(int w)
Definition: qtextstream.h:331
double fMaxWeight
max flux neutrino weight in input file
Definition: GNuMIFlux.h:401
g3numi * fG3NuMI
g3numi ntuple
Definition: GNuMIFlux.h:391
std::string pattern
Definition: regex_t.cc:35
void ScanForMaxWeight(void)
scan for max flux weight (before generating unweighted flux neutrinos)
Definition: GNuMIFlux.cxx:680
#define pNOTICE
Definition: Messenger.h:61
TRandom3 & RndFlux(void) const
rnd number generator used by flux drivers
Definition: RandomGen.h:71
string fNuFluxGen
"g3numi" "g4numi" or "flugg"
Definition: GNuMIFlux.h:390
void AddFile(TTree *tree, string fname)
Definition: GNuMIFlux.cxx:1133
long int fNNeutrinos
number of flux neutrinos thrown so far
Definition: GNuMIFlux.h:409
void CalcEffPOTsPerNu(void)
Definition: GNuMIFlux.cxx:439
bool GNuMIFlux::LoadConfig ( string  cfg)

load a named configuration

Definition at line 2324 of file GNuMIFlux.cxx.

2325 {
2326  const char* altxml = gSystem->Getenv("GNUMIFLUXXML");
2327  if ( altxml ) {
2328  SetXMLFileBase(altxml);
2329  }
2330  genie::flux::GNuMIFluxXMLHelper helper(this);
2331  return helper.LoadConfig(cfg);
2332 }
virtual void SetXMLFileBase(std::string xmlbasename="")
double genie::flux::GNuMIFlux::MaxEnergy ( void  )
inlinevirtual

declare the max flux neutrino energy that can be generated (for init. purposes)

Implements genie::GFluxI.

Definition at line 231 of file GNuMIFlux.h.

231 { return fMaxEv; }
double fMaxEv
maximum energy
Definition: GNuMIFlux.h:384
const TLorentzVector& genie::flux::GNuMIFlux::Momentum ( void  )
inlinevirtual

returns the flux neutrino 4-momentum

Implements genie::GFluxI.

Definition at line 235 of file GNuMIFlux.h.

235 { return fCurEntry->fgP4User; }
GNuMIFluxPassThroughInfo * fCurEntry
copy of current ntuple entry info (owned structure)
Definition: GNuMIFlux.h:436
TLorentzVector fgP4User
generated nu 4-momentum user coord
Definition: GNuMIFlux.h:104
void GNuMIFlux::MoveToZ0 ( double  z0)

move ray origin to user coord Z0

Definition at line 403 of file GNuMIFlux.cxx.

404 {
405  // move ray origin to specified user z0
406  // move beam coord entry correspondingly
407 
408 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
409  LOG("Flux", pNOTICE)
410  << "MoveToZ0 (z0usr=" << z0usr << ") before:"
411  << "\n p4 user: " << utils::print::P4AsShortString(&(fCurEntry->fgP4User))
412  << "\n x4 user: " << utils::print::X4AsString(&(fCurEntry->fgX4User));
413 #endif
414 
415  double pzusr = fCurEntry->fgP4User.Pz();
416  if ( TMath::Abs(pzusr) < 1.0e-30 ) {
417  // neutrino is moving almost entirely in x-y plane
418  LOG("Flux", pWARN)
419  << "MoveToZ0(" << z0usr << ") not possible due to pz_usr (" << pzusr << ")";
420  return;
421  }
422 
423  double scale = (z0usr - fCurEntry->fgX4User.Z()) / pzusr;
424  fCurEntry->fgX4User += (scale*fCurEntry->fgP4User);
426  // this scaling works for distances, but not the time component
427  fCurEntry->fgX4.SetT(0);
428  fCurEntry->fgX4User.SetT(0);
429 
430 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
431  LOG("Flux", pNOTICE)
432  << "MoveToZ0 (" << z0usr << ") after:"
433  << "\n x4 user: " << utils::print::X4AsString(&(fCurEntry->fgX4User));
434 #endif
435 
436 }
TLorentzVector fgP4
generated nu 4-momentum beam coord
Definition: GNuMIFlux.h:102
string P4AsShortString(const TLorentzVector *p)
Definition: PrintUtils.cxx:45
TLorentzVector fgX4User
generated nu 4-position user coord
Definition: GNuMIFlux.h:105
TLorentzVector fgX4
generated nu 4-position beam coord
Definition: GNuMIFlux.h:103
const double e
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
GNuMIFluxPassThroughInfo * fCurEntry
copy of current ntuple entry info (owned structure)
Definition: GNuMIFlux.h:436
#define pWARN
Definition: Messenger.h:60
double fLengthScaleU2B
scale factor beam user –> (cm)
Definition: GNuMIFlux.h:420
#define pNOTICE
Definition: Messenger.h:61
TLorentzVector fgP4User
generated nu 4-momentum user coord
Definition: GNuMIFlux.h:104
string X4AsString(const TLorentzVector *x)
Definition: PrintUtils.cxx:57
long int GNuMIFlux::NFluxNeutrinos ( void  ) const
virtual

of rays generated

number of flux neutrinos ray generated so far

Implements genie::flux::GFluxExposureI.

Definition at line 132 of file GNuMIFlux.cxx.

133 {
134  /// number of flux neutrinos ray generated so far
135  return fNNeutrinos;
136 }
long int fNNeutrinos
number of flux neutrinos thrown so far
Definition: GNuMIFlux.h:409
const GNuMIFluxPassThroughInfo& genie::flux::GNuMIFlux::PassThroughInfo ( void  )
inline

GNuMIFluxPassThroughInfo.

Definition at line 252 of file GNuMIFlux.h.

int genie::flux::GNuMIFlux::PdgCode ( void  )
inlinevirtual

returns the flux neutrino pdg code

Implements genie::GFluxI.

Definition at line 233 of file GNuMIFlux.h.

233 { return fCurEntry->fgPdgC; }
GNuMIFluxPassThroughInfo * fCurEntry
copy of current ntuple entry info (owned structure)
Definition: GNuMIFlux.h:436
int fgPdgC
generated nu pdg-code
Definition: GNuMIFlux.h:98
const TLorentzVector& genie::flux::GNuMIFlux::Position ( void  )
inlinevirtual

returns the flux neutrino 4-position (note: expect SI rather than physical units)

Implements genie::GFluxI.

Definition at line 236 of file GNuMIFlux.h.

236 { return fCurEntry->fgX4User; }
TLorentzVector fgX4User
generated nu 4-position user coord
Definition: GNuMIFlux.h:105
GNuMIFluxPassThroughInfo * fCurEntry
copy of current ntuple entry info (owned structure)
Definition: GNuMIFlux.h:436
double GNuMIFlux::POT_curr ( void  )

current average POT (RWH?)

Definition at line 477 of file GNuMIFlux.cxx.

477  {
478  // RWH: Not sure what POT_curr is supposed to represent I'll guess for
479  // now that that it means what I mean by UsedPOTs().
480  return UsedPOTs();
481 }
double UsedPOTs(void) const
of protons-on-target used
Definition: GNuMIFlux.cxx:464
void GNuMIFlux::PrintConfig ( void  )
virtual

print the current configuration

Implements genie::flux::GFluxFileConfigI.

Definition at line 2336 of file GNuMIFlux.cxx.

2337 {
2338 
2339  std::ostringstream s;
2340  PDGCodeList::const_iterator itr = fPdgCList->begin();
2341  for ( ; itr != fPdgCList->end(); ++itr) s << (*itr) << " ";
2342  s << "[rejected: ";
2343  itr = fPdgCListRej->begin();
2344  for ( ; itr != fPdgCListRej->end(); ++itr) s << (*itr) << " ";
2345  s << " ] ";
2346 
2347  std::ostringstream fpattout;
2348  for (size_t i = 0; i < fNuFluxFilePatterns.size(); ++i)
2349  fpattout << "\n [" << std::setw(3) << i << "] " << fNuFluxFilePatterns[i];
2350 
2351  std::ostringstream flistout;
2352  std::vector<std::string> flist = GetFileList();
2353  for (size_t i = 0; i < flist.size(); ++i)
2354  flistout << "\n [" << std::setw(3) << i << "] " << flist[i];
2355 
2356  TLorentzVector usr0(0,0,0,0);
2357  TLorentzVector usr0asbeam;
2358  User2BeamPos(usr0,usr0asbeam);
2359 
2360  const int w=10, p=6;
2361  std::ostringstream beamrot_str, beamrotinv_str;
2362  beamrot_str
2363  << "fBeamRot: " << std::setprecision(p) << "\n"
2364  << " [ "
2365  << std::setw(w) << fBeamRot.XX() << " "
2366  << std::setw(w) << fBeamRot.XY() << " "
2367  << std::setw(w) << fBeamRot.XZ() << " ]\n"
2368  << " [ "
2369  << std::setw(w) << fBeamRot.YX() << " "
2370  << std::setw(w) << fBeamRot.YY() << " "
2371  << std::setw(w) << fBeamRot.YZ() << " ]\n"
2372  << " [ "
2373  << std::setw(w) << fBeamRot.ZX() << " "
2374  << std::setw(w) << fBeamRot.ZY() << " "
2375  << std::setw(w) << fBeamRot.ZZ() << " ]";
2376  beamrotinv_str
2377  << "fBeamRotInv: " << std::setprecision(p) << "\n"
2378  << " [ "
2379  << std::setw(w) << fBeamRotInv.XX() << " "
2380  << std::setw(w) << fBeamRotInv.XY() << " "
2381  << std::setw(w) << fBeamRotInv.XZ() << " ]\n"
2382  << " [ "
2383  << std::setw(w) << fBeamRotInv.YX() << " "
2384  << std::setw(w) << fBeamRotInv.YY() << " "
2385  << std::setw(w) << fBeamRotInv.YZ() << " ]\n"
2386  << " [ "
2387  << std::setw(w) << fBeamRotInv.ZX() << " "
2388  << std::setw(w) << fBeamRotInv.ZY() << " "
2389  << std::setw(w) << fBeamRotInv.ZZ() << " ]";
2390 
2391  LOG("Flux", pNOTICE)
2392  << "GNuMIFlux Config:"
2393  << "\n Enu_max " << fMaxEv
2394  << "\n pdg-codes: " << s.str() << "\n "
2395  << (fG3NuMI?"g3numi":"")
2396  << (fG4NuMI?"g4numi":"")
2397  << (fFlugg?"flugg":"")
2398  << "/" << fNuFluxGen << " "
2399  << "(" << fNuFluxTreeName << "), " << fNEntries << " entries"
2400  << " (FilePOTs " << fFilePOTs << ") "
2401  << "in " << fNFiles << " files: "
2402  << flistout.str()
2403  << "\n from file patterns:"
2404  << fpattout.str()
2405  << "\n wgt max=" << fMaxWeight << " fudge=" << fMaxWgtFudge << " using "
2406  << fMaxWgtEntries << " entries"
2407  << "\n Z0 pushback " << fZ0
2408  << "\n used entry " << fIEntry << " " << fIUse << "/" << fNUse
2409  << " times, in " << fICycle << "/" << fNCycles << " cycles"
2410  << "\n SumWeight " << fSumWeight << " for " << fNNeutrinos << " neutrinos"
2411  << "\n EffPOTsPerNu " << fEffPOTsPerNu << " AccumPOTs " << fAccumPOTs
2412  << "\n GenWeighted \"" << (fGenWeighted?"true":"false") << ", "
2413  << "\", Detector location set \"" << (fDetLocIsSet?"true":"false") << "\""
2414  << "\n LengthUnits " << fLengthUnits << ", scale b2u " << fLengthScaleB2U
2415  << ", scale u2b " << fLengthScaleU2B
2416  << "\n User Flux Window: "
2420  << "\n Flux Window (cm, beam coord): "
2421  << "\n base " << utils::print::X4AsString(&fFluxWindowBase)
2422  << "\n dir1 " << utils::print::X4AsString(&fFluxWindowDir1) << " len " << fFluxWindowLen1
2423  << "\n dir2 " << utils::print::X4AsString(&fFluxWindowDir2) << " len " << fFluxWindowLen2
2424  << "\n normal " << utils::print::Vec3AsString(&(fWindowNormal))
2425  << "\n User Beam Origin: "
2426  << "\n base " << utils::print::X4AsString(&fBeamZero)
2427  << "\n " << beamrot_str.str() << " "
2428  << "\n Detector Origin (beam coord): "
2429  << "\n base " << utils::print::X4AsString(&usr0asbeam)
2430  << "\n " << beamrotinv_str.str() << " "
2431  << "\n UseFluxAtDetCenter " << fUseFluxAtDetCenter;
2432 
2433 }
double fSumWeight
sum of weights for nus thrown so far
Definition: GNuMIFlux.h:408
double fLengthUnits
units for coord in user exchanges
Definition: GNuMIFlux.h:418
long int fIUse
current # of times an entry has been used
Definition: GNuMIFlux.h:407
Long64_t fFilePOTs
of protons-on-target represented by all files
Definition: GNuMIFlux.h:398
Long64_t fIEntry
current flux ntuple entry
Definition: GNuMIFlux.h:396
bool fGenWeighted
does GenerateNext() give weights?
Definition: GNuMIFlux.h:413
TVector3 fFluxWindowPtUser[3]
user points of flux window
Definition: GNuMIFlux.h:426
TLorentzVector fFluxWindowDir2
extent for flux window (direction 2)
Definition: GNuMIFlux.h:429
Long64_t fNEntries
number of flux ntuple entries
Definition: GNuMIFlux.h:395
intermediate_table::const_iterator const_iterator
TLorentzVector fBeamZero
beam origin in user coords
Definition: GNuMIFlux.h:422
Q_EXPORT QTSManip setprecision(int p)
Definition: qtextstream.h:343
double fMaxWgtFudge
fudge factor for estimating max wgt
Definition: GNuMIFlux.h:402
string fNuFluxTreeName
Tree name "h10" (g3) or "nudata" (g4)
Definition: GNuMIFlux.h:388
void User2BeamPos(const TLorentzVector &usrxyz, TLorentzVector &beamxyz) const
Definition: GNuMIFlux.cxx:987
g4numi * fG4NuMI
g4numi ntuple
Definition: GNuMIFlux.h:392
std::vector< string > fNuFluxFilePatterns
(potentially wildcarded) path(s)
Definition: GNuMIFlux.h:387
PDGCodeList * fPdgCList
list of neutrino pdg-codes to generate
std::vector< std::string > GetFileList()
list of files currently part of chain
Definition: GNuMIFlux.cxx:2436
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
TLorentzVector fFluxWindowDir1
extent for flux window (direction 1)
Definition: GNuMIFlux.h:428
p
Definition: test.py:223
bool fDetLocIsSet
is a flux location (near/far) set?
Definition: GNuMIFlux.h:415
int fNFiles
number of files in chain
Definition: GNuMIFlux.h:394
long int fNCycles
times to cycle through the ntuple(s)
flugg * fFlugg
flugg ntuple
Definition: GNuMIFlux.h:393
double fAccumPOTs
POTs used so far.
Definition: GNuMIFlux.h:411
Q_EXPORT QTSManip setw(int w)
Definition: qtextstream.h:331
PDGCodeList * fPdgCListRej
list of nu pdg-codes seen but rejected
double fMaxWeight
max flux neutrino weight in input file
Definition: GNuMIFlux.h:401
double fLengthScaleB2U
scale factor beam (cm) –> user
Definition: GNuMIFlux.h:419
long int fMaxWgtEntries
of entries in estimating max wgt
Definition: GNuMIFlux.h:403
g3numi * fG3NuMI
g3numi ntuple
Definition: GNuMIFlux.h:391
TVector3 fWindowNormal
normal direction for flux window
Definition: GNuMIFlux.h:432
double fLengthScaleU2B
scale factor beam user –> (cm)
Definition: GNuMIFlux.h:420
int fUseFluxAtDetCenter
use flux at near (-1) or far (+1) det center from ntuple?
Definition: GNuMIFlux.h:416
double fEffPOTsPerNu
what a entry is worth ...
Definition: GNuMIFlux.h:410
#define pNOTICE
Definition: Messenger.h:61
TLorentzRotation fBeamRotInv
Definition: GNuMIFlux.h:424
string fNuFluxGen
"g3numi" "g4numi" or "flugg"
Definition: GNuMIFlux.h:390
string X4AsString(const TLorentzVector *x)
Definition: PrintUtils.cxx:57
TLorentzRotation fBeamRot
rotation applied beam –> user coord
Definition: GNuMIFlux.h:423
long int fNUse
how often to use same entry in a row
Definition: GNuMIFlux.h:406
string Vec3AsString(const TVector3 *vec)
Definition: PrintUtils.cxx:80
double fMaxEv
maximum energy
Definition: GNuMIFlux.h:384
static QCString * s
Definition: config.cpp:1042
TLorentzVector fFluxWindowBase
base point for flux window - beam coord
Definition: GNuMIFlux.h:427
long int fNNeutrinos
number of flux neutrinos thrown so far
Definition: GNuMIFlux.h:409
void GNuMIFlux::PrintCurrent ( void  )

print current entry from leaves

Definition at line 1004 of file GNuMIFlux.cxx.

1005 {
1006  LOG("Flux", pNOTICE) << "CurrentEntry:" << *fCurEntry;
1007 }
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
GNuMIFluxPassThroughInfo * fCurEntry
copy of current ntuple entry info (owned structure)
Definition: GNuMIFlux.h:436
#define pNOTICE
Definition: Messenger.h:61
void GNuMIFlux::ResetCurrent ( void  )
private

Definition at line 1106 of file GNuMIFlux.cxx.

1107 {
1108 // reset running values of neutrino pdg-code, 4-position & 4-momentum
1109 // and the input ntuple leaves
1110 
1112  fCurEntry->ResetCopy();
1113 }
GNuMIFluxPassThroughInfo * fCurEntry
copy of current ntuple entry info (owned structure)
Definition: GNuMIFlux.h:436
void GNuMIFlux::ScanForMaxWeight ( void  )

scan for max flux weight (before generating unweighted flux neutrinos)

Definition at line 680 of file GNuMIFlux.cxx.

681 {
682  if (!fDetLocIsSet) {
683  LOG("Flux", pERROR)
684  << "Specify a detector location before scanning for max weight";
685  return;
686  }
687 
688  // scan for the maximum weight
689  int ipos_estimator = fUseFluxAtDetCenter;
690  if ( ipos_estimator == 0 ) {
691  // within 100m of a known point?
692  double zbase = fFluxWindowBase.Z();
693  if ( TMath::Abs(zbase-103648.837) < 10000. ) ipos_estimator = -1; // use NearDet
694  if ( TMath::Abs(zbase-73534000. ) < 10000. ) ipos_estimator = +1; // use FarDet
695  }
696  if ( ipos_estimator != 0 ) {
697 
698  //// one can't really be sure which Nwtfar/Nwtnear this refers to
699  //// some gnumi files have "NOvA" weights
700  const char* ntwgtstrv[4] = { "Nimpwt*Nwtnear",
701  "Nimpwt*Nwtfar",
702  "Nimpwt*NWtNear[0]",
703  "Nimpwt*NWtFar[0]" };
704  int strindx = 0;
705  if ( ipos_estimator > 0 ) strindx = 1;
706  if ( fG4NuMI ) strindx += 2;
707  // set upper limit on how many entries to scan
708  Long64_t nscan = TMath::Min(fNEntries,200000LL);
709 
710  fNuFluxTree->Draw(ntwgtstrv[strindx],"","goff",nscan);
711  //std::cout << " Draw \"" << ntwgtstrv[strindx] << "\"" << std::endl;
712  //std::cout << " SelectedRows " << fNuFluxTree->GetSelectedRows()
713  // << " V1 " << fNuFluxTree->GetV1() << std::endl;
714 
715  Long64_t idx = TMath::LocMax(fNuFluxTree->GetSelectedRows(),
716  fNuFluxTree->GetV1());
717  //std::cout << "idx " << idx << " of " << fNuFluxTree->GetSelectedRows() << std::endl;
718  fMaxWeight = fNuFluxTree->GetV1()[idx];
719  LOG("Flux", pNOTICE) << "Maximum flux weight from Nwt in ntuple = "
720  << fMaxWeight;
721  if ( fMaxWeight <= 0 ) {
722  LOG("Flux", pFATAL) << "Non-positive maximum flux weight!";
723  exit(1);
724  }
725  }
726  // the above works only for things close to the MINOS stored weight
727  // values. otherwise we need to work out our own estimate.
728  double wgtgenmx = 0, enumx = 0;
729  TStopwatch t;
730  t.Start();
731  for (int itry=0; itry < fMaxWgtEntries; ++itry) {
732  this->GenerateNext_weighted();
733  double wgt = this->Weight();
734  if ( wgt > wgtgenmx ) wgtgenmx = wgt;
735  double enu = fCurEntry->fgP4.Energy();
736  if ( enu > enumx ) enumx = enu;
737  }
738  t.Stop();
739  t.Print("u");
740  LOG("Flux", pNOTICE) << "Maximum flux weight for spin = "
741  << wgtgenmx << ", energy = " << enumx
742  << " (" << fMaxWgtEntries << ")";
743 
744  if (wgtgenmx > fMaxWeight ) fMaxWeight = wgtgenmx;
745  // apply a fudge factor to estimated weight
746  fMaxWeight *= fMaxWgtFudge;
747  // adjust max energy?
748  if ( enumx*fMaxEFudge > fMaxEv ) {
749  LOG("Flux", pNOTICE) << "Adjust max: was=" << fMaxEv
750  << " now " << enumx << "*" << fMaxEFudge
751  << " = " << enumx*fMaxEFudge;
752  fMaxEv = enumx * fMaxEFudge;
753  }
754 
755  LOG("Flux", pNOTICE) << "Maximum flux weight = " << fMaxWeight
756  << ", energy = " << fMaxEv;
757 
758 }
TLorentzVector fgP4
generated nu 4-momentum beam coord
Definition: GNuMIFlux.h:102
#define pERROR
Definition: Messenger.h:59
double fMaxEFudge
fudge factor for estmating max enu (0=> use fixed 120GeV)
Definition: GNuMIFlux.h:404
#define pFATAL
Definition: Messenger.h:56
Long64_t fNEntries
number of flux ntuple entries
Definition: GNuMIFlux.h:395
double fMaxWgtFudge
fudge factor for estimating max wgt
Definition: GNuMIFlux.h:402
TChain * fNuFluxTree
TTree in g3numi or g4numi // REF ONLY!
Definition: GNuMIFlux.h:389
double Weight(void)
returns the flux neutrino weight (if any)
Definition: GNuMIFlux.h:234
g4numi * fG4NuMI
g4numi ntuple
Definition: GNuMIFlux.h:392
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
bool GenerateNext_weighted(void)
Definition: GNuMIFlux.cxx:204
bool fDetLocIsSet
is a flux location (near/far) set?
Definition: GNuMIFlux.h:415
GNuMIFluxPassThroughInfo * fCurEntry
copy of current ntuple entry info (owned structure)
Definition: GNuMIFlux.h:436
double fMaxWeight
max flux neutrino weight in input file
Definition: GNuMIFlux.h:401
long int fMaxWgtEntries
of entries in estimating max wgt
Definition: GNuMIFlux.h:403
int fUseFluxAtDetCenter
use flux at near (-1) or far (+1) det center from ntuple?
Definition: GNuMIFlux.h:416
#define pNOTICE
Definition: Messenger.h:61
double fMaxEv
maximum energy
Definition: GNuMIFlux.h:384
TLorentzVector fFluxWindowBase
base point for flux window - beam coord
Definition: GNuMIFlux.h:427
void genie::flux::GNuMIFlux::SetApplyWindowTiltWeight ( bool  apply = true)
inline
Parameters
applyapply wgt due to tilt of flux window relative to beam

Definition at line 302 of file GNuMIFlux.h.

303  { fApplyTiltWeight = apply; }
bool fApplyTiltWeight
wgt due to window normal not || beam
Definition: GNuMIFlux.h:414
void GNuMIFlux::SetBeamCenter ( TVector3  beam0)

Definition at line 920 of file GNuMIFlux.cxx.

921 {
922  // set coord transform between detector and beam
923  // NOTE: internally these are in "cm", but user might have set a preference
924  fBeamZero = TLorentzVector(beam0,0); // no time shift
925 }
TLorentzVector fBeamZero
beam origin in user coords
Definition: GNuMIFlux.h:422
void GNuMIFlux::SetBeamRotation ( TRotation  beamrot)

< beam (0,0,0) relative to user frame, beam direction in user frame

Definition at line 913 of file GNuMIFlux.cxx.

914 {
915  // rotation is really only 3-d vector, but we'll be operating on LorentzV's
916  fBeamRot = TLorentzRotation(beamrot);
917  fBeamRotInv = fBeamRot.Inverse();
918 }
TLorentzRotation fBeamRotInv
Definition: GNuMIFlux.h:424
TLorentzRotation fBeamRot
rotation applied beam –> user coord
Definition: GNuMIFlux.h:423
void GNuMIFlux::SetDefaults ( void  )
private

Definition at line 1077 of file GNuMIFlux.cxx.

1078 {
1079 // - Set default neutrino species list (nue, nuebar, numu, numubar) and
1080 // maximum energy (120 GeV).
1081 // These defaults can be overwritten by user calls (at the driver init) to
1082 // GNuMIlux::SetMaxEnergy(double Ev) and
1083 // GNuMIFlux::SetFluxParticles(const PDGCodeList & particles)
1084 // - Set the default file normalization to 1E+21 POT
1085 // - Set the default flux neutrino start z position at -5m (z=0 is the
1086 // detector centre).
1087 // - Set number of cycles to 1
1088 
1089  LOG("Flux", pNOTICE) << "Setting default GNuMIFlux driver options";
1090 
1091  PDGCodeList particles;
1092  particles.push_back(kPdgNuMu);
1093  particles.push_back(kPdgAntiNuMu);
1094  particles.push_back(kPdgNuE);
1095  particles.push_back(kPdgAntiNuE);
1096 
1097  this->SetFluxParticles (particles);
1098  this->SetMaxEnergy (120./*GeV*/); // was 200, but that would be wasteful
1099  this->SetUpstreamZ (-3.4e38); // way upstream ==> use flux window
1100  this->SetNumOfCycles (0);
1101  this->SetEntryReuse (1);
1102 
1103  this->SetXMLFileBase("GNuMIFlux.xml");
1104 }
void SetEntryReuse(long int nuse=1)
of times to use entry before moving to next
Definition: GNuMIFlux.cxx:768
const int kPdgNuE
Definition: PDGCodes.h:28
const int kPdgAntiNuE
Definition: PDGCodes.h:29
const int kPdgNuMu
Definition: PDGCodes.h:30
virtual void SetNumOfCycles(long int ncycle)
limit cycling through input files
A list of PDG codes.
Definition: PDGCodeList.h:32
virtual void SetXMLFileBase(std::string xmlbasename="")
virtual void SetUpstreamZ(double z0)
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
const int kPdgAntiNuMu
Definition: PDGCodes.h:31
virtual void SetFluxParticles(const PDGCodeList &particles)
specify list of flux neutrino species
void SetMaxEnergy(double Ev)
specify maximum flx neutrino energy
Definition: GNuMIFlux.cxx:760
#define pNOTICE
Definition: Messenger.h:61
void push_back(int pdg_code)
Definition: PDGCodeList.cxx:58
void GNuMIFlux::SetEntryReuse ( long int  nuse = 1)

of times to use entry before moving to next

Definition at line 768 of file GNuMIFlux.cxx.

769 {
770 // With nuse > 1 then the same entry in the file is used "nuse" times
771 // before moving on to the next entry in the ntuple
772 
773  fNUse = TMath::Max(1L, nuse);
774 }
long int fNUse
how often to use same entry in a row
Definition: GNuMIFlux.h:406
bool GNuMIFlux::SetFluxWindow ( GNuMIFlux::StdFluxWindow_t  stdwindow,
double  padding = 0 
)

return false if unhandled

Definition at line 801 of file GNuMIFlux.cxx.

802 {
803  // set some standard flux windows
804  // rwh: should also set detector coord transform
805  // rwh: padding allows add constant padding to pre-existing set
806  double padbeam = padding / fLengthScaleB2U; // user might set different units
807  LOG("Flux",pNOTICE)
808  << "SetBeamFluxWindow " << (int)stdwindow << " padding " << padbeam << " cm";
809 
810 
811  switch ( stdwindow ) {
812 #ifdef THIS_IS_NOT_YET_IMPLEMENTED
813  case kMinosNearDet:
814  SetBeamFluxWindow(103648.837);
815  break;
816  case kMinosFarDear:
817  SetBeamFluxWindow(73534000.);
818  break;
819  case kMinosNearRock:
820  SetBeamFluxWindow();
821  break;
822  case kMinosFarRock:
823  SetBeamFluxWindow();
824  break;
825 #endif
826  case kMinosNearCenter:
827  {
829  fFluxWindowBase = kPosCenterNearBeam;
830  fFluxWindowDir1 = TLorentzVector(); // no extent
831  fFluxWindowDir2 = TLorentzVector();
832  TLorentzVector usrpos;
833  Beam2UserPos(fFluxWindowBase, usrpos);
834  fFluxWindowPtUser[0] = usrpos.Vect();
837  fFluxWindowLen1 = 0;
838  fFluxWindowLen2 = 0;
839  break;
840  }
841  case kMinosFarCenter:
842  {
845  fFluxWindowDir1 = TLorentzVector(); // no extent
846  fFluxWindowDir2 = TLorentzVector();
847  TLorentzVector usrpos;
848  Beam2UserPos(fFluxWindowBase, usrpos);
849  fFluxWindowPtUser[0] = usrpos.Vect();
852  fFluxWindowLen1 = 0;
853  fFluxWindowLen2 = 0;
854  break;
855  }
856  default:
857  LOG("Flux", pERROR)
858  << "SetBeamFluxWindow - StdFluxWindow " << stdwindow
859  << " not yet implemented";
860  return false;
861  }
862  LOG("Flux",pNOTICE) << "about to CalcEffPOTsPerNu";
863  this->CalcEffPOTsPerNu();
864  return true;
865 }
#define pERROR
Definition: Messenger.h:59
void SetLengthUnits(double user_units)
Set units assumed by user.
Definition: GNuMIFlux.cxx:1206
TVector3 fFluxWindowPtUser[3]
user points of flux window
Definition: GNuMIFlux.h:426
TLorentzVector fFluxWindowDir2
extent for flux window (direction 2)
Definition: GNuMIFlux.h:429
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
TLorentzVector fFluxWindowDir1
extent for flux window (direction 1)
Definition: GNuMIFlux.h:428
double UnitFromString(string u)
Definition: UnitUtils.cxx:18
void Beam2UserPos(const TLorentzVector &beamxyz, TLorentzVector &usrxyz) const
Definition: GNuMIFlux.cxx:971
const TLorentzVector kPosCenterFarBeam(0., 0., 735340.00, 0.)
double fLengthScaleB2U
scale factor beam (cm) –> user
Definition: GNuMIFlux.h:419
#define pNOTICE
Definition: Messenger.h:61
TLorentzVector fFluxWindowBase
base point for flux window - beam coord
Definition: GNuMIFlux.h:427
void CalcEffPOTsPerNu(void)
Definition: GNuMIFlux.cxx:439
void GNuMIFlux::SetFluxWindow ( TVector3  p1,
TVector3  p2,
TVector3  p3 
)

3 points define a plane (by default in user coordinates)

Definition at line 867 of file GNuMIFlux.cxx.

869 {
870  // set flux window
871  // NOTE: internally these are in "cm", but user might have set a preference
873  fDetLocIsSet = true;
874 
875  fFluxWindowPtUser[0] = p0;
876  fFluxWindowPtUser[1] = p1;
877  fFluxWindowPtUser[2] = p2;
878 
879  // convert from user to beam coord and from 3 points to base + 2 directions
880  // apply units conversion
881  TLorentzVector ptbm0, ptbm1, ptbm2;
882  User2BeamPos(TLorentzVector(fFluxWindowPtUser[0],0),ptbm0);
883  User2BeamPos(TLorentzVector(fFluxWindowPtUser[1],0),ptbm1);
884  User2BeamPos(TLorentzVector(fFluxWindowPtUser[2],0),ptbm2);
885 
886  fFluxWindowBase = ptbm0;
887  fFluxWindowDir1 = ptbm1 - ptbm0;
888  fFluxWindowDir2 = ptbm2 - ptbm0;
889 
892  fWindowNormal = fFluxWindowDir1.Vect().Cross(fFluxWindowDir2.Vect()).Unit();
893 
894  double dot = fFluxWindowDir1.Dot(fFluxWindowDir2);
895  if ( TMath::Abs(dot) > 1.0e-8 )
896  LOG("Flux",pWARN) << "Dot product between window direction vectors was "
897  << dot << "; please check for orthoganality";
898 
899  LOG("Flux",pNOTICE) << "about to CalcEffPOTsPerNu";
900  this->CalcEffPOTsPerNu();
901 }
TVector3 fFluxWindowPtUser[3]
user points of flux window
Definition: GNuMIFlux.h:426
TLorentzVector fFluxWindowDir2
extent for flux window (direction 2)
Definition: GNuMIFlux.h:429
void User2BeamPos(const TLorentzVector &usrxyz, TLorentzVector &beamxyz) const
Definition: GNuMIFlux.cxx:987
const double e
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
TLorentzVector fFluxWindowDir1
extent for flux window (direction 1)
Definition: GNuMIFlux.h:428
bool fDetLocIsSet
is a flux location (near/far) set?
Definition: GNuMIFlux.h:415
#define pWARN
Definition: Messenger.h:60
TVector3 fWindowNormal
normal direction for flux window
Definition: GNuMIFlux.h:432
int fUseFluxAtDetCenter
use flux at near (-1) or far (+1) det center from ntuple?
Definition: GNuMIFlux.h:416
#define pNOTICE
Definition: Messenger.h:61
TLorentzVector fFluxWindowBase
base point for flux window - beam coord
Definition: GNuMIFlux.h:427
void CalcEffPOTsPerNu(void)
Definition: GNuMIFlux.cxx:439
void genie::flux::GNuMIFlux::SetGenWeighted ( bool  genwgt = false)
inline

toggle whether GenerateNext() returns weight=1 flux (initial default false)

Definition at line 292 of file GNuMIFlux.h.

void GNuMIFlux::SetLengthUnits ( double  user_units)

Set units assumed by user.

Definition at line 1206 of file GNuMIFlux.cxx.

1207 {
1208  // Set the scale factor for lengths going from beam (cm) to user coordinates
1209 
1210  // GNuMIFlux uses "cm" as the length unit consistently internally (this is
1211  // the length units used by both the g3 and g4 ntuples). User interactions
1212  // setting the beam-to-detector coordinate transform, flux window, and the
1213  // returned position might need to be in other units. Use:
1214  // double scale = genie::utils::units::UnitFromString("cm");
1215  // ( #include "Utils/UnitUtils.h for declaration )
1216  // to get the correct scale factor to pass in.
1217 
1218  double rescale = fLengthUnits / user_units;
1219  fLengthUnits = user_units;
1220  double cm = genie::utils::units::UnitFromString("cm");
1221  fLengthScaleB2U = cm / user_units;
1222  fLengthScaleU2B = user_units / cm;
1223 
1224  // in case we're changing units without resetting transform/window
1225  // not recommended, but should work
1227  fBeamZero *= rescale;
1228  fFluxWindowPtUser[0] *= rescale;
1229  fFluxWindowPtUser[1] *= rescale;
1230  fFluxWindowPtUser[2] *= rescale;
1231 
1232  // case GNuMIFlux::kmeter: fLengthScaleB2U = 0.01 ; break;
1233  // case GNuMIFlux::kcm: fLengthScaleB2U = 1. ; break;
1234  // case GNuMIFlux::kmm: fLengthScaleB2U = 10. ; break;
1235  // case GNuMIFlux::kfm: fLengthScaleB2U = 1.e13 ; break; // 10e-2m -> 10e-15m
1236 
1237 }
static constexpr double cm
Definition: Units.h:68
double fLengthUnits
units for coord in user exchanges
Definition: GNuMIFlux.h:418
TVector3 fFluxWindowPtUser[3]
user points of flux window
Definition: GNuMIFlux.h:426
TLorentzVector fgX4User
generated nu 4-position user coord
Definition: GNuMIFlux.h:105
TLorentzVector fBeamZero
beam origin in user coords
Definition: GNuMIFlux.h:422
double UnitFromString(string u)
Definition: UnitUtils.cxx:18
GNuMIFluxPassThroughInfo * fCurEntry
copy of current ntuple entry info (owned structure)
Definition: GNuMIFlux.h:436
double fLengthScaleB2U
scale factor beam (cm) –> user
Definition: GNuMIFlux.h:419
Quantity< ScaledUnit< typename Q::unit_t::baseunit_t, R >, T > rescale
Type of a quantity like Q, but with a different unit scale R.
Definition: quantities.h:927
double fLengthScaleU2B
scale factor beam user –> (cm)
Definition: GNuMIFlux.h:420
void genie::flux::GNuMIFlux::SetMaxEFudge ( double  fudge = 1.05)
inline
Parameters
fudgeextra fudge factor in estimating maximum energy

Definition at line 300 of file GNuMIFlux.h.

301  { fMaxEFudge = fudge; }
double fMaxEFudge
fudge factor for estmating max enu (0=> use fixed 120GeV)
Definition: GNuMIFlux.h:404
void GNuMIFlux::SetMaxEnergy ( double  Ev)

specify maximum flx neutrino energy

Definition at line 760 of file GNuMIFlux.cxx.

761 {
762  fMaxEv = TMath::Max(0.,Ev);
763 
764  LOG("Flux", pINFO)
765  << "Declared maximum flux neutrino energy: " << fMaxEv;
766 }
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
#define pINFO
Definition: Messenger.h:62
double fMaxEv
maximum energy
Definition: GNuMIFlux.h:384
void genie::flux::GNuMIFlux::SetMaxWgtScan ( double  fudge = 1.05,
long int  nentries = 2500000 
)
inline
Parameters
nentriesconfiguration when estimating max weight

Definition at line 298 of file GNuMIFlux.h.

299  { fMaxWgtFudge = fudge; fMaxWgtEntries = nentries; }
double fMaxWgtFudge
fudge factor for estimating max wgt
Definition: GNuMIFlux.h:402
long int fMaxWgtEntries
of entries in estimating max wgt
Definition: GNuMIFlux.h:403
void GNuMIFlux::SetTreeName ( string  name)

set input tree name (default: "h10")

Definition at line 776 of file GNuMIFlux.cxx.

777 {
779 }
static QCString name
Definition: declinfo.cpp:673
string fNuFluxTreeName
Tree name "h10" (g3) or "nudata" (g4)
Definition: GNuMIFlux.h:388
double genie::flux::GNuMIFlux::SumWeight ( void  ) const
inline

integrated weight for flux neutrinos looped so far

Definition at line 267 of file GNuMIFlux.h.

double GNuMIFlux::UsedPOTs ( void  ) const

of protons-on-target used

Definition at line 464 of file GNuMIFlux.cxx.

465 {
466 // Compute current number of flux POTs
467 
468  if (!fNuFluxTree) {
469  LOG("Flux", pWARN)
470  << "The flux driver has not been properly configured";
471  return 0;
472  }
473  return fAccumPOTs;
474 }
TChain * fNuFluxTree
TTree in g3numi or g4numi // REF ONLY!
Definition: GNuMIFlux.h:389
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
double fAccumPOTs
POTs used so far.
Definition: GNuMIFlux.h:411
#define pWARN
Definition: Messenger.h:60
void GNuMIFlux::UseFluxAtFarDetCenter ( void  )

Definition at line 791 of file GNuMIFlux.cxx.

792 {
795  fFluxWindowDir1 = TLorentzVector(); // no extent
796  fFluxWindowDir2 = TLorentzVector();
797  fUseFluxAtDetCenter = +1;
798  fDetLocIsSet = true;
799 }
void SetLengthUnits(double user_units)
Set units assumed by user.
Definition: GNuMIFlux.cxx:1206
TLorentzVector fFluxWindowDir2
extent for flux window (direction 2)
Definition: GNuMIFlux.h:429
TLorentzVector fFluxWindowDir1
extent for flux window (direction 1)
Definition: GNuMIFlux.h:428
double UnitFromString(string u)
Definition: UnitUtils.cxx:18
bool fDetLocIsSet
is a flux location (near/far) set?
Definition: GNuMIFlux.h:415
const TLorentzVector kPosCenterFarBeam(0., 0., 735340.00, 0.)
int fUseFluxAtDetCenter
use flux at near (-1) or far (+1) det center from ntuple?
Definition: GNuMIFlux.h:416
TLorentzVector fFluxWindowBase
base point for flux window - beam coord
Definition: GNuMIFlux.h:427
void GNuMIFlux::UseFluxAtNearDetCenter ( void  )

force weights at MINOS detector "center" as found in ntuple

Definition at line 781 of file GNuMIFlux.cxx.

782 {
784  fFluxWindowBase = kPosCenterNearBeam;
785  fFluxWindowDir1 = TLorentzVector(); // no extent
786  fFluxWindowDir2 = TLorentzVector();
787  fUseFluxAtDetCenter = -1;
788  fDetLocIsSet = true;
789 }
void SetLengthUnits(double user_units)
Set units assumed by user.
Definition: GNuMIFlux.cxx:1206
TLorentzVector fFluxWindowDir2
extent for flux window (direction 2)
Definition: GNuMIFlux.h:429
TLorentzVector fFluxWindowDir1
extent for flux window (direction 1)
Definition: GNuMIFlux.h:428
double UnitFromString(string u)
Definition: UnitUtils.cxx:18
bool fDetLocIsSet
is a flux location (near/far) set?
Definition: GNuMIFlux.h:415
int fUseFluxAtDetCenter
use flux at near (-1) or far (+1) det center from ntuple?
Definition: GNuMIFlux.h:416
TLorentzVector fFluxWindowBase
base point for flux window - beam coord
Definition: GNuMIFlux.h:427
void GNuMIFlux::User2BeamDir ( const TLorentzVector &  usrdir,
TLorentzVector &  beamdir 
) const

Definition at line 992 of file GNuMIFlux.cxx.

994 {
995  beamdir = fLengthScaleU2B*(fBeamRotInv*usrdir);
996 }
double fLengthScaleU2B
scale factor beam user –> (cm)
Definition: GNuMIFlux.h:420
TLorentzRotation fBeamRotInv
Definition: GNuMIFlux.h:424
void GNuMIFlux::User2BeamP4 ( const TLorentzVector &  usrp4,
TLorentzVector &  beamp4 
) const

Definition at line 997 of file GNuMIFlux.cxx.

999 {
1000  beamp4 = fBeamRotInv*usrp4;
1001 }
TLorentzRotation fBeamRotInv
Definition: GNuMIFlux.h:424
void GNuMIFlux::User2BeamPos ( const TLorentzVector &  usrxyz,
TLorentzVector &  beamxyz 
) const

Definition at line 987 of file GNuMIFlux.cxx.

989 {
990  beamxyz = fLengthScaleU2B*(fBeamRotInv*(usrxyz-fBeamZero));
991 }
TLorentzVector fBeamZero
beam origin in user coords
Definition: GNuMIFlux.h:422
double fLengthScaleU2B
scale factor beam user –> (cm)
Definition: GNuMIFlux.h:420
TLorentzRotation fBeamRotInv
Definition: GNuMIFlux.h:424
double genie::flux::GNuMIFlux::Weight ( void  )
inlinevirtual

returns the flux neutrino weight (if any)

Implements genie::GFluxI.

Definition at line 234 of file GNuMIFlux.h.

234 { return fWeight; }
double fWeight
current neutrino weight, =1 if generating unweighted entries
Definition: GNuMIFlux.h:400

Member Data Documentation

double genie::flux::GNuMIFlux::fAccumPOTs
private

POTs used so far.

Definition at line 411 of file GNuMIFlux.h.

bool genie::flux::GNuMIFlux::fApplyTiltWeight
private

wgt due to window normal not || beam

Definition at line 414 of file GNuMIFlux.h.

TLorentzRotation genie::flux::GNuMIFlux::fBeamRot
private

rotation applied beam –> user coord

Definition at line 423 of file GNuMIFlux.h.

TLorentzRotation genie::flux::GNuMIFlux::fBeamRotInv
private

Definition at line 424 of file GNuMIFlux.h.

TLorentzVector genie::flux::GNuMIFlux::fBeamZero
private

beam origin in user coords

Definition at line 422 of file GNuMIFlux.h.

GNuMIFluxPassThroughInfo* genie::flux::GNuMIFlux::fCurEntry
private

copy of current ntuple entry info (owned structure)

Definition at line 436 of file GNuMIFlux.h.

bool genie::flux::GNuMIFlux::fDetLocIsSet
private

is a flux location (near/far) set?

Definition at line 415 of file GNuMIFlux.h.

double genie::flux::GNuMIFlux::fEffPOTsPerNu
private

what a entry is worth ...

Definition at line 410 of file GNuMIFlux.h.

bool genie::flux::GNuMIFlux::fEnd
private

end condition reached

Definition at line 385 of file GNuMIFlux.h.

Long64_t genie::flux::GNuMIFlux::fFilePOTs
private

of protons-on-target represented by all files

Definition at line 398 of file GNuMIFlux.h.

flugg* genie::flux::GNuMIFlux::fFlugg
private

flugg ntuple

Definition at line 393 of file GNuMIFlux.h.

TLorentzVector genie::flux::GNuMIFlux::fFluxWindowBase
private

base point for flux window - beam coord

Definition at line 427 of file GNuMIFlux.h.

TLorentzVector genie::flux::GNuMIFlux::fFluxWindowDir1
private

extent for flux window (direction 1)

Definition at line 428 of file GNuMIFlux.h.

TLorentzVector genie::flux::GNuMIFlux::fFluxWindowDir2
private

extent for flux window (direction 2)

Definition at line 429 of file GNuMIFlux.h.

double genie::flux::GNuMIFlux::fFluxWindowLen1
private

Definition at line 430 of file GNuMIFlux.h.

double genie::flux::GNuMIFlux::fFluxWindowLen2
private

Definition at line 431 of file GNuMIFlux.h.

TVector3 genie::flux::GNuMIFlux::fFluxWindowPtUser[3]
private

user points of flux window

Definition at line 426 of file GNuMIFlux.h.

g3numi* genie::flux::GNuMIFlux::fG3NuMI
private

g3numi ntuple

Definition at line 391 of file GNuMIFlux.h.

g4numi* genie::flux::GNuMIFlux::fG4NuMI
private

g4numi ntuple

Definition at line 392 of file GNuMIFlux.h.

bool genie::flux::GNuMIFlux::fGenWeighted
private

does GenerateNext() give weights?

Definition at line 413 of file GNuMIFlux.h.

TLorentzVector genie::flux::GNuMIFlux::fgX4dkvtx
private

decay 4-position beam coord

Definition at line 434 of file GNuMIFlux.h.

Long64_t genie::flux::GNuMIFlux::fIEntry
private

current flux ntuple entry

Definition at line 396 of file GNuMIFlux.h.

long int genie::flux::GNuMIFlux::fIUse
private

current # of times an entry has been used

Definition at line 407 of file GNuMIFlux.h.

double genie::flux::GNuMIFlux::fLengthScaleB2U
private

scale factor beam (cm) –> user

Definition at line 419 of file GNuMIFlux.h.

double genie::flux::GNuMIFlux::fLengthScaleU2B
private

scale factor beam user –> (cm)

Definition at line 420 of file GNuMIFlux.h.

double genie::flux::GNuMIFlux::fLengthUnits
private

units for coord in user exchanges

Definition at line 418 of file GNuMIFlux.h.

double genie::flux::GNuMIFlux::fMaxEFudge
private

fudge factor for estmating max enu (0=> use fixed 120GeV)

Definition at line 404 of file GNuMIFlux.h.

double genie::flux::GNuMIFlux::fMaxEv
private

maximum energy

Definition at line 384 of file GNuMIFlux.h.

double genie::flux::GNuMIFlux::fMaxWeight
private

max flux neutrino weight in input file

Definition at line 401 of file GNuMIFlux.h.

long int genie::flux::GNuMIFlux::fMaxWgtEntries
private

of entries in estimating max wgt

Definition at line 403 of file GNuMIFlux.h.

double genie::flux::GNuMIFlux::fMaxWgtFudge
private

fudge factor for estimating max wgt

Definition at line 402 of file GNuMIFlux.h.

Long64_t genie::flux::GNuMIFlux::fNEntries
private

number of flux ntuple entries

Definition at line 395 of file GNuMIFlux.h.

int genie::flux::GNuMIFlux::fNFiles
private

number of files in chain

Definition at line 394 of file GNuMIFlux.h.

long int genie::flux::GNuMIFlux::fNNeutrinos
private

number of flux neutrinos thrown so far

Definition at line 409 of file GNuMIFlux.h.

std::vector<string> genie::flux::GNuMIFlux::fNuFluxFilePatterns
private

(potentially wildcarded) path(s)

Definition at line 387 of file GNuMIFlux.h.

string genie::flux::GNuMIFlux::fNuFluxGen
private

"g3numi" "g4numi" or "flugg"

Definition at line 390 of file GNuMIFlux.h.

TChain* genie::flux::GNuMIFlux::fNuFluxTree
private

TTree in g3numi or g4numi // REF ONLY!

Definition at line 389 of file GNuMIFlux.h.

string genie::flux::GNuMIFlux::fNuFluxTreeName
private

Tree name "h10" (g3) or "nudata" (g4)

Definition at line 388 of file GNuMIFlux.h.

long int genie::flux::GNuMIFlux::fNUse
private

how often to use same entry in a row

Definition at line 406 of file GNuMIFlux.h.

Long64_t genie::flux::GNuMIFlux::fNuTot
private

cummulative # of entries (=fNEntries)

Definition at line 397 of file GNuMIFlux.h.

double genie::flux::GNuMIFlux::fSumWeight
private

sum of weights for nus thrown so far

Definition at line 408 of file GNuMIFlux.h.

int genie::flux::GNuMIFlux::fUseFluxAtDetCenter
private

use flux at near (-1) or far (+1) det center from ntuple?

Definition at line 416 of file GNuMIFlux.h.

double genie::flux::GNuMIFlux::fWeight
private

current neutrino weight, =1 if generating unweighted entries

Definition at line 400 of file GNuMIFlux.h.

TVector3 genie::flux::GNuMIFlux::fWindowNormal
private

normal direction for flux window

Definition at line 432 of file GNuMIFlux.h.


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