A generic GENIE flux driver. Generates a 'cylindrical' neutrino beam along the input direction, with the input transverse radius and centered at the input position. The energies are generated from the input energy spectrum (TH1D). Multiple neutrino species can be generated (you will need to supply an energy spectrum for each). More...

#include <GCylindTH1Flux.h>

Inheritance diagram for genie::flux::GCylindTH1Flux:

[legend]

Collaboration diagram for genie::flux::GCylindTH1Flux:

[legend]

Public Member Functions
	GCylindTH1Flux ()
	~GCylindTH1Flux ()
void	SetNuDirection (const TVector3 &direction)
void	SetBeamSpot (const TVector3 &spot)
void	SetTransverseRadius (double Rt)
void	AddEnergySpectrum (int nu_pdgc, TH1D *spectrum)
void	SetRtDependence (string rdep)
const PDGCodeList &	FluxParticles (void)
	declare list of flux neutrinos that can be generated (for init. purposes)
double	MaxEnergy (void)
	declare the max flux neutrino energy that can be generated (for init. purposes)
bool	GenerateNext (void)
	generate the next flux neutrino (return false in err)
int	PdgCode (void)
	returns the flux neutrino pdg code
double	Weight (void)
	returns the flux neutrino weight (if any)
const TLorentzVector &	Momentum (void)
	returns the flux neutrino 4-momentum
const TLorentzVector &	Position (void)
	returns the flux neutrino 4-position (note: expect SI rather than physical units)
bool	End (void)
	true if no more flux nu's can be thrown (eg reaching end of beam sim ntuples)
long int	Index (void)
	returns corresponding index for current flux neutrino (e.g. for a flux ntuple returns the current entry number)
void	Clear (Option_t *opt)
	reset state variables based on opt
void	GenerateWeighted (bool gen_weighted)
	set whether to generate weighted or unweighted neutrinos
Public Member Functions inherited from genie::GFluxI
virtual	~GFluxI ()

Private Member Functions
void	Initialize (void)
void	CleanUp (void)
void	ResetSelection (void)
void	AddAllFluxes (void)
int	SelectNeutrino (double Ev)
double	GeneratePhi (void) const
double	GenerateRt (void) const

Private Attributes
double	fMaxEv
	maximum energy
PDGCodeList *	fPdgCList
	list of neutrino pdg-codes
int	fgPdgC
	running generated nu pdg-code
TLorentzVector	fgP4
	running generated nu 4-momentum
TLorentzVector	fgX4
	running generated nu 4-position
vector< TH1D * >	fSpectrum
	flux = f(Ev), 1/neutrino species
TH1D *	fTotSpectrum
	combined flux = f(Ev)
TVector3 *	fDirVec
	neutrino direction
TVector3 *	fBeamSpot
	beam spot position
double	fRt
	transverse size of neutrino beam
TF1 *	fRtDep
	transverse radius dependence

Additional Inherited Members
Protected Member Functions inherited from genie::GFluxI
	GFluxI ()

Detailed Description

A generic GENIE flux driver. Generates a 'cylindrical' neutrino beam along the input direction, with the input transverse radius and centered at the input position. The energies are generated from the input energy spectrum (TH1D). Multiple neutrino species can be generated (you will need to supply an energy spectrum for each).

Author: Costas Andreopoulos <c.andreopoulos \at cern.ch> University of Liverpool

Created:\n July 4, 2005

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

Definition at line 43 of file GCylindTH1Flux.h.

Constructor & Destructor Documentation

◆ GCylindTH1Flux()

GCylindTH1Flux::GCylindTH1Flux ( )

Definition at line 35 of file GCylindTH1Flux.cxx.

{
  this->Initialize();
}

References Initialize().

◆ ~GCylindTH1Flux()

GCylindTH1Flux::~GCylindTH1Flux ( )

Definition at line 40 of file GCylindTH1Flux.cxx.

{
  this->CleanUp();
}

References CleanUp().

Member Function Documentation

◆ AddAllFluxes()

void GCylindTH1Flux::AddAllFluxes ( void )

private

Definition at line 218 of file GCylindTH1Flux.cxx.

{
  LOG("Flux", pNOTICE) << "Computing combined flux";
 
  if(fTotSpectrum) delete fTotSpectrum;
 
  vector<TH1D *>::const_iterator spectrum_iter;
 
  unsigned int inu=0;
  for(spectrum_iter = fSpectrum.begin();
                       spectrum_iter != fSpectrum.end(); ++spectrum_iter) {
     TH1D * spectrum = *spectrum_iter;
 
     if(inu==0) { fTotSpectrum = new TH1D(*spectrum); }
     else       { fTotSpectrum->Add(spectrum);        }
     inu++;
  }
}

References fSpectrum, fTotSpectrum, LOG, and pNOTICE.

Referenced by AddEnergySpectrum().

◆ AddEnergySpectrum()

void GCylindTH1Flux::AddEnergySpectrum	(	int	nu_pdgc,
		TH1D *	spectrum )

Definition at line 183 of file GCylindTH1Flux.cxx.

{
  LOG("Flux", pNOTICE) << "Adding flux spectrum for pdg = " << nu_pdgc;
 
  fPdgCList->push_back(nu_pdgc);
 
  bool accepted = (count(fPdgCList->begin(),fPdgCList->end(),nu_pdgc) == 1);
  if(!accepted) {
     LOG ("Flux", pWARN)
            << "The pdg-code isn't recognized and the spectrum was ignored";
  } else {
     fSpectrum.push_back(spectrum);
 
     int    nb  = spectrum->GetNbinsX();
     Axis_t max = spectrum->GetBinLowEdge(nb)+spectrum->GetBinWidth(nb);
     fMaxEv = TMath::Max(fMaxEv, (double)max);
 
     LOG("Flux", pNOTICE)
          << "Updating maximum energy of flux particles to: " << fMaxEv;
 
     this->AddAllFluxes(); // update combined flux
  }
}

References AddAllFluxes(), fMaxEv, fPdgCList, fSpectrum, LOG, pNOTICE, and pWARN.

Referenced by main().

◆ CleanUp()

void GCylindTH1Flux::CleanUp ( void )

private

Definition at line 145 of file GCylindTH1Flux.cxx.

{
  LOG("Flux", pNOTICE) << "Cleaning up...";
 
  if (fDirVec     ) delete fDirVec;
  if (fBeamSpot   ) delete fBeamSpot;
  if (fPdgCList   ) delete fPdgCList;
  if (fTotSpectrum) delete fTotSpectrum;
  if (fRtDep      ) delete fRtDep;
 
  unsigned int nspectra = fSpectrum.size();
  for(unsigned int i = 0; i < nspectra; i++) {
     TH1D * spectrum = fSpectrum[i];
     delete spectrum;
     spectrum = 0;
  }
}

References fBeamSpot, fDirVec, fPdgCList, fRtDep, fSpectrum, fTotSpectrum, LOG, and pNOTICE.

Referenced by ~GCylindTH1Flux().

◆ Clear()

void GCylindTH1Flux::Clear ( Option_t * opt )

virtual

reset state variables based on opt

Implements genie::GFluxI.

Definition at line 103 of file GCylindTH1Flux.cxx.

{
// Dummy clear method needed to conform to GFluxI interface
//
  LOG("Flux", pERROR) <<
      "No Clear(Option_t * opt) method implemented for opt: "<< opt;
}

References LOG, and pERROR.

◆ End()

bool genie::flux::GCylindTH1Flux::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 64 of file GCylindTH1Flux.h.

64{ return false; }

◆ FluxParticles()

const PDGCodeList & genie::flux::GCylindTH1Flux::FluxParticles ( void )

inlinevirtual

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

Implements genie::GFluxI.

Definition at line 57 of file GCylindTH1Flux.h.

57{ return *fPdgCList; }

References fPdgCList.

◆ GenerateNext()

bool GCylindTH1Flux::GenerateNext ( void )

virtual

generate the next flux neutrino (return false in err)

Implements genie::GFluxI.

Definition at line 45 of file GCylindTH1Flux.cxx.

{
  //-- Reset previously generated neutrino code / 4-p / 4-x
  this->ResetSelection();
 
  //-- Generate an energy from the 'combined' spectrum histogram
  //   and compute the momentum vector
  double Ev = (double) fTotSpectrum->GetRandom();
 
  TVector3 p3(*fDirVec); // momentum along the neutrino direction
  p3.SetMag(Ev);         // with |p|=Ev
  // EDIT: Check if we're running with DM beam
  if (fPdgCList->ExistsInPDGCodeList(kPdgDarkMatter) || fPdgCList->ExistsInPDGCodeList(kPdgAntiDarkMatter)) {
    double Md = PDGLibrary::Instance()->Find(kPdgDarkMatter)->Mass();
    p3.SetMag(TMath::Sqrt(Ev*Ev - Md*Md));
  }
 
  fgP4.SetPxPyPzE(p3.Px(), p3.Py(), p3.Pz(), Ev);
 
  //-- Select a neutrino species from the flux fractions at the
  //   selected energy
  fgPdgC = (*fPdgCList)[this->SelectNeutrino(Ev)];
 
  //-- Compute neutrino 4-x
 
  if(fRt <= 0) {
    fgX4.SetXYZT(0.,0.,0.,0.);
  }
  else {
    // Create a vector (vec) that points to a random position at a disk
    // of radius Rt passing through the origin, perpendicular to the
    // input direction.
    TVector3 vec0(*fDirVec);           // vector along the input direction
    TVector3 vec = vec0.Orthogonal();  // orthogonal vector
 
    double psi = this->GeneratePhi();  // rndm angle [0,2pi]
    double Rt  = this->GenerateRt();   // rndm R [0,Rtransverse]
 
    vec.Rotate(psi,vec0); // rotate around original vector
    vec.SetMag(Rt);       // set new norm
 
    // Set the neutrino position as beam_spot + vec
    double x = fBeamSpot->X() + vec.X();
    double y = fBeamSpot->Y() + vec.Y();
    double z = fBeamSpot->Z() + vec.Z();
 
    fgX4.SetXYZT(x,y,z,0.);
  }
 
  LOG("Flux", pINFO) << "Generated neutrino pdg-code: " << fgPdgC;
  LOG("Flux", pINFO)
        << "Generated neutrino p4: " << utils::print::P4AsShortString(&fgP4);
  LOG("Flux", pINFO)
             << "Generated neutrino x4: " << utils::print::X4AsString(&fgX4);
 
  return true;
}

References fBeamSpot, fDirVec, fgP4, fgPdgC, fgX4, genie::PDGLibrary::Find(), fPdgCList, fRt, fTotSpectrum, GeneratePhi(), GenerateRt(), genie::PDGLibrary::Instance(), genie::kPdgAntiDarkMatter, genie::kPdgDarkMatter, LOG, genie::utils::print::P4AsShortString(), pINFO, ResetSelection(), SelectNeutrino(), and genie::utils::print::X4AsString().

◆ GeneratePhi()

double GCylindTH1Flux::GeneratePhi ( void ) const

private

Definition at line 271 of file GCylindTH1Flux.cxx.

{
  RandomGen * rnd = RandomGen::Instance();
  double phi = 2.*kPi * rnd->RndFlux().Rndm(); // [0,2pi]
  return phi;
}

References genie::RandomGen::Instance(), genie::constants::kPi, and genie::RandomGen::RndFlux().

Referenced by GenerateNext().

◆ GenerateRt()

double GCylindTH1Flux::GenerateRt ( void ) const

private

Definition at line 278 of file GCylindTH1Flux.cxx.

{
  double Rt = fRtDep->GetRandom(); // rndm R [0,Rtransverse]
  return Rt;
}

References fRtDep.

Referenced by GenerateNext().

◆ GenerateWeighted()

void GCylindTH1Flux::GenerateWeighted ( bool gen_weighted )

virtual

set whether to generate weighted or unweighted neutrinos

Implements genie::GFluxI.

Definition at line 111 of file GCylindTH1Flux.cxx.

{
// Dummy implementation needed to conform to GFluxI interface
//
  LOG("Flux", pERROR) <<
      "No GenerateWeighted(bool gen_weighted) method implemented for " <<
      "gen_weighted: " << gen_weighted;
}

References LOG, and pERROR.

◆ Index()

long int genie::flux::GCylindTH1Flux::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 65 of file GCylindTH1Flux.h.

65{ return -1; }

◆ Initialize()

void GCylindTH1Flux::Initialize ( void )

private

Definition at line 120 of file GCylindTH1Flux.cxx.

{
  LOG("Flux", pNOTICE) << "Initializing GCylindTH1Flux driver";
 
  fMaxEv       = 0;
  fPdgCList    = new PDGCodeList;
  fTotSpectrum = 0;
  fDirVec      = 0;
  fBeamSpot    = 0;
  fRt          =-1;
  fRtDep       = 0;
 
  this->ResetSelection();
  this->SetRtDependence("x");
  //eg, other example: this->SetRtDependence("pow(x,2)");
}

References fBeamSpot, fDirVec, fMaxEv, fPdgCList, fRt, fRtDep, fTotSpectrum, LOG, pNOTICE, ResetSelection(), and SetRtDependence().

Referenced by GCylindTH1Flux().

◆ MaxEnergy()

double genie::flux::GCylindTH1Flux::MaxEnergy ( void )

inlinevirtual

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

Implements genie::GFluxI.

Definition at line 58 of file GCylindTH1Flux.h.

58{ return fMaxEv; }

References fMaxEv.

◆ Momentum()

const TLorentzVector & genie::flux::GCylindTH1Flux::Momentum ( void )

inlinevirtual

returns the flux neutrino 4-momentum

Implements genie::GFluxI.

Definition at line 62 of file GCylindTH1Flux.h.

62{ return fgP4; }

References fgP4.

◆ PdgCode()

int genie::flux::GCylindTH1Flux::PdgCode ( void )

inlinevirtual

returns the flux neutrino pdg code

Implements genie::GFluxI.

Definition at line 60 of file GCylindTH1Flux.h.

60{ return fgPdgC; }

References fgPdgC.

◆ Position()

const TLorentzVector & genie::flux::GCylindTH1Flux::Position ( void )

inlinevirtual

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

Implements genie::GFluxI.

Definition at line 63 of file GCylindTH1Flux.h.

63{ return fgX4; }

References fgX4.

◆ ResetSelection()

void GCylindTH1Flux::ResetSelection ( void )

private

Definition at line 137 of file GCylindTH1Flux.cxx.

{
// initializing running neutrino pdg-code, 4-position, 4-momentum
  fgPdgC = 0;
  fgP4.SetPxPyPzE (0.,0.,0.,0.);
  fgX4.SetXYZT    (0.,0.,0.,0.);
}

References fgP4, fgPdgC, and fgX4.

Referenced by GenerateNext(), and Initialize().

◆ SelectNeutrino()

int GCylindTH1Flux::SelectNeutrino ( double Ev )

private

Definition at line 237 of file GCylindTH1Flux.cxx.

{
  const unsigned int n = fPdgCList->size();
  double fraction[n];
 
  vector<TH1D *>::const_iterator spectrum_iter;
 
  unsigned int inu=0;
  for(spectrum_iter = fSpectrum.begin();
                       spectrum_iter != fSpectrum.end(); ++spectrum_iter) {
     TH1D * spectrum = *spectrum_iter;
     fraction[inu++] = spectrum->GetBinContent(spectrum->FindBin(Ev));
  }
 
  double sum = 0;
  for(inu = 0; inu < n; inu++) {
     sum += fraction[inu];
     fraction[inu] = sum;
     LOG("Flux", pDEBUG) << "SUM-FRACTION(0->" << inu <<") = " << sum;
  }
 
  RandomGen * rnd = RandomGen::Instance();
  double R = sum * rnd->RndFlux().Rndm();
 
  LOG("Flux", pDEBUG) << "R e [0,SUM] = " << R;
 
  for(inu = 0; inu < n; inu++) {if ( R < fraction[inu] ) return inu;}
 
  LOG("Flux", pERROR) << "Could not select a neutrino species";
  assert(false);
 
  return -1;
}

References fPdgCList, fSpectrum, genie::RandomGen::Instance(), LOG, pDEBUG, pERROR, and genie::RandomGen::RndFlux().

Referenced by GenerateNext().

◆ SetBeamSpot()

void GCylindTH1Flux::SetBeamSpot ( const TVector3 & spot )

Definition at line 169 of file GCylindTH1Flux.cxx.

{
  if(fBeamSpot) delete fBeamSpot;
  fBeamSpot = new TVector3(spot);
}

References fBeamSpot.

Referenced by main().

◆ SetNuDirection()

void GCylindTH1Flux::SetNuDirection ( const TVector3 & direction )

Definition at line 163 of file GCylindTH1Flux.cxx.

{
  if(fDirVec) delete fDirVec;
  fDirVec = new TVector3(direction);
}

References fDirVec.

Referenced by main().

◆ SetRtDependence()

void GCylindTH1Flux::SetRtDependence ( string rdep )

Definition at line 207 of file GCylindTH1Flux.cxx.

{
// Set the (functional form of) Rt dependence as string, eg "x*x+sin(x)"
// You do not need to set this method. The default behaviour is to generate
// flux neutrinos uniformly over the area of the cylinder's cross section.
 
  if(fRtDep) delete fRtDep;
 
  fRtDep = new TF1("rdep", rdep.c_str(), 0,fRt);
}

References fRt, and fRtDep.

Referenced by Initialize().

◆ SetTransverseRadius()

void GCylindTH1Flux::SetTransverseRadius ( double Rt )

Definition at line 175 of file GCylindTH1Flux.cxx.

{
  LOG ("Flux", pNOTICE) << "Setting R[transverse] = " << Rt;
  fRt = Rt;
 
  if(fRtDep) fRtDep->SetRange(0,Rt);
}

References fRt, fRtDep, LOG, and pNOTICE.

Referenced by main().

◆ Weight()

double genie::flux::GCylindTH1Flux::Weight ( void )

inlinevirtual

returns the flux neutrino weight (if any)

Implements genie::GFluxI.

Definition at line 61 of file GCylindTH1Flux.h.

61{ return 1.0; }

Member Data Documentation

◆ fBeamSpot

TVector3* genie::flux::GCylindTH1Flux::fBeamSpot

private

beam spot position

Definition at line 89 of file GCylindTH1Flux.h.

Referenced by CleanUp(), GenerateNext(), Initialize(), and SetBeamSpot().

◆ fDirVec

TVector3* genie::flux::GCylindTH1Flux::fDirVec

private

neutrino direction

Definition at line 88 of file GCylindTH1Flux.h.

Referenced by CleanUp(), GenerateNext(), Initialize(), and SetNuDirection().

◆ fgP4

TLorentzVector genie::flux::GCylindTH1Flux::fgP4

private

running generated nu 4-momentum

Definition at line 84 of file GCylindTH1Flux.h.

Referenced by GenerateNext(), Momentum(), and ResetSelection().

◆ fgPdgC

int genie::flux::GCylindTH1Flux::fgPdgC

private

running generated nu pdg-code

Definition at line 83 of file GCylindTH1Flux.h.

Referenced by GenerateNext(), PdgCode(), and ResetSelection().

◆ fgX4

TLorentzVector genie::flux::GCylindTH1Flux::fgX4

private

running generated nu 4-position

Definition at line 85 of file GCylindTH1Flux.h.

Referenced by GenerateNext(), Position(), and ResetSelection().

◆ fMaxEv

double genie::flux::GCylindTH1Flux::fMaxEv

private

maximum energy

Definition at line 81 of file GCylindTH1Flux.h.

Referenced by AddEnergySpectrum(), Initialize(), and MaxEnergy().

◆ fPdgCList

PDGCodeList* genie::flux::GCylindTH1Flux::fPdgCList

private

list of neutrino pdg-codes

Definition at line 82 of file GCylindTH1Flux.h.

Referenced by AddEnergySpectrum(), CleanUp(), FluxParticles(), GenerateNext(), Initialize(), and SelectNeutrino().

◆ fRt

double genie::flux::GCylindTH1Flux::fRt

private

transverse size of neutrino beam

Definition at line 90 of file GCylindTH1Flux.h.

Referenced by GenerateNext(), Initialize(), SetRtDependence(), and SetTransverseRadius().

◆ fRtDep

TF1* genie::flux::GCylindTH1Flux::fRtDep

private

transverse radius dependence

Definition at line 91 of file GCylindTH1Flux.h.

Referenced by CleanUp(), GenerateRt(), Initialize(), SetRtDependence(), and SetTransverseRadius().

◆ fSpectrum

vector<TH1D *> genie::flux::GCylindTH1Flux::fSpectrum

private

flux = f(Ev), 1/neutrino species

Definition at line 86 of file GCylindTH1Flux.h.

Referenced by AddAllFluxes(), AddEnergySpectrum(), CleanUp(), and SelectNeutrino().

◆ fTotSpectrum

TH1D* genie::flux::GCylindTH1Flux::fTotSpectrum

private

combined flux = f(Ev)

Definition at line 87 of file GCylindTH1Flux.h.

Referenced by AddAllFluxes(), CleanUp(), GenerateNext(), and Initialize().

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

Public Member Functions

Private Member Functions

Private Attributes

Additional Inherited Members

Detailed Description

Constructor & Destructor Documentation

◆ GCylindTH1Flux()

◆ ~GCylindTH1Flux()

Member Function Documentation

◆ AddAllFluxes()

◆ AddEnergySpectrum()

◆ CleanUp()

◆ Clear()

◆ End()

◆ FluxParticles()

◆ GenerateNext()

◆ GeneratePhi()

◆ GenerateRt()

◆ GenerateWeighted()

◆ Index()

◆ Initialize()

◆ MaxEnergy()

◆ Momentum()

◆ PdgCode()

◆ Position()

◆ ResetSelection()

◆ SelectNeutrino()

◆ SetBeamSpot()

◆ SetNuDirection()

◆ SetRtDependence()

◆ SetTransverseRadius()

◆ Weight()

Member Data Documentation

◆ fBeamSpot

◆ fDirVec

◆ fgP4

◆ fgPdgC

◆ fgX4

◆ fMaxEv

◆ fPdgCList

◆ fRt

◆ fRtDep

◆ fSpectrum

◆ fTotSpectrum