/search.css" rel="stylesheet" type="text/css"/> /search.js">
| Classes | Job Modules | Data Objects | Services | Algorithms | Tools | Packages | Directories | Tracs |

In This Package:

Public Member Functions | Static Public Member Functions | Private Member Functions | Private Attributes
GtInverseBeta Class Reference

Convert old standalone generator into GenTool style class qinghe@princeton.edu Sep 17, 2009. More...

#include <GtInverseBeta.h>

Inheritance diagram for GtInverseBeta:
Inheritance graph
[legend]
Collaboration diagram for GtInverseBeta:
Collaboration graph
[legend]

List of all members.

Public Member Functions

 GtInverseBeta (const std::string &type, const std::string &name, const IInterface *parent)
virtual ~GtInverseBeta ()
virtual StatusCode initialize ()
virtual StatusCode finalize ()
virtual StatusCode mutate (HepMC::GenEvent &event)
void GenerateInverseBetaPrimaries (TF2 *, const double, Hep3Vector &, Hep3Vector &, Hep3Vector &)

Static Public Member Functions

static const InterfaceID & interfaceID ()

Private Member Functions

HepMC::GenParticlemake_particle (Hep3Vector &, double, int)

Private Attributes

TF2 * funcInvBetaProb
double m_neutrino_angle_in_deg
bool m_positron_only
bool m_neutron_only

Detailed Description

Convert old standalone generator into GenTool style class qinghe@princeton.edu Sep 17, 2009.

Definition at line 21 of file GtInverseBeta.h.


Constructor & Destructor Documentation

GtInverseBeta::GtInverseBeta ( const std::string &  type,
const std::string &  name,
const IInterface *  parent 
)

Definition at line 33 of file GtInverseBeta.cc.

  : GaudiTool(type,name,parent)
{
  declareInterface<IHepMCEventMutator>(this);
  
  declareProperty("NeutrinoAngle", m_neutrino_angle_in_deg=0, 
                  "neutrino angle in degree");
  declareProperty("PositronOnly",m_positron_only=0,"Only positron");
  declareProperty("NeutronOnly",m_neutron_only=0,"Only neutron");
}
GtInverseBeta::~GtInverseBeta ( ) [virtual]

Definition at line 46 of file GtInverseBeta.cc.

{
}

Member Function Documentation

StatusCode GtInverseBeta::initialize ( ) [virtual]

Definition at line 173 of file GtInverseBeta.cc.

{
  static int count = 0;
  count++;
  char itoa[100];
  sprintf(itoa,"%i",count);

  string filename("InvBeta");
  filename.append(itoa);  
  funcInvBetaProb = new TF2(filename.c_str(),inverse_beta_prob,
                            0.01,60.0,-1.,1., 0); //"0" parameters
  //set the bins on x and y to help generate pseudo random numbers
  //according to the reaction probabilities
  funcInvBetaProb->SetNpx(1000);//energy
  funcInvBetaProb->SetNpy(100);//costh

  return StatusCode::SUCCESS;
}
StatusCode GtInverseBeta::finalize ( ) [virtual]

Definition at line 192 of file GtInverseBeta.cc.

{
  delete funcInvBetaProb;

  return StatusCode::SUCCESS;
}
StatusCode GtInverseBeta::mutate ( HepMC::GenEvent event) [virtual]

Implements IHepMCEventMutator.

Definition at line 60 of file GtInverseBeta.cc.

{
  //  gReactorFlux = new KRLReactorFlux(); // defaults to Petr's fit
  //  gInverseBeta = new KRLInverseBeta();  
    
  Hep3Vector pNeutrino; //incoming neutrino momentum
  Hep3Vector pPositron, pNeutron; // the e+ and n momentum vectors. Unit GeV/c
  
  GenerateInverseBetaPrimaries(funcInvBetaProb, m_neutrino_angle_in_deg,
                               pNeutrino, 
                               pPositron, pNeutron);
  
  //the new GenVertex will be deleted by GenEvent destructor
  HepMC::GenVertex* vertex = new HepMC::GenVertex(HepMC::FourVector(0,0,0,0));
  event.add_vertex(vertex);
  event.set_signal_process_vertex(vertex);
  
  int pid = 0;
  HepMC::GenParticle* particle1 = 0;
  HepMC::GenParticle* particle2 = 0;
  if(!m_positron_only && !m_neutron_only){
    pid = -11;
    particle1 = make_particle(pPositron, electron_mass_c2, pid);
    vertex->add_particle_out(particle1);
    
    pid = 2112;
    particle2 = make_particle(pNeutron, neutron_mass_c2, pid);
    vertex->add_particle_out(particle2);
  }  
  else if(m_positron_only){
    pid = -11;
    particle1 = make_particle(pPositron, electron_mass_c2, pid);
    vertex->add_particle_out(particle1);
  }
  else if(m_neutron_only){
    pid = 2112;
    particle1 = make_particle(pPositron, neutron_mass_c2, pid);
    vertex->add_particle_out(particle1);
  }
  
  int mo_pid = -12;
  HepMC::FourVector fourmom(pNeutrino.x()/MeV,pNeutrino.y()/MeV,
                              pNeutrino.z()/MeV,pNeutrino.mag()/MeV);
  HepMC::GenParticle* particle_in = 
    new HepMC::GenParticle(fourmom, mo_pid, 2);
  // Save mother particle as info particle (default status=2)
  vertex->add_particle_in(particle_in);
  
  return StatusCode::SUCCESS;
}
void GtInverseBeta::GenerateInverseBetaPrimaries ( TF2 *  invBetaProb,
const double  nu_angle_wrt_y_in_deg,
Hep3Vector &  pnu,
Hep3Vector &  p1,
Hep3Vector &  p2 
)

Definition at line 111 of file GtInverseBeta.cc.

{
  //now create a TF2 object to calculate the reaction probability
  //given the neutrino energy and the inverse beta probability
  //x: eNu in MeV
  //y: cos(th_e)
//  TF2 *invBetaProb = new TF2("InvBeta",inverse_beta_prob,
//                               0.01,60.0,-1.,1., 0); //"0" parameters
//  //set the bins on x and y to help generate pseudo random numbers
//  //according to the reaction probabilities
//  invBetaProb->SetNpx(1000);//energy
//  invBetaProb->SetNpy(100);//costh

  //generate neutrino energy and outgoing electron cos angle
  // djaffe: 12may09 modified to avoid rare case when e+ energy
  //        is returned as identically zero which indicates that below-threshold
  //        production according to KRLInverseBeta.
  Double_t aEngNu(0), aCosThEl(0);
  Double_t Ee(0);
  while ( Ee == 0 ){
    invBetaProb->GetRandom2(aEngNu, aCosThEl);
    
    //cout<<aEngNu<<" : "<<acos(aCosThEl)*360./M_2PI<<endl;
    
    //at this point, assume the incident nu momentum along z
    
    //calculate positron energy in MeV
    //    Ee = gInverseBeta->Ee1(aEngNu, aCosThEl);
    KRLInverseBeta gInverseBeta;
    Ee = gInverseBeta.Ee1(aEngNu, aCosThEl);

  }
  //generate positron phi angle
  Double_t aPhiEl = gRandom->Uniform(0.,2.0*M_PI);
  Double_t kMe = electron_mass_c2/MeV;
  Double_t Pe = sqrt(Ee*Ee-kMe*kMe);
  p1.setRThetaPhi(Pe*MeV, acos(aCosThEl), aPhiEl);
  
  //now use the 3-momentum conservation to set neutron recoil 
  //momentum
  pnu.setRThetaPhi(aEngNu*MeV,0,0);
  p2 = pnu-p1;
  
  //now make the vector rotation. first rotate around x by 90 degrees
  //so that +z axis is now +y axis
  pnu.rotateX(M_PI/2.);
  p1.rotateX(M_PI/2.);
  p2.rotateX(M_PI/2.);
  
  //if nu_angle_wrt_y is non-zero, rotate again
  if(nu_angle_wrt_y_in_deg!=0) {
    pnu.rotateX(2*M_PI/360.*nu_angle_wrt_y_in_deg);
    p1.rotateX(2*M_PI/360.*nu_angle_wrt_y_in_deg);
    p2.rotateX(2*M_PI/360.*nu_angle_wrt_y_in_deg);  
  }

  return;
}
HepMC::GenParticle * GtInverseBeta::make_particle ( Hep3Vector &  vect,
double  mass,
int  pid 
) [private]

Definition at line 50 of file GtInverseBeta.cc.

{
  double energy = sqrt(vect.mag2()+mass*mass);
  HepMC::FourVector fourmom(vect.x()/MeV,vect.y()/MeV,vect.z()/MeV,
                            energy/MeV);
  //the new GenParticle will be deleted by GenVertex destructor
  return new HepMC::GenParticle(fourmom,pid,1/*=status*/);
}

Member Data Documentation

Definition at line 41 of file GtInverseBeta.h.

Definition at line 43 of file GtInverseBeta.h.

Definition at line 44 of file GtInverseBeta.h.

Definition at line 45 of file GtInverseBeta.h.


The documentation for this class was generated from the following files:
| Classes | Job Modules | Data Objects | Services | Algorithms | Tools | Packages | Directories | Tracs |

Generated on Fri May 16 2014 10:22:51 for InvBetaDecay by doxygen 1.7.4