EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
SmearBeAST_0_0.cxx
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file SmearBeAST_0_0.cxx
1 #include "eicsmear/erhic/VirtualParticle.h"
2 #include "eicsmear/smear/Acceptance.h"
3 #include "eicsmear/smear/Device.h"
4 #include "eicsmear/smear/Detector.h"
5 #include "eicsmear/smear/Smearer.h"
6 #include "eicsmear/smear/ParticleMCS.h"
7 #include "eicsmear/smear/PerfectID.h"
8 #include <eicsmear/smear/Smear.h>
9 #include <eicsmear/erhic/ParticleMC.h>
10 #include "Math/Vector4D.h"
11 
12 Smear::Detector BuildBeAST_0_0() {
13 
14  gSystem->Load("libeicsmear");
15 
16  // Calorimeter resolution usually given as sigma_E/E = const% + stocastic%/Sqrt{E}
17  // EIC Smear needs absolute sigma: sigma_E = Sqrt{const*const*E*E + stoc*stoc*E}
18 
19  // Create the EM Calorimeter
20  Smear::Device emcalBck(Smear::kE, "sqrt(0.01*0.01*E*E + 0.015*0.015*E)");
21  Smear::Device emcalMidBck(Smear::kE, "sqrt(0.01*0.01*E*E + 0.07*0.07*E)");
22  Smear::Device emcalMid(Smear::kE, "sqrt(0.01*0.01*E*E + 0.10*0.10*E)");
23  Smear::Device emcalFwd(Smear::kE, "sqrt(0.01*0.01*E*E + 0.07*0.07*E)");
24 
25 
26  // Create the Forward/Backward Hadron Calorimeter
27  Smear::Device hcalFwd(Smear::kE, "sqrt(0.015*0.015*E*E + 0.50*0.50*E)");
28  Smear::Device hcalBck(Smear::kE, "sqrt(0.015*0.015*E*E + 0.50*0.50*E)");
29 
30  // Create the Hypothetical Mid Rap Calorimeter
31  //Smear::Device hcalMid(Smear::kE, "sqrt(0.07*0.07*E*E + 0.85*0.85*E)"); // ~CMS
32  Smear::Device hcalMid(Smear::kE, "sqrt(0.02*0.02*E*E + 0.35*0.35*E)"); // ~Zeus
33 
34  // Create Forward/Backward Hadron Calorimeter to Measure Charged Hadrons after the Tracker
35  Smear::Device hcalTrkFwd(Smear::kE, "sqrt(0.015*0.015*E*E + 0.50*0.50*E)");
36  Smear::Device hcalTrkBck(Smear::kE, "sqrt(0.015*0.015*E*E + 0.50*0.50*E)");
37 
38  // Create our tracking capabilities, by a combination of mometum, theta and phi Devices.
39  // The momentum parametrization (a*p + b) gives sigma_P/P in percent.
40  // So Multiply through by P and divide by 100 to get absolute sigma_P
41  // Theta and Phi parametrizations give absolute sigma in miliradians
42 
43  // Track Momentum
44  //Smear::Device momentum(Smear::kP, "(P*P*(0.0182031 + 0.00921047*pow((-log(tan(theta/2.0))), 2) - 0.00291243*pow((-log(tan(theta/2.0))), 4) + 0.000264353*pow((-log(tan(theta/2.0))), 6)) + (0.209681 + 0.275144*pow((-log(tan(theta/2.0))), 2) - 0.0436536*pow((-log(tan(theta/2.0))), 4) + 0.00367412*pow((-log(tan(theta/2.0))), 6)))*0.01");
45  Smear::Device momentum(Smear::kP, "(P*P*(0.0182031 + 0.00921047*pow((-log(tan(theta/2.0))), 2) - 0.00291243*pow((-log(tan(theta/2.0))), 4) + 0.000264353*pow((-log(tan(theta/2.0))), 6)) + P*(0.209681 + 0.275144*pow((-log(tan(theta/2.0))), 2) - 0.0436536*pow((-log(tan(theta/2.0))), 4) + 0.00367412*pow((-log(tan(theta/2.0))), 6)))*0.01");
46  Smear::Device trackTheta(Smear::kTheta, "((1.0/(1.0*P))*(0.752935 + 0.280370*pow((-log(tan(theta/2.0))), 2) - 0.0359713*pow((-log(tan(theta/2.0))), 4) + 0.00200623*pow((-log(tan(theta/2.0))), 6)) + 0.0282315 - 0.00998623*pow((-log(tan(theta/2.0))), 2) + 0.00117487*pow((-log(tan(theta/2.0))), 4) - 0.0000443918*pow((-log(tan(theta/2.0))), 6))*0.001");
47  Smear::Device trackPhi(Smear::kPhi, "((1.0/(1.0*P))*(0.743977 + 0.753393*pow((-log(tan(theta/2.0))), 2) + 0.0634184*pow((-log(tan(theta/2.0))), 4) + 0.0128001*pow((-log(tan(theta/2.0))), 6)) + 0.0308753 + 0.0480770*pow((-log(tan(theta/2.0))), 2) - 0.0129859*pow((-log(tan(theta/2.0))), 4) + 0.00109374*pow((-log(tan(theta/2.0))), 6))*0.001");
48 
49 
50  // Need these to keep the component not smeared
51  // Momentum for EM
52  Smear::Device momentumEM(Smear::kP, "0");
53  Smear::Device trackThetaEM(Smear::kTheta, "0");
54  Smear::Device trackPhiEM(Smear::kPhi, "0");
55 
56  // Momentum for Neutral Hadrons
57  Smear::Device momentumHad(Smear::kP, "0");
58  Smear::Device trackThetaHad(Smear::kTheta, "0");
59  Smear::Device trackPhiHad(Smear::kPhi, "0");
60 
61  // Momentum for Charged Hadrons After the Tracker
62  Smear::Device momentumHadTrkBck(Smear::kP, "0");
63  Smear::Device trackThetaHadTrkBck(Smear::kTheta, "0");
64  Smear::Device trackPhiHadTrkBck(Smear::kPhi, "0");
65 
66  Smear::Device momentumHadTrkFwd(Smear::kP, "0");
67  Smear::Device trackThetaHadTrkFwd(Smear::kTheta, "0");
68  Smear::Device trackPhiHadTrkFwd(Smear::kPhi, "0");
69 
70  // Energy for Tracks
71  Smear::Device trackEnergy(Smear::kE, "0");
72 
73  // Create a smearer for momentum
74  //Smear::Device momentum(Smear::kP, "0.01 * P");
75  //Smear::Device theta(Smear::kTheta, "0.05 * P");
76  //Smear::Device phi(Smear::kPhi, "0"); // "0" indicates perf585ect performance i.e. sigma(phi) = 0
77 
78  // Create a based on Hermes RICH.
79  //Smear::ParticleID rich("PIDMatrix.dat");
80 
81 
82  // Set Spatial Extent of Devices (in radians)
83  // eta = -4.5 -> theta = 3.1194
84  // eta = -4.0 -> theta = 3.1050
85  // eta = -3.5 -> theta = 3.0812
86  // eta = -3.0 -> theta = 3.0421
87  // eta = -2.0 -> theta = 2.8726
88  // eta = -1.0 -> theta = 2.4366
89  // eta = +1.0 -> theta = 0.7050
90  // eta = +2.0 -> theta = 0.2690
91  // eta = +3.0 -> theta = 0.0995
92  // eta = +3.5 -> theta = 0.0604
93  // eta = +4.0 -> theta = 0.0366
94  // eta = +4.5 -> theta = 0.0222
95 
96  // The BeAST design has calorimetry extending to +-4, here extend to +-4.5 to match the extent of my particle level jets
97  // Can place a cut on jet thrust axis in the analysis if need be
98 
99  // Set Up EMCal Zones
100  Smear::Acceptance::Zone emBck(2.8726,3.1194,0.,TMath::TwoPi(),0.,TMath::Infinity(),0.,TMath::Infinity(),0.,TMath::Infinity(),-TMath::Infinity(),TMath::Infinity());
101  Smear::Acceptance::Zone emMidBck(2.4366,2.8726,0.,TMath::TwoPi(),0.,TMath::Infinity(),0.,TMath::Infinity(),0.,TMath::Infinity(),-TMath::Infinity(),TMath::Infinity());
102  Smear::Acceptance::Zone emMid(0.7050,2.4366,0.,TMath::TwoPi(),0.,TMath::Infinity(),0.,TMath::Infinity(),0.,TMath::Infinity(),-TMath::Infinity(),TMath::Infinity());
103  Smear::Acceptance::Zone emFwd(0.0222,0.7050,0.,TMath::TwoPi(),0.,TMath::Infinity(),0.,TMath::Infinity(),0.,TMath::Infinity(),-TMath::Infinity(),TMath::Infinity());
104  Smear::Acceptance::Zone emTot(0.0222,3.1194,0.,TMath::TwoPi(),0.,TMath::Infinity(),0.,TMath::Infinity(),0.,TMath::Infinity(),-TMath::Infinity(),TMath::Infinity());
105 
106  // Set Up HCal Zone
107  Smear::Acceptance::Zone hBck(2.4366,3.1194,0.,TMath::TwoPi(),0.,TMath::Infinity(),0.,TMath::Infinity(),0.,TMath::Infinity(),-TMath::Infinity(),TMath::Infinity());
108  Smear::Acceptance::Zone hFwd(0.0222,0.7050,0.,TMath::TwoPi(),0.,TMath::Infinity(),0.,TMath::Infinity(),0.,TMath::Infinity(),-TMath::Infinity(),TMath::Infinity());
109  Smear::Acceptance::Zone hMid(0.7050,2.4366,0.,TMath::TwoPi(),0.,TMath::Infinity(),0.,TMath::Infinity(),0.,TMath::Infinity(),-TMath::Infinity(),TMath::Infinity());
110  Smear::Acceptance::Zone hTot(0.0222,3.1194,0.,TMath::TwoPi(),0.,TMath::Infinity(),0.,TMath::Infinity(),0.,TMath::Infinity(),-TMath::Infinity(),TMath::Infinity());
111  // Do HCal Over Full Acceptance (-4,4) and just kill neutrals at midrapidity in the jet finder for standard BeAST
112 
113  // Set Up HCal After Tracker Zone
114  Smear::Acceptance::Zone hTrkBck(3.0812,3.1194,0.,TMath::TwoPi(),0.,TMath::Infinity(),0.,TMath::Infinity(),0.,TMath::Infinity(),-TMath::Infinity(),TMath::Infinity());
115  Smear::Acceptance::Zone hTrkFwd(0.0222,0.0604,0.,TMath::TwoPi(),0.,TMath::Infinity(),0.,TMath::Infinity(),0.,TMath::Infinity(),-TMath::Infinity(),TMath::Infinity());
116 
117  // Set Up Tracking Zone
118  Smear::Acceptance::Zone trk(0.0604,3.0812,0.,TMath::TwoPi(),0.,TMath::Infinity(),0.,TMath::Infinity(),0.,TMath::Infinity(),-TMath::Infinity(),TMath::Infinity());
119 
120 
121  // Assign acceptance to calorimeters
122  emcalBck.Accept.SetGenre(Smear::kElectromagnetic);
123  emcalMidBck.Accept.SetGenre(Smear::kElectromagnetic);
124  emcalMid.Accept.SetGenre(Smear::kElectromagnetic);
125  emcalFwd.Accept.SetGenre(Smear::kElectromagnetic);
126 
127  emcalBck.Accept.AddZone(emBck);
128  emcalMidBck.Accept.AddZone(emMidBck);
129  emcalMid.Accept.AddZone(emMid);
130  emcalFwd.Accept.AddZone(emFwd);
131 
132  hcalBck.Accept.AddParticle(2112);
133  hcalBck.Accept.AddParticle(-2112);
134  hcalBck.Accept.AddParticle(130);
135  hcalFwd.Accept.AddParticle(2112);
136  hcalFwd.Accept.AddParticle(-2112);
137  hcalFwd.Accept.AddParticle(130);
138  hcalMid.Accept.AddParticle(2112);
139  hcalMid.Accept.AddParticle(-2112);
140  hcalMid.Accept.AddParticle(130);
141 
142  hcalBck.Accept.AddZone(hBck);
143  hcalFwd.Accept.AddZone(hFwd);
144  hcalMid.Accept.AddZone(hMid);
145 
146  // Assign acceptance to calorimeters for charged hadrons past tracker
147  hcalTrkBck.Accept.SetGenre(Smear::kHadronic);
148  hcalTrkFwd.Accept.SetGenre(Smear::kHadronic);
149 
150  hcalTrkBck.Accept.SetCharge(Smear::kCharged);
151  hcalTrkFwd.Accept.SetCharge(Smear::kCharged);
152 
153  hcalTrkBck.Accept.AddZone(hTrkBck);
154  hcalTrkFwd.Accept.AddZone(hTrkFwd);
155 
156  // Assign acceptance to tracker
157  momentum.Accept.SetGenre(Smear::kHadronic);
158  trackTheta.Accept.SetGenre(Smear::kHadronic);
159  trackPhi.Accept.SetGenre(Smear::kHadronic);
160 
161  momentum.Accept.SetCharge(Smear::kCharged);
162  trackTheta.Accept.SetCharge(Smear::kCharged);
163  trackPhi.Accept.SetCharge(Smear::kCharged);
164 
165  momentum.Accept.AddZone(trk);
166  trackTheta.Accept.AddZone(trk);
167  trackPhi.Accept.AddZone(trk);
168 
169  // Assign acceptance for calorimeter momentum
170  momentumEM.Accept.SetGenre(Smear::kElectromagnetic);
171  trackThetaEM.Accept.SetGenre(Smear::kElectromagnetic);
172  trackPhiEM.Accept.SetGenre(Smear::kElectromagnetic);
173 
174  momentumEM.Accept.AddZone(emTot);
175  trackThetaEM.Accept.AddZone(emTot);
176  trackPhiEM.Accept.AddZone(emTot);
177 
178  momentumHad.Accept.AddParticle(2112);
179  trackThetaHad.Accept.AddParticle(2112);
180  trackPhiHad.Accept.AddParticle(2112);
181 
182  momentumHad.Accept.AddParticle(-2112);
183  trackThetaHad.Accept.AddParticle(-2112);
184  trackPhiHad.Accept.AddParticle(-2112);
185 
186  momentumHad.Accept.AddParticle(130);
187  trackThetaHad.Accept.AddParticle(130);
188  trackPhiHad.Accept.AddParticle(130);
189 
190  momentumHad.Accept.AddZone(hTot);
191  trackThetaHad.Accept.AddZone(hTot);
192  trackPhiHad.Accept.AddZone(hTot);
193 
194  // Assign acceptance for calorimeter momentum for charged hadrons past tracker
195  momentumHadTrkBck.Accept.SetGenre(Smear::kHadronic);
196  trackThetaHadTrkBck.Accept.SetGenre(Smear::kHadronic);
197  trackPhiHadTrkBck.Accept.SetGenre(Smear::kHadronic);
198 
199  momentumHadTrkBck.Accept.SetCharge(Smear::kCharged);
200  trackThetaHadTrkBck.Accept.SetCharge(Smear::kCharged);
201  trackPhiHadTrkBck.Accept.SetCharge(Smear::kCharged);
202 
203  momentumHadTrkFwd.Accept.SetGenre(Smear::kHadronic);
204  trackThetaHadTrkFwd.Accept.SetGenre(Smear::kHadronic);
205  trackPhiHadTrkFwd.Accept.SetGenre(Smear::kHadronic);
206 
207  momentumHadTrkFwd.Accept.SetCharge(Smear::kCharged);
208  trackThetaHadTrkFwd.Accept.SetCharge(Smear::kCharged);
209  trackPhiHadTrkFwd.Accept.SetCharge(Smear::kCharged);
210 
211  momentumHadTrkBck.Accept.AddZone(hTrkBck);
212  trackThetaHadTrkBck.Accept.AddZone(hTrkBck);
213  trackPhiHadTrkBck.Accept.AddZone(hTrkBck);
214 
215  momentumHadTrkFwd.Accept.AddZone(hTrkFwd);
216  trackThetaHadTrkFwd.Accept.AddZone(hTrkFwd);
217  trackPhiHadTrkFwd.Accept.AddZone(hTrkFwd);
218 
219  // Assign acceptance for track energy
220  trackEnergy.Accept.SetGenre(Smear::kHadronic);
221 
222  trackEnergy.Accept.SetCharge(Smear::kCharged);
223 
224  trackEnergy.Accept.AddZone(trk);
225 
226 
227  //emcal.Accept.AddZone(central);
228  //momentum.Accept.AddZone(central);
229  //theta.Accept.AddZone(central);
230  //phi.Accept.AddZone(central);
232 
233  // Create a detector and add the devices
234  Smear::Detector det;
235  det.AddDevice(emcalBck);
236  det.AddDevice(emcalMidBck);
237  det.AddDevice(emcalMid);
238  det.AddDevice(emcalFwd);
239  det.AddDevice(hcalBck);
240  det.AddDevice(hcalFwd);
241  det.AddDevice(hcalMid);
242  det.AddDevice(hcalTrkBck);
243  det.AddDevice(hcalTrkFwd);
244  det.AddDevice(momentum);
245  det.AddDevice(trackTheta);
246  det.AddDevice(trackPhi);
247  // det.AddDevice(momentumEM);
248  // det.AddDevice(trackThetaEM);
249  // det.AddDevice(trackPhiEM);
250  det.AddDevice(momentumHad);
251  det.AddDevice(trackThetaHad);
252  det.AddDevice(trackPhiHad);
253  det.AddDevice(momentumHadTrkBck);
254  det.AddDevice(trackThetaHadTrkBck);
255  det.AddDevice(trackPhiHadTrkBck);
256  det.AddDevice(momentumHadTrkFwd);
257  det.AddDevice(trackThetaHadTrkFwd);
258  det.AddDevice(trackPhiHadTrkFwd);
259  det.AddDevice(trackEnergy);
260  det.SetEventKinematicsCalculator("NM JB DA"); // The detector will calculate event kinematics from smeared values
261 
262  //det.AddDevice(emcal);
263  //det.AddDevice(momentum);
264  //det.AddDevice(theta);
265  //det.AddDevice(phi);
267 
268 
269 
270 
271  return det;
272 }