EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
SmearMatrixDetector_0_2_B3T.cxx
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file SmearMatrixDetector_0_2_B3T.cxx
1 // Based on
2 // https://physdiv.jlab.org/DetectorMatrix/
3 // From November 21, 2020
4 //
5 // Tracking assumes a 3T field
6 //
7 // Important Notes:
8 // - Where ranges are given, the more conservative number is chosen.
9 // - Without available specifications, angular resolution is assumed to be perfect.
10 // - Momentum and energy acceptance specifications are NOT implemented
11 // -- a) "90% acceptance" is not specific enough to introduce efficiency
12 // -- b) It is not a priori clear whether cuts should apply to truth or smeared value
13 // -- A "looper" may not reach the detector.
14 // -- On the other hand, a calorimetry deposit may fluctuate above threshold
15 // - PID for pi/K/p is implemented with perfect PID in the acceptance (and momentum range)
16 // - Electron PID is not implemented. Very high purity is important and not handled well with the above approach
17 
20 #include "eicsmear/smear/Device.h"
22 #include "eicsmear/smear/Smearer.h"
25 #include <eicsmear/smear/Smear.h>
27 #include "Math/Vector4D.h"
28 
29 // declare static --> local to this file, won't clash with others
30 static double ThetaFromEta( const double eta );
31 
33  gSystem->Load("libeicsmear");
34 
35  // Create the detector object to hold all devices
36  Smear::Detector det;
37 
38  // The framework provides implementations of three kinematic calculation methods
39  // from smeared values
40  // NM - "Null method" uses the scattered lepton.
41  // DA - Double Angle method
42  // JB - Jacquet-Blondel method
43  // Users should be mindful of the limitations and assumptions in
44  // the calculation of smeared kinematics. Sophisticated calculations are best redone in user code
45  // from the smeared particles directly.
46  det.SetEventKinematicsCalculator("NM DA JB");
47 
48  // Perfect phi and theta for all particles
49  // ---------------------------------------
50  // total coverage of the handbook for tracker, ecal, and hcal is -3.5 < eta < 3.5
51  Smear::Acceptance::Zone AngleZoneCommon(ThetaFromEta ( 3.5 ),ThetaFromEta ( -3.5 ));
52  Smear::Device SmearThetaCommon(Smear::kTheta, "0.0");
53  SmearThetaCommon.Accept.AddZone(AngleZoneCommon);
54  SmearThetaCommon.Accept.SetGenre(Smear::kAll);
55  det.AddDevice(SmearThetaCommon);
56 
57  Smear::Device SmearPhiCommon(Smear::kPhi, "0.0");
58  SmearPhiCommon.Accept.AddZone(AngleZoneCommon);
59  SmearPhiCommon.Accept.SetGenre(Smear::kAll);
60  det.AddDevice(SmearPhiCommon);
61 
62  // Tracking (B = 3 T)
63  // --------
64  // Note: Smear::kCharged checks pdg charge, so includes muons (good)
65  // NOT implemented
66  // Minimum pT for B = 3 T:
67  // 150 MeV/c for -3.0 < eta < -2.5
68  // 220 MeV/c for -2.5 < eta < -2.0
69  // 160 MeV/c for -2.0 < eta < -1.5
70  // 300 MeV/c for -1.5 < eta < -1.0
71 
72  // eta = -3.5 -- -2.5
73  // sigma_p/p ~ 0.1% p + 2%
74  Smear::Acceptance::Zone TrackBack1Zone(ThetaFromEta ( -2.5 ),ThetaFromEta ( -3.5 ));
75  Smear::Device TrackBack1P(Smear::kP, "sqrt( pow ( 0.001*P*P, 2) + pow ( 0.02*P, 2) )");
76  TrackBack1P.Accept.AddZone(TrackBack1Zone);
77  TrackBack1P.Accept.SetCharge(Smear::kCharged);
78  det.AddDevice(TrackBack1P);
79 
80  // eta = -2.5 -- -1
81  // sigma_p/p ~ 0.02% p + 1%
82  Smear::Acceptance::Zone TrackBack2Zone(ThetaFromEta ( -1 ),ThetaFromEta ( -2.5 ));
83  Smear::Device TrackBack2P(Smear::kP, "sqrt( pow ( 0.0002*P*P, 2) + pow ( 0.01*P, 2) )");
84  TrackBack2P.Accept.AddZone(TrackBack2Zone);
85  TrackBack2P.Accept.SetCharge(Smear::kCharged);
86  det.AddDevice(TrackBack2P);
87 
88  // eta = -1 -- +1
89  // sigma_p/p ~ 0.02% p + 0.05%
90  Smear::Acceptance::Zone TrackBarrelZone(ThetaFromEta ( 1 ),ThetaFromEta ( -1 ));
91  Smear::Device TrackBarrelP(Smear::kP, "sqrt( pow ( 0.0002*P*P, 2) + pow ( 0.005*P, 2) )");
92  TrackBarrelP.Accept.AddZone(TrackBarrelZone);
93  TrackBarrelP.Accept.SetCharge(Smear::kCharged);
94  det.AddDevice(TrackBarrelP);
95 
96  // eta = 1 -- 2.5
97  // sigma_p/p ~ 0.02% p+1%
98  Smear::Acceptance::Zone TrackFwd2Zone(ThetaFromEta ( 2.5 ),ThetaFromEta ( 1 ));
99  Smear::Device TrackFwd2P(Smear::kP, "sqrt( pow ( 0.0002*P*P, 2) + pow ( 0.01*P, 2) )");
100  TrackFwd2P.Accept.AddZone(TrackFwd2Zone);
101  TrackFwd2P.Accept.SetCharge(Smear::kCharged);
102  det.AddDevice(TrackFwd2P);
103 
104  // eta = 2.5 -- 3.5
105  // sigma_p/p ~ 0.1% p+2%
106  Smear::Acceptance::Zone TrackFwd1Zone(ThetaFromEta ( 3.5 ),ThetaFromEta ( 2.5 ));
107  Smear::Device TrackFwd1P(Smear::kP, "sqrt( pow ( 0.001*P*P, 2) + pow ( 0.02*P, 2) )");
108  TrackFwd1P.Accept.AddZone(TrackFwd1Zone);
109  TrackFwd1P.Accept.SetCharge(Smear::kCharged);
110  det.AddDevice(TrackFwd1P);
111 
112  // PID
113  // ---
114  // Perfect pi/K/P PID in the region indicating >3 sigma separation
115  // Make sure to not cover more than is covered by the other detectors.
116  // No minimum momentum is used.
117  // Accept charged hadrons (excludes e, mu)
118 
119  // Back
120  // eta = -3.5 -- -1
121  // p < 10 GeV
122  Smear::Acceptance::Zone PidBackZone(ThetaFromEta(-1),ThetaFromEta(-3.5),
123  0., TMath::TwoPi(), // phi
124  0., TMath::Infinity(), // E
125  0., 10 ); // p
126  Smear::PerfectID PidBack;
127  PidBack.Accept.AddZone(PidBackZone);
130  det.AddDevice( PidBack );
131 
132  // Barrel
133  // eta = -1 -- 1
134  // p < 6 GeV
136  0., TMath::TwoPi(), // phi
137  0., TMath::Infinity(), // E
138  0., 6 ); // p
139  Smear::PerfectID PidBarrel;
140  PidBarrel.Accept.AddZone(PidBarrelZone);
141  PidBarrel.Accept.SetCharge(Smear::kCharged);
142  PidBarrel.Accept.SetGenre(Smear::kHadronic);
143  det.AddDevice( PidBarrel );
144 
145  // Forward
146  // eta = -1 -- 1
147  // p < 50 GeV
148  // not sure what to make of "(worse approaching 3.5)"
150  0., TMath::TwoPi(), // phi
151  0., TMath::Infinity(), // E
152  0., 50 ); // p
153  Smear::PerfectID PidFwd;
154  PidFwd.Accept.AddZone(PidFwdZone);
157  det.AddDevice( PidFwd );
158 
159  // EM Calorimeters
160  // ---------------
161  // Note: Smear::kElectromagnetic == gamma + e. Does not include muons (good)
162  // Calorimeter resolution usually given as sigma_E/E = A% / E + B%/Sqrt{E} + C%
163  // EIC Smear needs absolute sigma: sigma_E = Sqrt{A*A + B*B*E + C*C*E*E}
164 
165  // Back
166  // eta = -3.5 -- -2
167  // E> 20 MeV not used because it suppresses smeared-up particles
168  // 1%/E + 2.5%/sqrtE + 1%
169  // A B C
170  Smear::Acceptance::Zone EmcalBackZone(ThetaFromEta ( -2 ),ThetaFromEta ( -3.5 ));
171  Smear::Device EmcalBack(Smear::kE, "sqrt( pow ( 0.01,2 ) + pow( 0.025,2)*E + pow ( 0.01*E,2 ) )");
172  EmcalBack.Accept.AddZone(EmcalBackZone);
174  det.AddDevice(EmcalBack);
175 
176  // MidBack
177  // eta = -2 -- -1
178  // E> 50 MeV not used because it suppresses smeared-up particles
179  // 2%/E + (4-8)%/sqrtE + 2% (Upper limit achievable with 50 cm space. Better resolution requires ~65 cm space allocated )
180  // A B choose 8% C
181  Smear::Acceptance::Zone EmcalMidBackZone(ThetaFromEta ( -1 ),ThetaFromEta ( -2 ));
182  Smear::Device EmcalMidBack(Smear::kE, "sqrt( pow ( 0.02,2 ) + pow( 0.08,2)*E + pow ( 0.02*E,2 ) )");
183  EmcalMidBack.Accept.AddZone(EmcalMidBackZone);
185  det.AddDevice(EmcalMidBack);
186 
187  // Barrel
188  // eta = -1 -- 1
189  // E > 100 MeV (50 MeV if higher resolution)
190  // not used because it suppresses smeared-up particles
191  // 2%/EāŠ•(12-14)%/āˆšEāŠ•(2-3)% for 30 cm space
192  // A better stochastic term can be achieved with more space:
193  // 2.5% with crystals 35cm
194  // 10% sampling 40cm
195  // 4% SciGlass 65cm
196  // (*Better resolution requires ~65 cm space allocated )"
197  // choosing
198  // 2%/E + 14%/sqrtE + 3%
199  // A B C
200  Smear::Acceptance::Zone EmcalBarrelZone(ThetaFromEta ( 1 ),ThetaFromEta ( -1 ));
201  Smear::Device EmcalBarrel(Smear::kE, "sqrt( pow ( 0.02,2 ) + pow( 0.14,2)*E + pow ( 0.03*E,2 ) )");
202  EmcalBarrel.Accept.AddZone(EmcalBarrelZone);
204  det.AddDevice(EmcalBarrel);
205 
206  // Forward
207  // eta = 1 -- 3.5
208  // E > 50 MeV not used because it suppresses smeared-up particles
209  // "2%/E + (4*-12)%/sqrtE + 2% Upper limit achievable with 40cm space
210  // *Better resolution requires ~65 cm space allocated"
211  // choosing
212  // 2%/E + 12%/sqrtE + 2%
213  // A B C
214  Smear::Acceptance::Zone EmcalFwdZone(ThetaFromEta ( 3.5 ),ThetaFromEta ( 1 )); // E
215  Smear::Device EmcalFwd(Smear::kE, "sqrt( pow ( 0.02,2 ) + pow( 0.12,2)*E + pow ( 0.02*E,2 ) )");
216  EmcalFwd.Accept.AddZone(EmcalFwdZone);
218  det.AddDevice(EmcalFwd);
219 
220  // Hadronic Calorimeters
221  // ----------------------
222  // Note: kHadronic == |pdg|>110.
223 
224  // Back
225  // eta = -3.5 -- -1
226  // E>500 MeV not used because it suppresses smeared-up particles
227  // (Better resolution required more space and R&D)
228  // stoch. = 50%, const = 10%
229  Smear::Acceptance::Zone HcalBackZone(ThetaFromEta ( -1 ),ThetaFromEta ( -3.5 ) );
230  Smear::Device HcalBack(Smear::kE, "sqrt(pow( 0.1*E, 2) + pow ( 0.5,2) *E)");
231  HcalBack.Accept.AddZone(HcalBackZone);
232  HcalBack.Accept.SetGenre(Smear::kHadronic);
233  det.AddDevice(HcalBack);
234 
235  // Barrel
236  // eta = -1 -- 1
237  // E>500 MeV not used because it suppresses smeared-up particles
238  // stoch. = 100%, const = 10%
239  Smear::Acceptance::Zone HcalBarrelZone(ThetaFromEta ( 1 ),ThetaFromEta ( -1 ) ); // E
240  Smear::Device HcalBarrel(Smear::kE, "sqrt( pow( 0.1*E, 2) + pow( 1.0, 2)*E)");
241  HcalBarrel.Accept.AddZone(HcalBarrelZone);
242  HcalBarrel.Accept.SetGenre(Smear::kHadronic);
243  det.AddDevice(HcalBarrel);
244 
245  // Forward
246  // eta = 1 -- 3.5
247  // E>500 MeV not used because it suppresses smeared-up particles
248  // stoch. = 50%, const = 10%
249  // "(35%/sqrtE not achievable)"
250  Smear::Acceptance::Zone HcalFwdZone(ThetaFromEta ( 3.5 ),ThetaFromEta ( 1 ));
251  Smear::Device HcalFwd(Smear::kE, "sqrt(pow( 0.1*E, 2) + pow ( 0.5,2) *E)");
252  HcalFwd.Accept.AddZone(HcalFwdZone);
254  det.AddDevice(HcalFwd);
255 
256  // Not covered:
257  // Tracker material budget X/X0 <~5%
258  // Tracker pointing resolution
259  // Low-Q^2 tagger: -6.9<eta<-5.8: Delta_theta/theta < 1.5%; 10^-6 < Q2 < 10^-2 GeV2
260  // Proton spectrometer: eta>6.2: sigma_intrinsic(|t|)/|t| < 1%; Acceptance: 0.2 < pT < 1.2 GeV/c
261  // Barrel vertexing: sigma_xyz ~ 20 microns, d0(z) ~ d0(r phi) ~ (20 microns)/(pT [GeV]) + 5 microns
262  // Central muon detection
263 
264  return det;
265 }
266 
267 // -------------------------------------------------------------------
268 double ThetaFromEta( const double eta ) {
269  if ( !std::isnan(eta) && !std::isinf(eta) ) {
270  return 2.0 * atan( exp( -eta ));
271  }
272  throw std::runtime_error("ThetaFromEta called with NaN or Inf");
273  return -1;
274 }