ATLAS Offline Software
InducedChargeModel.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2022 CERN for the benefit of the ATLAS collaboration
3 */
4 
5 //-----------------------------------------------
6 // 2020.8.12 Implementation in Athena by Manabu Togawa
7 // 2017.7.24 update for the case of negative VD (type-P)
8 // 2017.7.17 updated
9 // 2016.4.3 for ICM (induced charge model) by Taka Kondo (KEK)
10 //-----------------------------------------------
11 
13 
14 // ROOT
15 #include "TH2.h"
16 #include "TFile.h"
17 
18 // C++ Standard Library
19 #include <algorithm>
20 #include <cmath>
21 #include <iostream>
22 
24 const double InducedChargeModel::s_e = Gaudi::Units::e_SI; // [Coulomb]
25 // Sorted depletion voltage values
26 const std::vector<float> InducedChargeModel::s_VFD0 =
27  { -180., -150., -120., -90., -60., -30., 0., 30., 70.};
28 const std::vector<std::string> InducedChargeModel::s_VFD0str =
29  {"M180", "M150", "M120", "M90", "M60", "M30", "0", "30", "70"};
30 
32  AthMessaging("InducedChargeModel"),
33  m_EFieldModel(model)
34 {
36  m_data.resize(maxHash);
37 }
38 
41  const float vbias,
43  const AtlasFieldCacheCondObj* fieldCondObj,
44  const ToolHandle<ISiliconConditionsTool>& siConditionsTool,
45  CLHEP::HepRandomEngine* rndmEngine,
46  const EventContext& ctx) const {
47  std::lock_guard<std::mutex> lock(m_mutex);
48 
49  // If cache exists, cache is used.
50  SCT_InducedChargeModelData* p_data = m_data[element->identifyHash()].get();
51  if (p_data) return p_data;
52 
53  ATH_MSG_DEBUG("--- Induced Charged Model Paramter (Begin) --------");
54 
55  HepGeom::Point3D<double> localpos(0., 0., 0.); // The center of wafer
56  HepGeom::Point3D<double> globalpos = element->globalPositionHit(localpos);
57  double point[3] = {globalpos.x(), globalpos.y(), globalpos.z()};
58  double field[3] = {0., 0., 0.};
59  MagField::AtlasFieldCache fieldCache;
60  fieldCondObj->getInitializedCache(fieldCache);
61  fieldCache.getField(point, field);
62  const Amg::Vector3D magneticField(field[0], field[1], field[2]); // in kTesla
63 
64  //---------------setting basic parameters---------------------------
65  std::unique_ptr<SCT_InducedChargeModelData> data =
66  std::make_unique<SCT_InducedChargeModelData>(vdepl,
67  vbias,
68  element,
69  magneticField,
72  siConditionsTool,
73  rndmEngine,
74  ctx);
75 
76  //--- set electric fields ---
77  setEField(*data);
78 
79  //--------------------------------------------------------
80 
81  ATH_MSG_DEBUG(" EFieldModel 0 (Flat Diode Model), 1 (FEM solusions) 2 (uniform E)\t"<< data->m_EFieldModel);
82  ATH_MSG_DEBUG(" DepletionVoltage_VD\t"<< data->m_VD <<" V");
83  ATH_MSG_DEBUG(" BiasVoltage_VB \t"<< data->m_VB <<" V");
84  ATH_MSG_DEBUG(" MagneticField_B \t"
85  << data->m_magneticField[Amg::x] << " "
86  << data->m_magneticField[Amg::y] << " "
87  << data->m_magneticField[Amg::z] <<" kTesla");
88  ATH_MSG_DEBUG(" Temperature_T \t"<< data->m_T <<" K");
89  ATH_MSG_DEBUG(" TransportTimeStep \t"<< m_transportTimeStep <<" ns");
90  ATH_MSG_DEBUG(" TransportTimeMax\t"<< m_transportTimeMax <<" ns");
91  ATH_MSG_DEBUG(" BulkDepth\t\t"<< m_bulk_depth << " cm");
92  ATH_MSG_DEBUG(" DepletionDepth\t\t"<< data->m_depletion_depth << " cm");
93  ATH_MSG_DEBUG(" StripPitch\t\t"<< m_strip_pitch << " cm");
94  ATH_MSG_DEBUG("--- Induced Charged Model Paramter (End) ---------");
95 
96  m_data[element->identifyHash()] = std::move(data);
97 
98  return m_data[element->identifyHash()].get();
99 }
100 
101 //---------------------------------------------------------------------
102 // transport
103 //---------------------------------------------------------------------
105  const bool isElectron,
106  const double x0, const double y0,
107  double* Q_m2, double* Q_m1, double* Q_00, double* Q_p1, double* Q_p2,
108  const IdentifierHash hashId,
109  const ToolHandle<ISiPropertiesTool>& siPropertiesTool,
110  const EventContext& ctx) const {
111  // transport electrons and holes in the bulk
112  // T. Kondo, 2010.9.9, 2017.7.20
113  // External parameters to be specified
114  // m_transportTimeMax [nsec]
115  // m_transportTimeStep [nsec]
116  // m_bulk_depth [cm]
117  // Induced currents are added to every m_transportTimeStep (defailt is 0.5 ns):
118  // Q_m2[100],Q_m1[100],Q_00[100],Q_p1[100],Q_p2[100] NTransportSteps=100
119  // means 50 ns storages of charges in each strips.
120 
121  double x = x0; // original electron/hole position [cm]
122  double y = y0; // original electron/hole position [cm]
123  bool isInBulk = true;
124  double t_current = 0.;
125  double qstrip[NStrips]; // induced charges on strips due to an electron or a hole
126  double vx, vy, D;
127  for (int istrip=StartStrip; istrip<=EndStrip; istrip++) {
128  qstrip[istrip+Offset] = (isElectron ? -1. : +1.) * induced(data, istrip, x, y); // original induced charge by electron or hole
129  }
130  while (t_current < m_transportTimeMax) {
131  if (!isInBulk) break;
132  if (!getVxVyD(data, isElectron, x, y, vx, vy, D, hashId, siPropertiesTool, ctx)) break;
133  double delta_y = vy * m_transportTimeStep / Gaudi::Units::second; // ns -> s
134  y += delta_y;
135  double dt = m_transportTimeStep; // [nsec]
136  if (isElectron) {
137  if (y < m_y_origin_min) {
138  isInBulk = false;
139  dt = (m_y_origin_min - (y-delta_y))/delta_y * m_transportTimeStep;
140  y = m_y_origin_min;
141  }
142  } else {
143  if (y > m_bulk_depth) {
144  isInBulk = false; // outside bulk region
145  dt = (m_bulk_depth - (y-delta_y))/delta_y * m_transportTimeStep;
146  y = m_bulk_depth;
147  }
148  }
149  t_current += dt;
150  x += vx * dt / Gaudi::Units::second; // ns -> s
151  double diffusion = std::sqrt(2.* D * dt / Gaudi::Units::second); // ns -> s
152  y += diffusion * CLHEP::RandGaussZiggurat::shoot(data.m_rndmEngine, 0.0, 1.0);
153  x += diffusion * CLHEP::RandGaussZiggurat::shoot(data.m_rndmEngine, 0.0, 1.0);
154  if (isElectron) {
155  if (y < m_y_origin_min) {
156  y = m_y_origin_min;
157  isInBulk = false;
158  }
159  } else {
160  if (y > m_bulk_depth) {
161  y = m_bulk_depth;
162  isInBulk = false; // outside bulk region
163  }
164  }
165 
166  // get induced current by subtracting induced charges
167  for (int istrip=StartStrip; istrip<=EndStrip; istrip++) {
168  double qnew = (isElectron ? -1. : +1.) * induced(data, istrip, x, y);
169  int jj = istrip + Offset;
170  double dq = qnew - qstrip[jj];
171  qstrip[jj] = qnew;
172 
173  int jt = static_cast<int>((t_current+0.001) / m_transportTimeStep);
174  if (jt < NTransportSteps) {
175  switch (istrip) {
176  case -2: Q_m2[jt] += dq; break;
177  case -1: Q_m1[jt] += dq; break;
178  case 0: Q_00[jt] += dq; break;
179  case +1: Q_p1[jt] += dq; break;
180  case +2: Q_p2[jt] += dq; break;
181  }
182  }
183 
184  }
185  } // end of electron or hole tracing
186 }
187 
188 //---------------------------------------------------------------------
189 // holeTransport
190 //---------------------------------------------------------------------
192  const double x0, const double y0,
193  double* Q_m2, double* Q_m1, double* Q_00, double* Q_p1, double* Q_p2,
194  const IdentifierHash hashId,
195  const ToolHandle<ISiPropertiesTool>& siPropertiesTool,
196  const EventContext& ctx) const {
197  // Q_m2[100],Q_m1[100],Q_00[100],Q_p1[100],Q_p2[100] NTransportSteps=100
198  const bool isElectron = false;
199  transport(data,
200  isElectron,
201  x0, y0,
202  Q_m2, Q_m1, Q_00, Q_p1, Q_p2,
203  hashId,
204  siPropertiesTool,
205  ctx);
206 }
207 
208 //---------------------------------------------------------------------
209 // electronTransport
210 //---------------------------------------------------------------------
212  const double x0, const double y0,
213  double* Q_m2, double* Q_m1, double* Q_00, double* Q_p1, double* Q_p2,
214  const IdentifierHash hashId,
215  const ToolHandle<ISiPropertiesTool>& siPropertiesTool,
216  const EventContext& ctx) const {
217  // Q_m2[100],Q_m1[100],Q_00[100],Q_p1[100],Q_p2[100] NTransportSteps=100
218  const bool isElectron = true;
219  transport(data,
220  isElectron,
221  x0, y0,
222  Q_m2, Q_m1, Q_00, Q_p1, Q_p2,
223  hashId,
224  siPropertiesTool,
225  ctx);
226 }
227 
228 //---------------------------------------------------------------
229 // Get parameters for electron or hole transport
230 //---------------------------------------------------------------
232  const bool isElectron,
233  const double x, const double y, double& vx, double& vy, double& D,
234  const IdentifierHash hashId,
235  const ToolHandle<ISiPropertiesTool>& siPropertiesTool,
236  const EventContext& ctx) const {
237  double Ex, Ey;
238  getEField(data, x, y, Ex, Ey); // [V/cm]
239  if (Ey > 0.) {
240  if (isElectron) {
241  // because electron has negative charge
242  Ex *= -1.;
243  Ey *= -1.;
244  }
245  const double E = std::sqrt(Ex*Ex+Ey*Ey);
246  const double mu = (isElectron ?
247  siPropertiesTool->getSiProperties(hashId, ctx).calcElectronDriftMobility(data.m_T, E*CLHEP::volt/CLHEP::cm) :
248  siPropertiesTool->getSiProperties(hashId, ctx).calcHoleDriftMobility (data.m_T, E*CLHEP::volt/CLHEP::cm))
250  const double v = mu * E;
251  const double r = (isElectron ?
252  siPropertiesTool->getSiProperties(hashId, ctx).calcElectronHallFactor(data.m_T) :
253  siPropertiesTool->getSiProperties(hashId, ctx).calcHoleHallFactor(data.m_T));
254 
255  const double tanLA = data.m_element->design().readoutSide()
256  * r * mu * (isElectron ? -1. : +1.) // because e has negative charge.
257  * data.m_element->hitDepthDirection()
258  * data.m_element->hitPhiDirection()
259  * (data.m_element->normal().cross(data.m_magneticField)).dot(data.m_element->phiAxis())
260  / Gaudi::Units::tesla * Gaudi::Units::gauss * 1000.; // 1000. is to convert kTesla to Tesla.
261 
262  const double secLA = std::sqrt(1.+tanLA*tanLA);
263  const double cosLA = 1./secLA;
264  const double sinLA = tanLA/secLA;
265  vy = v * (Ey*cosLA - Ex*sinLA)/E;
266  vx = v * (Ex*cosLA + Ey*sinLA)/E;
267  D = s_kB * data.m_T * mu/ s_e;
268  return true;
269  }
270  return false;
271 }
272 
273 //-------------------------------------------------------------------
274 // calculation of induced charge using Weighting (Ramo) function
275 //-------------------------------------------------------------------
277  const int istrip, const double x, const double y) const
278 {
279  // x and y are the location of charge (e or hole)
280  // induced charge on the strip "istrip" situated at the height y = d
281  // the center of the strip (istrip=0) is x = 0.004 [cm]
282 
283  // Ramo potential : 80 divisions (81 points) with 5 um intervals from 40-440 um.
284  // In sensor depth directions 114 divisions (115 points) with 2.5 nm intervals for 285 um.
285  const double deltax = 0.0005, deltay = 0.00025;
286 
287  if (y < 0. || y > m_bulk_depth) return 0.;
288  double xc = m_strip_pitch * (istrip + 0.5);
289  double dx = std::abs(x-xc);
290  int ix = static_cast<int>(dx / deltax);
291  if (ix >= NRamoPoints) return 0.;
292  int iy = static_cast<int>(y / deltay);
293  double fx = (dx - ix*deltax) / deltax;
294  double fy = ( y - iy*deltay) / deltay;
295  int ix1 = ix + 1;
296  int iy1 = iy + 1;
297  double P = m_PotentialValue[data.m_EFieldModel][ix ][iy ] *(1.-fx)*(1.-fy)
298  + m_PotentialValue[data.m_EFieldModel][ix1][iy ] * fx *(1.-fy)
299  + m_PotentialValue[data.m_EFieldModel][ix ][iy1] *(1.-fx)* fy
300  + m_PotentialValue[data.m_EFieldModel][ix1][iy1] * fx* fy;
301 
302  return P;
303 }
304 
305 //--------------------------------------------------------------
306 // Electric Field Ex and Ey
307 //--------------------------------------------------------------
309  const double x, const double y, double& Ex, double& Ey) const {
310  // m_EFieldModel == FlatDiodeModel 0; flat diode model
311  // m_EFieldModel == FEMsolutions 1; FEM solusions
312  // m_EFieldModel == UniformE 2; uniform E-field model
313  // x == 0.0040 [cm] : at the center of strip
314  // x must be within 0 and 0.008 [cm]
315 
316  // Electric field : 16 divisions (17 points) with 2.5 um intervals from 0-40 um.
317  // In sensor depth directions 114 divisions (115 points) with 2.5 nm intervals for 285 um.
318  const double deltay=0.00025, deltax=0.00025; // [cm]
319 
320  Ex = 0.;
321  Ey = 0.;
322  if (y < 0. || y > m_bulk_depth) return;
323  //---------- case for FEM analysis solution -------------------------
324  if (data.m_EFieldModel==FEMsolutions) {
325  int iy = static_cast<int>(y/deltay);
326  double fy = (y-iy*deltay) / deltay;
327  double xhalfpitch = m_strip_pitch/2.;
328  double x999 = 999.*m_strip_pitch;
329  double xx = fmod(x+x999, m_strip_pitch);
330  double xxx = xx;
331  if (xx > xhalfpitch) xxx = m_strip_pitch - xx;
332  int ix = static_cast<int>(xxx/deltax);
333  double fx = (xxx - ix*deltax) / deltax;
334 
335  int ix1 = ix+1;
336  int iy1 = iy+1;
337  double Ex00 = 0., Ex10 = 0., Ex01 = 0., Ex11 = 0.;
338  double Ey00 = 0., Ey10 = 0., Ey01 = 0., Ey11 = 0.;
339  //-------- pick up appropriate data file for m_VD-----------------------
340 
341  Ex00 = data.m_ExValue[ix][iy]; Ex10 = data.m_ExValue[ix1][iy];
342  Ex01 = data.m_ExValue[ix][iy1]; Ex11 = data.m_ExValue[ix1][iy1];
343  Ey00 = data.m_EyValue[ix][iy]; Ey10 = data.m_EyValue[ix1][iy];
344  Ey01 = data.m_EyValue[ix][iy1]; Ey11 = data.m_EyValue[ix1][iy1];
345 
346  //---------------- end of data bank search---
347  Ex = Ex00*(1.-fx)*(1.-fy) + Ex10*fx*(1.-fy)
348  + Ex01*(1.-fx)* fy + Ex11*fx* fy;
349  if (xx > xhalfpitch) Ex = -Ex;
350  Ey = Ey00*(1.-fx)*(1.-fy) + Ey10*fx*(1.-fy)
351  + Ey01*(1.-fx)* fy + Ey11*fx* fy;
352  }
353  //---------- case for uniform electric field ------------------------
354  else if (data.m_EFieldModel==UniformE) {
355  if (m_bulk_depth - y < data.m_depletion_depth && data.m_depletion_depth >0.) {
356  Ey = data.m_VB / data.m_depletion_depth;
357  } else {
358  Ey = 0.;
359  }
360  }
361  //---------- case for flat diode model ------------------------------
362  else if (data.m_EFieldModel==FlatDiodeModel) {
363  if (data.m_VB > std::abs(data.m_VD)) {
364  Ey = (data.m_VB+data.m_VD)/m_bulk_depth - 2.*data.m_VD*(m_bulk_depth-y)/(m_bulk_depth*m_bulk_depth);
365  } else {
366  if (m_bulk_depth - y < data.m_depletion_depth && data.m_depletion_depth >0.) {
367  double Emax = 2.* data.m_depletion_depth * std::abs(data.m_VD) / (m_bulk_depth*m_bulk_depth);
368  Ey = Emax*(1.-(m_bulk_depth-y)/data.m_depletion_depth);
369  } else {
370  Ey = 0.;
371  }
372  }
373  }
374 }
375 
378  // Return index for s_VFD0 and EFieldModel.
379  // If vdepl is out of range of s_VFD0, EFieldModel of FlatDiodeModel is returned.
380  size_t iVFD = 0;
381  if (data.m_VD<s_VFD0[0] || data.m_VD>s_VFD0[s_VFD0.size()-1]) { // Out of range
382  data.m_EFieldModel = FlatDiodeModel;
383  } else {
384  float dVFD2 = 0.;
385  for (iVFD=0; iVFD<s_VFD0.size()-1; iVFD++) {
386  dVFD2 = s_VFD0[iVFD+1] - data.m_VD;
387  if (dVFD2 >= 0.) break;
388  }
389  }
390  return iVFD;
391 }
392 
394 
396  // Loading Ramo potential and Electric field.
397  // In strip pitch directions :
398  // Ramo potential : 80 divisions (81 points) with 5 um intervals from 40-440 um.
399  // Electric field : 16 divisions (17 points) with 2.5 um intervals from 0-40 um.
400  // In sensor depth directions (for both potential and Electric field):
401  // 114 divisions (115 points) with 2.5 nm intervals for 285 um.
402 
403  TFile *hfile = new TFile(PathResolverFindCalibFile("SCT_Digitization/SCT_InducedChargedModel.root").c_str());
404 
405  // For Ramo Potential
406  TH2F* h_Potential_FDM = static_cast<TH2F*>(hfile->Get("Potential_FDM"));
407  TH2F* h_Potential_FEM = static_cast<TH2F*>(hfile->Get("Potential_FEM"));
408  // FDM case keeps only FDM potential.
409  // FEM case keeps both FDM and FEM potentials.
412  else m_PotentialValue.resize(FEMsolutions+1);
413  for (int ix=0; ix<NRamoPoints; ix++) {
414  for (int iy=0; iy<NDepthPoints; iy++) {
415  m_PotentialValue[FlatDiodeModel][ix][iy] = h_Potential_FDM->GetBinContent(ix+1, iy+1);
417  m_PotentialValue[FEMsolutions][ix][iy] = h_Potential_FEM->GetBinContent(ix+1, iy+1);
418  }
419  }
420  }
421  }
422  h_Potential_FDM->Delete();
423  h_Potential_FEM->Delete();
424 
426  size_t i=0;
427  m_ExValue.resize(s_VFD0str.size()+1);
428  m_EyValue.resize(s_VFD0str.size()+1);
429  for (const std::string& str : s_VFD0str) {
430  TH2F* h_Ex = static_cast<TH2F*>(hfile->Get(("Ex_FEM_"+str+"_0").c_str()));
431  TH2F* h_Ey = static_cast<TH2F*>(hfile->Get(("Ey_FEM_"+str+"_0").c_str()));
432  for (int ix=0; ix <NEFieldPoints; ix++){
433  for (int iy=0; iy<NDepthPoints; iy++) {
434  m_ExValue[i][ix][iy] = h_Ex->GetBinContent(ix+1, iy+1);
435  m_EyValue[i][ix][iy] = h_Ey->GetBinContent(ix+1, iy+1);
436  }
437  }
438  h_Ex->Delete();
439  h_Ey->Delete();
440  i++;
441  }
442 
443  TH2F* h_ExHV = static_cast<TH2F*>(hfile->Get("Ex_FEM_0_100"));
444  TH2F* h_EyHV = static_cast<TH2F*>(hfile->Get("Ey_FEM_0_100"));
445  for (int ix=0; ix <NEFieldPoints; ix++){
446  for (int iy=0; iy<NDepthPoints; iy++) {
447  m_ExValue[i][ix][iy] = h_ExHV->GetBinContent(ix+1, iy+1);
448  m_EyValue[i][ix][iy] = h_EyHV->GetBinContent(ix+1, iy+1);
449  }
450  }
451  h_ExHV->Delete();
452  h_EyHV->Delete();
453  }
454 
455  hfile->Close();
456 }
457 
459  if (data.m_EFieldModel!=FEMsolutions) return;
460 
461  size_t iVFD = getFEMIndex(data);
462  if (data.m_EFieldModel==FlatDiodeModel) {
463  ATH_MSG_INFO("Changed to Flat Diode Model since deplettion volage is out of range. ("
464  << s_VFD0[0] << " < m_VD < " << s_VFD0[s_VFD0.size()-1] << " is allowed.)");
465  return;
466  }
467 
468  const float dVFD1 = data.m_VD - s_VFD0[iVFD];
469  const float dVFD2 = s_VFD0[iVFD+1] - data.m_VD;
470 
471  const float scalingFactor = data.m_VB / 100.; // Reference voltage is 100 V
472 
473  const size_t iHV = s_VFD0.size();
474 
475  // For Electric field
476  for (int ix=0; ix <NEFieldPoints; ix++){
477  for (int iy=0; iy<NDepthPoints; iy++) {
478  float Ex1 = m_ExValue[iVFD ][ix][iy];
479  float Ex2 = m_ExValue[iVFD+1][ix][iy];
480  float ExHV = m_ExValue[iHV ][ix][iy];
481  data.m_ExValue[ix][iy] = (Ex1*dVFD2+Ex2*dVFD1)/(dVFD1+dVFD2);
482  data.m_ExValue[ix][iy] += ExHV*scalingFactor;
483 
484  float Ey1 = m_EyValue[iVFD ][ix][iy];
485  float Ey2 = m_EyValue[iVFD+1][ix][iy];
486  float EyHV = m_EyValue[iHV ][ix][iy];
487  data.m_EyValue[ix][iy] = (Ey1*dVFD2+Ey2*dVFD1)/(dVFD1+dVFD2);
488  data.m_EyValue[ix][iy] += EyHV*scalingFactor;
489  }
490  }
491 }
InducedChargeModel::EndStrip
@ EndStrip
Definition: InducedChargeModel.h:49
InducedChargeModel::NTransportSteps
@ NTransportSteps
Definition: InducedChargeModel.h:47
beamspotman.r
def r
Definition: beamspotman.py:676
InducedChargeModel::UniformE
@ UniformE
Definition: InducedChargeModel.h:46
python.SystemOfUnits.second
int second
Definition: SystemOfUnits.py:120
data
char data[hepevt_bytes_allocation_ATLAS]
Definition: HepEvt.cxx:11
python.SystemOfUnits.joule
int joule
Definition: SystemOfUnits.py:151
python.PhysicalConstants.k_Boltzmann
float k_Boltzmann
Definition: PhysicalConstants.py:114
python.SystemOfUnits.s
int s
Definition: SystemOfUnits.py:131
InducedChargeModel::SCT_InducedChargeModelData
Definition: InducedChargeModel.h:51
ATH_MSG_INFO
#define ATH_MSG_INFO(x)
Definition: AthMsgStreamMacros.h:31
AtlasFieldCacheCondObj
Definition: AtlasFieldCacheCondObj.h:19
InDetDD::SolidStateDetectorElementBase
Definition: SolidStateDetectorElementBase.h:132
InducedChargeModel::getEField
void getEField(const SCT_InducedChargeModelData &data, const double x, const double y, double &Ex, double &Ey) const
Definition: InducedChargeModel.cxx:308
python.SystemOfUnits.e_SI
float e_SI
Definition: SystemOfUnits.py:138
InducedChargeModel::m_EyValue
std::vector< std::array< std::array< double, NDepthPoints >, NEFieldPoints > > m_EyValue
Definition: InducedChargeModel.h:167
InducedChargeModel::s_VFD0
static const std::vector< float > s_VFD0
Definition: InducedChargeModel.h:176
InducedChargeModel::setWaferData
SCT_InducedChargeModelData * setWaferData(const float vdepl, const float vbias, const InDetDD::SolidStateDetectorElementBase *element, const AtlasFieldCacheCondObj *fieldCondObj, const ToolHandle< ISiliconConditionsTool > &siConditionsTool, CLHEP::HepRandomEngine *rndmEngine, const EventContext &ctx) const
Definition: InducedChargeModel.cxx:40
DMTest::P
P_v1 P
Definition: P.h:23
m_data
std::vector< T > m_data
Definition: TrackTruthMatchingBaseAlg.cxx:660
InducedChargeModel::EFieldModel
EFieldModel
Definition: InducedChargeModel.h:46
InducedChargeModel::m_y_origin_min
double m_y_origin_min
Definition: InducedChargeModel.h:155
InducedChargeModel.h
InducedChargeModel::setEField
void setEField(SCT_InducedChargeModelData &data) const
Definition: InducedChargeModel.cxx:458
Amg::y
@ y
Definition: GeoPrimitives.h:35
InducedChargeModel::m_mutex
std::mutex m_mutex
Definition: InducedChargeModel.h:174
InducedChargeModel::NEFieldPoints
@ NEFieldPoints
Definition: InducedChargeModel.h:48
InducedChargeModel::getVxVyD
bool getVxVyD(const SCT_InducedChargeModelData &data, const bool isElectron, const double x, const double y, double &vx, double &vy, double &D, const IdentifierHash hashId, const ToolHandle< ISiPropertiesTool > &siPropertiesTool, const EventContext &ctx) const
Definition: InducedChargeModel.cxx:231
ReadOfcFromCool.field
field
Definition: ReadOfcFromCool.py:48
InducedChargeModel::m_strip_pitch
double m_strip_pitch
Definition: InducedChargeModel.h:154
python.TrigEgammaMonitorHelper.TH2F
def TH2F(name, title, nxbins, bins_par2, bins_par3, bins_par4, bins_par5=None, bins_par6=None, path='', **kwargs)
Definition: TrigEgammaMonitorHelper.py:45
x
#define x
InDetDD::SolidStateDetectorElementBase::identifyHash
virtual IdentifierHash identifyHash() const override final
identifier hash (inline)
InducedChargeModel::StartStrip
@ StartStrip
Definition: InducedChargeModel.h:49
cm
const double cm
Definition: Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimParametrization/tools/FCAL_ChannelMap.cxx:25
InducedChargeModel::Offset
@ Offset
Definition: InducedChargeModel.h:49
Amg::z
@ z
Definition: GeoPrimitives.h:36
Run3DQTestingDriver.dq
dq
Definition: Run3DQTestingDriver.py:108
AtlasFieldCacheCondObj::getInitializedCache
void getInitializedCache(MagField::AtlasFieldCache &cache) const
get B field cache for evaluation as a function of 2-d or 3-d position.
Definition: AtlasFieldCacheCondObj.h:32
InducedChargeModel::m_transportTimeStep
double m_transportTimeStep
Definition: InducedChargeModel.h:149
InducedChargeModel::getFEMIndex
static size_t getFEMIndex(SCT_InducedChargeModelData &data)
Definition: InducedChargeModel.cxx:377
lumiFormat.i
int i
Definition: lumiFormat.py:85
InducedChargeModel::transport
void transport(const SCT_InducedChargeModelData &data, const bool isElectron, const double x0, const double y0, double *Q_m2, double *Q_m1, double *Q_00, double *Q_p1, double *Q_p2, const IdentifierHash hashId, const ToolHandle< ISiPropertiesTool > &siPropertiesTool, const EventContext &ctx) const
Definition: InducedChargeModel.cxx:104
ATH_MSG_DEBUG
#define ATH_MSG_DEBUG(x)
Definition: AthMsgStreamMacros.h:29
CaloNoise_fillDB.dt
dt
Definition: CaloNoise_fillDB.py:58
xAOD::EgammaHelpers::isElectron
bool isElectron(const xAOD::Egamma *eg)
is the object an electron (not Fwd)
Definition: EgammaxAODHelpers.cxx:12
Amg::x
@ x
Definition: GeoPrimitives.h:34
InducedChargeModel::m_ExValue
std::vector< std::array< std::array< double, NDepthPoints >, NEFieldPoints > > m_ExValue
Definition: InducedChargeModel.h:166
AthMessaging
Class to provide easy MsgStream access and capabilities.
Definition: AthMessaging.h:55
InducedChargeModel::holeTransport
void holeTransport(const SCT_InducedChargeModelData &data, const double x0, const double y0, double *Q_m2, double *Q_m1, double *Q_00, double *Q_p1, double *Q_p2, const IdentifierHash hashId, const ToolHandle< ISiPropertiesTool > &siPropertiesTool, const EventContext &ctx) const
Definition: InducedChargeModel.cxx:191
python.SystemOfUnits.volt
int volt
Definition: SystemOfUnits.py:204
dot.dot
def dot(G, fn, nodesToHighlight=[])
Definition: dot.py:5
InducedChargeModel::NDepthPoints
@ NDepthPoints
Definition: InducedChargeModel.h:48
InducedChargeModel::m_EFieldModel
EFieldModel m_EFieldModel
Definition: InducedChargeModel.h:148
InDetDD::SolidStateDetectorElementBase::globalPositionHit
HepGeom::Point3D< double > globalPositionHit(const HepGeom::Point3D< double > &simulationLocalPos) const
transform a hit local position into a global position (inline):
InducedChargeModel::s_e
static const double s_e
Definition: InducedChargeModel.h:145
VP1PartSpect::E
@ E
Definition: VP1PartSpectFlags.h:21
python.SystemOfUnits.tesla
int tesla
Definition: SystemOfUnits.py:228
InducedChargeModel::NRamoPoints
@ NRamoPoints
Definition: InducedChargeModel.h:48
Amg::Vector3D
Eigen::Matrix< double, 3, 1 > Vector3D
Definition: GeoPrimitives.h:47
InducedChargeModel::NStrips
@ NStrips
Definition: InducedChargeModel.h:49
InducedChargeModel::s_VFD0str
static const std::vector< std::string > s_VFD0str
Definition: InducedChargeModel.h:177
PathResolverFindCalibFile
std::string PathResolverFindCalibFile(const std::string &logical_file_name)
Definition: PathResolver.cxx:431
InducedChargeModel::m_transportTimeMax
double m_transportTimeMax
Definition: InducedChargeModel.h:150
python.PyAthena.v
v
Definition: PyAthena.py:154
y
#define y
InducedChargeModel::electronTransport
void electronTransport(const SCT_InducedChargeModelData &data, const double x0, const double y0, double *Q_m2, double *Q_m1, double *Q_00, double *Q_p1, double *Q_p2, const IdentifierHash hashId, const ToolHandle< ISiPropertiesTool > &siPropertiesTool, const EventContext &ctx) const
Definition: InducedChargeModel.cxx:211
python.SystemOfUnits.gauss
int gauss
Definition: SystemOfUnits.py:230
GlobalVariables.Emax
Emax
Definition: GlobalVariables.py:185
correlationModel::model
model
Definition: AsgElectronEfficiencyCorrectionTool.cxx:46
makeTRTBarrelCans.dx
tuple dx
Definition: makeTRTBarrelCans.py:20
MagField::AtlasFieldCache
Local cache for magnetic field (based on MagFieldServices/AtlasFieldSvcTLS.h)
Definition: AtlasFieldCache.h:43
InducedChargeModel::loadICMParameters
void loadICMParameters()
Definition: InducedChargeModel.cxx:395
MagField::AtlasFieldCache::getField
void getField(const double *ATH_RESTRICT xyz, double *ATH_RESTRICT bxyz, double *ATH_RESTRICT deriv=nullptr)
get B field value at given position xyz[3] is in mm, bxyz[3] is in kT if deriv[9] is given,...
Definition: AtlasFieldCache.cxx:42
str
Definition: BTagTrackIpAccessor.cxx:11
InducedChargeModel::InducedChargeModel
InducedChargeModel(size_t maxHash, EFieldModel model=FEMsolutions)
Definition: InducedChargeModel.cxx:31
IdentifierHash
This is a "hash" representation of an Identifier. This encodes a 32 bit index which can be used to lo...
Definition: IdentifierHash.h:25
InducedChargeModel::induced
double induced(const SCT_InducedChargeModelData &data, const int istrip, const double x, const double y) const
Definition: InducedChargeModel.cxx:276
InducedChargeModel::m_bulk_depth
double m_bulk_depth
Definition: InducedChargeModel.h:153
CaloNoise_fillDB.mu
mu
Definition: CaloNoise_fillDB.py:53
InducedChargeModel::FEMsolutions
@ FEMsolutions
Definition: InducedChargeModel.h:46
InducedChargeModel::m_PotentialValue
std::vector< std::array< std::array< double, NDepthPoints >, NRamoPoints > > m_PotentialValue
Definition: InducedChargeModel.h:165
defineDB.hfile
string hfile
Definition: JetTagCalibration/share/defineDB.py:37
InducedChargeModel::s_kB
static const double s_kB
Definition: InducedChargeModel.h:144
InducedChargeModel::FlatDiodeModel
@ FlatDiodeModel
Definition: InducedChargeModel.h:46