G4UA::FluxRecorder Class Reference

#include <FluxRecorder.h>

Inheritance diagram for G4UA::FluxRecorder:

Collaboration diagram for G4UA::FluxRecorder:

Public Member Functions
	FluxRecorder ()
virtual void	BeginOfRunAction (const G4Run *) override
virtual void	EndOfRunAction (const G4Run *) override
virtual void	EndOfEventAction (const G4Event *) override
virtual void	UserSteppingAction (const G4Step *) override

Private Types
enum	scoringVolume {   RPCOlz, RPCOmz, RPCOhz, RPCMlz,   RPCMmz, RPCMhz, MDTIlz, MDTImz,   MDTIhz, LMDTo, LMDTm, LMDTi,   BMDTo, BMDTm, BMDTi, SWo,   SWm, SWi, SWc, SWt,   bflz, bfhz, bslz, bshz,   btlz, bthz, ffle, ffme,   ffhe, fsle, fsme, fshe,   ftle, ftme, fthe, CY1,   CY2, lastVol }

Private Member Functions
void	findVolume (const double, const double, const double, const double)

Private Attributes
TH1D *	m_flux [lastVol][9][2] = {}
TH1D *	m_fluxE [lastVol][9] = {}
double	m_nev
std::vector< int >	m_list

Detailed Description

Definition at line 22 of file FluxRecorder.h.

Member Enumeration Documentation

scoringVolume

enum G4UA::FluxRecorder::scoringVolume

private

Enumerator
RPCOlz
RPCOmz
RPCOhz
RPCMlz
RPCMmz
RPCMhz
MDTIlz
MDTImz
MDTIhz
LMDTo
LMDTm
LMDTi
BMDTo
BMDTm
BMDTi
SWo
SWm
SWi
SWc
SWt
bflz
bfhz
bslz
bshz
btlz
bthz
ffle
ffme
ffhe
fsle
fsme
fshe
ftle
ftme
fthe
CY1
CY2
lastVol

Definition at line 34 of file FluxRecorder.h.

                       {RPCOlz,RPCOmz,RPCOhz,RPCMlz,RPCMmz,RPCMhz,MDTIlz,MDTImz,MDTIhz,
                        LMDTo,LMDTm,LMDTi,BMDTo,BMDTm,BMDTi,SWo,SWm,SWi,SWc,SWt,
                        bflz,bfhz,bslz,bshz,btlz,bthz,ffle,ffme,ffhe,fsle,fsme,fshe,ftle,ftme,fthe,
                        CY1,CY2,lastVol};

Constructor & Destructor Documentation

FluxRecorder()

G4UA::FluxRecorder::FluxRecorder

(

)

Definition at line 26 of file FluxRecorder.cxx.

    : m_nev(0.0)
  {
  }

Member Function Documentation

BeginOfRunAction()

void G4UA::FluxRecorder::BeginOfRunAction ( const G4Run *  )

overridevirtual

Definition at line 31 of file FluxRecorder.cxx.

  {
    char nom[120];
    double timebins[101],ebins[101];
    for (int i=0;i<101;++i){
      timebins[i] = std::pow(10.,9.*double(i)/100.);
      ebins[i] = std::pow(10.,-11.+16.*double(i)/100.);
    }
    for (int i=0;i<lastVol;++i){
      for (int j=0;j<9;++j){
        for (int k=0;k<2;++k){
          sprintf(nom,"Flux_%i_%i_%i",i,j,k);
          m_flux[i][j][k] = new TH1D(nom,"",100,timebins);
        } // energy
        sprintf(nom,"FluxE_%i_%i",i,j);
        m_fluxE[i][j] = new TH1D(nom,"",100,ebins);
      } // PDGID
    } // Volume
 
  }

EndOfEventAction()

void G4UA::FluxRecorder::EndOfEventAction ( const G4Event *  )

overridevirtual

Definition at line 52 of file FluxRecorder.cxx.

  {
    m_nev+=1.;
  }

EndOfRunAction()

void G4UA::FluxRecorder::EndOfRunAction ( const G4Run *  )

overridevirtual

Definition at line 57 of file FluxRecorder.cxx.

                                               {
    TFile * f = new TFile("flux.root","RECREATE");
    f->cd();
    char nom[80];
    for (int i=0;i<lastVol;++i){
      for (int j=0;j<9;++j){
        for (int k=0;k<2;++k){
          sprintf(nom,"Flux_%i_%i_%i",i,j,k);
          m_flux[i][j][k]->Scale(1./m_nev);
          m_flux[i][j][k]->Write(nom);
        } // Energy
        sprintf(nom,"FluxE_%i_%i",i,j);
        m_fluxE[i][j]->Scale(1./m_nev);
        m_fluxE[i][j]->Write(nom);
      } // PDGID
    } // Volume
    f->Close();
  }

findVolume()

void G4UA::FluxRecorder::findVolume ( const double  r1, const double  z1, const double  r2, const double  z2  )

private

Definition at line 120 of file FluxRecorder.cxx.

  {
 
    // This horrible code should be fixed
    const static double dim[lastVol][4] = {
      {980.,1000.,0.,400.} , {980.,1000.,400.,800.} , {980.,1000.,800.,1200.} ,
      {750.,770.,0.,200.} , {750.,770.,200.,450.} , {750.,770.,450.,850.} ,
      {480.,540.,0.,200.} , {480.,540.,200.,400.} , {480.,540.,400.,600.} ,
      {500.,1000.,1320.,1400.} , {280.,500.,1320.,1400.} , {180.,280.,1320.,1400.} ,
      {600.,1200.,2120.,2180.} , {400.,600.,2120.,2180.} , {300.,400.,2120.,2180.} ,
      {450.,600.,700.,760.} , {320.,450.,700.,760.} , {220.,320.,700.,760.} , {100.,200.,720.,760.} , {220.,320.,680.,700.} ,
      {460.,460.1,0.,350.} , {460.,460.1,350.,655.} , {750.,750.1,0.,740.} , {750.,750.1,740.,1000.} ,
      {1020.,1020.1,0.,740.} , {1020.,1020.1,740.,1281.} ,
      {370.,532.,705.,705.1} , {194.,370.,705.,705.1} , {99.,190.,730.,690.} ,
      {683.,1020.,1301.,1301.1} , {264.,683.,1301.,1301.1} , {176.,264.,1301.,1301.1} ,
      {600.,1170.,2040.,2040.1} , {447.,700.,2207.,2207.1} , {298.,447.,2207.,2207.1} ,
      // Extras from Charlie
      {750.,770.,850.,950.} , {480.,540.,600.,720.} };
 
 
    for (int i=0;i<lastVol;++i){
      // Crossing outward in r over r1
      int hit = 0;
      if (i!=28){
        double myZ1 = z1+(r1!=r2?(z2-z1)/(r2-r1)*(dim[i][0]-r1):(z2-z1)*0.5);
        double myZ2 = z1+(r1!=r2?(z2-z1)/(r2-r1)*(dim[i][1]-r1):(z2-z1)*0.5);
        double myR1 = r1+(z1!=z2?(r2-r1)/(z2-z1)*(dim[i][2]-z1):(r2-r1)*0.5);
        double myR2 = r1+(z1!=z2?(r2-r1)/(z2-z1)*(dim[i][3]-z1):(r2-r1)*0.5);
 
        // Crossing outward in r over r1
        if      ( r1<dim[i][0] && r2>dim[i][0] && myZ1>dim[i][2] && myZ1<dim[i][3] ) hit = 1;
        // Crossing inward in r over r2
        else if ( r1>dim[i][1] && r2<dim[i][1] && myZ2>dim[i][2] && myZ2<dim[i][3] ) hit = 2;
        // Crossing rightward in z over z1
        else if ( z1<dim[i][2] && z2>dim[i][2] && myR1>dim[i][0] && myR1<dim[i][1] ) hit = 3;
        // Crossing leftward in z over z2
        else if ( z1>dim[i][3] && z2<dim[i][3] && myR2>dim[i][0] && myR2<dim[i][1] ) hit = 4;
 
        // Special handling for stupid FLUKA slanted volume
      } else {
        // corner test...
        if ( (r1==dim[i][0]&&z1==dim[i][2])||(r2==dim[i][0]&&r2==dim[i][2]) ) hit = 1;
        else {
 
          double denom = (r2 - r1)*(dim[i][3] - dim[i][2]) - (z2 - z1)*(dim[i][1] - dim[i][0]);
          double nume_a = (z2 - z1)*(dim[i][0] - r1) - (r2 - r1)*(dim[i][2] - z1);
          double nume_b = (dim[i][3] - dim[i][2])*(dim[i][0] - r1) - (dim[i][1] - dim[i][0])*(dim[i][2] - z1);
 
          if(!denom){
            if(!nume_a && !nume_b) hit=1; // they are the same line
          } else {
            double ua = nume_a / denom;
            double ub = nume_b / denom;
            if(ua >= 0. && ua <= 1. && ub >= 0. && ub <= 1.) hit=2; // They intersect
          }
        }
      }
 
      if (hit){
        /*
          if (myIndex!=-1)
          std::cout << "Particle moving from ( " << r1 << " , " << z1 << " ) to ( " << r2 << " , " << z2
          << " ) appears to cross boundaries " << myIndex << " and " << i << " with " << myZ1
          << " , " << myZ2 << " , " << myR1 << " , " << myR2 << " via " << hit << " of dim "
          << dim[i][0] << " , " << dim[i][1] << " , " << dim[i][2] << " , " << dim[i][3] << std::endl;
        */
        m_list.push_back(i);
      } // if we scored
    } // Loop over all volumes
  }

UserSteppingAction()

void G4UA::FluxRecorder::UserSteppingAction ( const G4Step *  aStep )

overridevirtual

Definition at line 76 of file FluxRecorder.cxx.

                                                          {
    int pdgid = 8, energy=(aStep->GetTrack()->GetKineticEnergy()>10.)?1:0;
    if (aStep->GetTrack()->GetDefinition()==G4Gamma::Definition()){
      pdgid=0;
      energy = (aStep->GetTrack()->GetKineticEnergy()>0.03)?1:0;
    } else if (aStep->GetTrack()->GetDefinition()==G4Proton::Definition()){
      pdgid=1;
      energy = (aStep->GetTrack()->GetKineticEnergy()>10.)?1:0;
    } else if (aStep->GetTrack()->GetDefinition()==G4PionPlus::Definition() ||
               aStep->GetTrack()->GetDefinition()==G4PionMinus::Definition()){
      pdgid=2;
      energy = (aStep->GetTrack()->GetKineticEnergy()>10.)?1:0;
    } else if(aStep->GetTrack()->GetDefinition()==G4MuonPlus::Definition() ||
              aStep->GetTrack()->GetDefinition()==G4MuonMinus::Definition()){
      pdgid=3;
      energy = (aStep->GetTrack()->GetKineticEnergy()>10.)?1:0;
    } else if(aStep->GetTrack()->GetDefinition()==G4Electron::Definition() ||
              aStep->GetTrack()->GetDefinition()==G4Positron::Definition()){
      if (aStep->GetTrack()->GetKineticEnergy()>0.5){
        pdgid=4; energy=1;
      } else {
        pdgid=5;
        energy=(aStep->GetTrack()->GetKineticEnergy()>0.01)?1:0;
      }
    } else if(aStep->GetTrack()->GetDefinition()==G4Neutron::Definition()){
      if (aStep->GetTrack()->GetKineticEnergy()>10.){
        pdgid=6; energy=1;
      } else {
        pdgid=7;
        energy=(aStep->GetTrack()->GetKineticEnergy()>0.1?1:0);
      }
    } else if (!aStep->GetTrack()->GetDefinition()->GetPDGCharge()) return; // Not one of those and not a primary?
 
    m_list.clear();
    findVolume( aStep->GetPreStepPoint()->GetPosition().rho()*0.1, std::fabs(aStep->GetPreStepPoint()->GetPosition().z()*0.1),
                aStep->GetPostStepPoint()->GetPosition().rho()*0.1, std::fabs(aStep->GetPostStepPoint()->GetPosition().z()*0.1) ); // units are cm
    if (m_list.size()==0) return;
 
    for (unsigned int i=0;i<m_list.size();++i){
      m_flux[m_list[i]][pdgid][energy]->Fill( aStep->GetTrack()->GetGlobalTime() );
      m_fluxE[m_list[i]][pdgid]->Fill( aStep->GetTrack()->GetKineticEnergy() );
    }
  }

Member Data Documentation

m_flux

TH1D* G4UA::FluxRecorder::m_flux[lastVol][9][2] = {}

private

Definition at line 38 of file FluxRecorder.h.

m_fluxE

TH1D* G4UA::FluxRecorder::m_fluxE[lastVol][9] = {}

private

Definition at line 39 of file FluxRecorder.h.

m_list

std::vector<int> G4UA::FluxRecorder::m_list

private

Definition at line 42 of file FluxRecorder.h.

m_nev

double G4UA::FluxRecorder::m_nev

private

Definition at line 41 of file FluxRecorder.h.

---

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

• FluxRecorder.h
• FluxRecorder.cxx