genie::flux::GCylindTH1Flux Class Reference

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:

Collaboration diagram for genie::flux::GCylindTH1Flux:

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

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 cern.ch> University of Liverpool

Created:

July 4, 2005

License:

Copyright (c) 2003-2024, 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

(

)

Definition at line 35 of file GCylindTH1Flux.cxx.

References Initialize().

{
  this->Initialize();
}

GCylindTH1Flux::~GCylindTH1Flux

(

)

Definition at line 40 of file GCylindTH1Flux.cxx.

References CleanUp().

{
  this->CleanUp();
}

Member Function Documentation

void GCylindTH1Flux::AddAllFluxes ( void  )

private

Definition at line 218 of file GCylindTH1Flux.cxx.

References fSpectrum, fTotSpectrum, LOG, and pNOTICE.

Referenced by AddEnergySpectrum().

{
  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++;
  }
}

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

Definition at line 183 of file GCylindTH1Flux.cxx.

References AddAllFluxes(), fMaxEv, fPdgCList, fSpectrum, LOG, genie::units::nb, pNOTICE, genie::PDGCodeList::push_back(), and pWARN.

Referenced by main().

{
  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
  }
}

void GCylindTH1Flux::CleanUp ( void  )

private

Definition at line 145 of file GCylindTH1Flux.cxx.

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

Referenced by ~GCylindTH1Flux().

{
  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;
  }
}

void GCylindTH1Flux::Clear ( Option_t *  opt )

virtual

reset state variables based on opt

Implements genie::GFluxI.

Definition at line 103 of file GCylindTH1Flux.cxx.

References LOG, and pERROR.

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

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.

{ return  false;     }

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.

References fPdgCList.

{ return *fPdgCList; }

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.

References genie::PDGCodeList::ExistsInPDGCodeList(), 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().

{
  //-- 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;
}

double GCylindTH1Flux::GeneratePhi ( void  ) const

private

Definition at line 271 of file GCylindTH1Flux.cxx.

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

Referenced by GenerateNext().

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

double GCylindTH1Flux::GenerateRt ( void  ) const

private

Definition at line 278 of file GCylindTH1Flux.cxx.

References fRtDep.

Referenced by GenerateNext().

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

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.

References LOG, and pERROR.

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

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.

{ return -1;         }

void GCylindTH1Flux::Initialize ( void  )

private

Definition at line 120 of file GCylindTH1Flux.cxx.

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

{
  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)");
}

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.

References fMaxEv.

{ return  fMaxEv;    }

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.

References fgP4.

{ return  fgP4;      }

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

inlinevirtual

returns the flux neutrino pdg code

Implements genie::GFluxI.

Definition at line 60 of file GCylindTH1Flux.h.

References fgPdgC.

{ return  fgPdgC;    }

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.

References fgX4.

{ return  fgX4;      }

void GCylindTH1Flux::ResetSelection ( void  )

private

Definition at line 137 of file GCylindTH1Flux.cxx.

References fgP4, fgPdgC, and fgX4.

Referenced by GenerateNext(), and Initialize().

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

int GCylindTH1Flux::SelectNeutrino ( double  Ev )

private

Definition at line 237 of file GCylindTH1Flux.cxx.

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

Referenced by GenerateNext().

{
  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;
}

void GCylindTH1Flux::SetBeamSpot

(

const TVector3 &

spot

)

Definition at line 169 of file GCylindTH1Flux.cxx.

References fBeamSpot.

Referenced by main().

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

void GCylindTH1Flux::SetNuDirection

(

const TVector3 &

direction

)

Definition at line 163 of file GCylindTH1Flux.cxx.

References fDirVec.

Referenced by main().

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

void GCylindTH1Flux::SetRtDependence

(

string

rdep

)

Definition at line 207 of file GCylindTH1Flux.cxx.

References fRt, and fRtDep.

Referenced by Initialize().

{
// 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);
}

void GCylindTH1Flux::SetTransverseRadius

(

double

Rt

)

Definition at line 175 of file GCylindTH1Flux.cxx.

References fRt, fRtDep, LOG, and pNOTICE.

Referenced by main().

{
  LOG ("Flux", pNOTICE) << "Setting R[transverse] = " << Rt;
  fRt = Rt;

  if(fRtDep) fRtDep->SetRange(0,Rt);
}

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.

{ return  1.0;       }

Member Data Documentation

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

private

beam spot position

Definition at line 89 of file GCylindTH1Flux.h.

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

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

private

neutrino direction

Definition at line 88 of file GCylindTH1Flux.h.

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

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().

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().

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().

double genie::flux::GCylindTH1Flux::fMaxEv

private

maximum energy

Definition at line 81 of file GCylindTH1Flux.h.

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

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().

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().

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().

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().

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:

• GCylindTH1Flux.h
• GCylindTH1Flux.cxx