ATLAS Offline Software
Loading...
Searching...
No Matches
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

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.

Constructor & Destructor Documentation

◆ FluxRecorder()

G4UA::FluxRecorder::FluxRecorder ( )

Definition at line 26 of file FluxRecorder.cxx.

27 : m_nev(0.0)
28 {
29 }

Member Function Documentation

◆ BeginOfRunAction()

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

Definition at line 31 of file FluxRecorder.cxx.

32 {
33 char nom[120];
34 double timebins[101],ebins[101];
35 for (int i=0;i<101;++i){
36 timebins[i] = std::pow(10.,9.*double(i)/100.);
37 ebins[i] = std::pow(10.,-11.+16.*double(i)/100.);
38 }
39 for (int i=0;i<lastVol;++i){
40 for (int j=0;j<9;++j){
41 for (int k=0;k<2;++k){
42 sprintf(nom,"Flux_%i_%i_%i",i,j,k);
43 m_flux[i][j][k] = new TH1D(nom,"",100,timebins);
44 } // energy
45 sprintf(nom,"FluxE_%i_%i",i,j);
46 m_fluxE[i][j] = new TH1D(nom,"",100,ebins);
47 } // PDGID
48 } // Volume
49
50 }
TH1D * m_flux[lastVol][9][2]
TH1D * m_fluxE[lastVol][9]

◆ EndOfEventAction()

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

Definition at line 52 of file FluxRecorder.cxx.

53 {
54 m_nev+=1.;
55 }

◆ EndOfRunAction()

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

Definition at line 57 of file FluxRecorder.cxx.

57 {
58 TFile * f = new TFile("flux.root","RECREATE");
59 f->cd();
60 char nom[80];
61 for (int i=0;i<lastVol;++i){
62 for (int j=0;j<9;++j){
63 for (int k=0;k<2;++k){
64 sprintf(nom,"Flux_%i_%i_%i",i,j,k);
65 m_flux[i][j][k]->Scale(1./m_nev);
66 m_flux[i][j][k]->Write(nom);
67 } // Energy
68 sprintf(nom,"FluxE_%i_%i",i,j);
69 m_fluxE[i][j]->Scale(1./m_nev);
70 m_fluxE[i][j]->Write(nom);
71 } // PDGID
72 } // Volume
73 f->Close();
74 }

◆ 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.

121 {
122
123 // This horrible code should be fixed
124 const static double dim[lastVol][4] = {
125 {980.,1000.,0.,400.} , {980.,1000.,400.,800.} , {980.,1000.,800.,1200.} ,
126 {750.,770.,0.,200.} , {750.,770.,200.,450.} , {750.,770.,450.,850.} ,
127 {480.,540.,0.,200.} , {480.,540.,200.,400.} , {480.,540.,400.,600.} ,
128 {500.,1000.,1320.,1400.} , {280.,500.,1320.,1400.} , {180.,280.,1320.,1400.} ,
129 {600.,1200.,2120.,2180.} , {400.,600.,2120.,2180.} , {300.,400.,2120.,2180.} ,
130 {450.,600.,700.,760.} , {320.,450.,700.,760.} , {220.,320.,700.,760.} , {100.,200.,720.,760.} , {220.,320.,680.,700.} ,
131 {460.,460.1,0.,350.} , {460.,460.1,350.,655.} , {750.,750.1,0.,740.} , {750.,750.1,740.,1000.} ,
132 {1020.,1020.1,0.,740.} , {1020.,1020.1,740.,1281.} ,
133 {370.,532.,705.,705.1} , {194.,370.,705.,705.1} , {99.,190.,730.,690.} ,
134 {683.,1020.,1301.,1301.1} , {264.,683.,1301.,1301.1} , {176.,264.,1301.,1301.1} ,
135 {600.,1170.,2040.,2040.1} , {447.,700.,2207.,2207.1} , {298.,447.,2207.,2207.1} ,
136 // Extras from Charlie
137 {750.,770.,850.,950.} , {480.,540.,600.,720.} };
138
139
140 for (int i=0;i<lastVol;++i){
141 // Crossing outward in r over r1
142 int hit = 0;
143 if (i!=28){
144 double myZ1 = z1+(r1!=r2?(z2-z1)/(r2-r1)*(dim[i][0]-r1):(z2-z1)*0.5);
145 double myZ2 = z1+(r1!=r2?(z2-z1)/(r2-r1)*(dim[i][1]-r1):(z2-z1)*0.5);
146 double myR1 = r1+(z1!=z2?(r2-r1)/(z2-z1)*(dim[i][2]-z1):(r2-r1)*0.5);
147 double myR2 = r1+(z1!=z2?(r2-r1)/(z2-z1)*(dim[i][3]-z1):(r2-r1)*0.5);
148
149 // Crossing outward in r over r1
150 if ( r1<dim[i][0] && r2>dim[i][0] && myZ1>dim[i][2] && myZ1<dim[i][3] ) hit = 1;
151 // Crossing inward in r over r2
152 else if ( r1>dim[i][1] && r2<dim[i][1] && myZ2>dim[i][2] && myZ2<dim[i][3] ) hit = 2;
153 // Crossing rightward in z over z1
154 else if ( z1<dim[i][2] && z2>dim[i][2] && myR1>dim[i][0] && myR1<dim[i][1] ) hit = 3;
155 // Crossing leftward in z over z2
156 else if ( z1>dim[i][3] && z2<dim[i][3] && myR2>dim[i][0] && myR2<dim[i][1] ) hit = 4;
157
158 // Special handling for stupid FLUKA slanted volume
159 } else {
160 // corner test...
161 if ( (r1==dim[i][0]&&z1==dim[i][2])||(r2==dim[i][0]&&r2==dim[i][2]) ) hit = 1;
162 else {
163
164 double denom = (r2 - r1)*(dim[i][3] - dim[i][2]) - (z2 - z1)*(dim[i][1] - dim[i][0]);
165 double nume_a = (z2 - z1)*(dim[i][0] - r1) - (r2 - r1)*(dim[i][2] - z1);
166 double nume_b = (dim[i][3] - dim[i][2])*(dim[i][0] - r1) - (dim[i][1] - dim[i][0])*(dim[i][2] - z1);
167
168 if(!denom){
169 if(!nume_a && !nume_b) hit=1; // they are the same line
170 } else {
171 double ua = nume_a / denom;
172 double ub = nume_b / denom;
173 if(ua >= 0. && ua <= 1. && ub >= 0. && ub <= 1.) hit=2; // They intersect
174 }
175 }
176 }
177
178 if (hit){
179 /*
180 if (myIndex!=-1)
181 std::cout << "Particle moving from ( " << r1 << " , " << z1 << " ) to ( " << r2 << " , " << z2
182 << " ) appears to cross boundaries " << myIndex << " and " << i << " with " << myZ1
183 << " , " << myZ2 << " , " << myR1 << " , " << myR2 << " via " << hit << " of dim "
184 << dim[i][0] << " , " << dim[i][1] << " , " << dim[i][2] << " , " << dim[i][3] << std::endl;
185 */
186 m_list.push_back(i);
187 } // if we scored
188 } // Loop over all volumes
189 }
std::vector< int > m_list

◆ UserSteppingAction()

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

Definition at line 76 of file FluxRecorder.cxx.

76 {
77 int pdgid = 8, energy=(aStep->GetTrack()->GetKineticEnergy()>10.)?1:0;
78 if (aStep->GetTrack()->GetDefinition()==G4Gamma::Definition()){
79 pdgid=0;
80 energy = (aStep->GetTrack()->GetKineticEnergy()>0.03)?1:0;
81 } else if (aStep->GetTrack()->GetDefinition()==G4Proton::Definition()){
82 pdgid=1;
83 energy = (aStep->GetTrack()->GetKineticEnergy()>10.)?1:0;
84 } else if (aStep->GetTrack()->GetDefinition()==G4PionPlus::Definition() ||
85 aStep->GetTrack()->GetDefinition()==G4PionMinus::Definition()){
86 pdgid=2;
87 energy = (aStep->GetTrack()->GetKineticEnergy()>10.)?1:0;
88 } else if(aStep->GetTrack()->GetDefinition()==G4MuonPlus::Definition() ||
89 aStep->GetTrack()->GetDefinition()==G4MuonMinus::Definition()){
90 pdgid=3;
91 energy = (aStep->GetTrack()->GetKineticEnergy()>10.)?1:0;
92 } else if(aStep->GetTrack()->GetDefinition()==G4Electron::Definition() ||
93 aStep->GetTrack()->GetDefinition()==G4Positron::Definition()){
94 if (aStep->GetTrack()->GetKineticEnergy()>0.5){
95 pdgid=4; energy=1;
96 } else {
97 pdgid=5;
98 energy=(aStep->GetTrack()->GetKineticEnergy()>0.01)?1:0;
99 }
100 } else if(aStep->GetTrack()->GetDefinition()==G4Neutron::Definition()){
101 if (aStep->GetTrack()->GetKineticEnergy()>10.){
102 pdgid=6; energy=1;
103 } else {
104 pdgid=7;
105 energy=(aStep->GetTrack()->GetKineticEnergy()>0.1?1:0);
106 }
107 } else if (!aStep->GetTrack()->GetDefinition()->GetPDGCharge()) return; // Not one of those and not a primary?
108
109 m_list.clear();
110 findVolume( aStep->GetPreStepPoint()->GetPosition().rho()*0.1, std::fabs(aStep->GetPreStepPoint()->GetPosition().z()*0.1),
111 aStep->GetPostStepPoint()->GetPosition().rho()*0.1, std::fabs(aStep->GetPostStepPoint()->GetPosition().z()*0.1) ); // units are cm
112 if (m_list.size()==0) return;
113
114 for (unsigned int i=0;i<m_list.size();++i){
115 m_flux[m_list[i]][pdgid][energy]->Fill( aStep->GetTrack()->GetGlobalTime() );
116 m_fluxE[m_list[i]][pdgid]->Fill( aStep->GetTrack()->GetKineticEnergy() );
117 }
118 }
void findVolume(const double, const double, const double, const double)

Member Data Documentation

◆ m_flux

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

Definition at line 38 of file FluxRecorder.h.

38{};

◆ m_fluxE

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

Definition at line 39 of file FluxRecorder.h.

39{};

◆ 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: