ATLAS Offline Software
Public Member Functions | Protected Attributes | Private Member Functions | Private Attributes | List of all members
LArBarrelCalculator Class Reference

#include <LArBarrelCalculator.h>

Inheritance diagram for LArBarrelCalculator:
Collaboration diagram for LArBarrelCalculator:

Public Member Functions

 LArBarrelCalculator (const std::string &name, ISvcLocator *pSvcLocator)
 
virtual StatusCode initialize () override final
 
virtual StatusCode finalize () override final
 
 LArBarrelCalculator (const LArBarrelCalculator &)=delete
 
LArBarrelCalculatoroperator= (const LArBarrelCalculator &)=delete
 
void initializeForSDCreation () override final
 
virtual G4bool Process (const G4Step *a_step, std::vector< LArHitData > &hdata) const override final
 
virtual G4bool isInTime (G4double hitTime) const override final
 
virtual G4float OOTcut () const override final
 

Protected Attributes

Gaudi::Property< bool > m_BirksLaw {this, "BirksLaw", true}
 
Gaudi::Property< double > m_Birksk {this, "Birksk", 0.05832}
 
Gaudi::Property< double > m_OOTcut {this, "OOTcut", 300*CLHEP::ns}
 

Private Member Functions

G4bool FiducialCuts (G4double, G4double, G4double) const
 
void InitHV ()
 
double ScaleHV (double, double, double, double) const
 

Private Attributes

ServiceHandle< ILArBarrelGeometrym_geometry {this, "GeometryCalculator", "LArBarrelGeometry"}
 
const AccMapm_accmap {nullptr}
 
std::unique_ptr< MapEtam_etamap1 {}
 
std::unique_ptr< MapEtam_etamap2 {}
 
std::unique_ptr< MapEtam_etamap3 {}
 
Gaudi::Property< bool > m_IflCur {this, "EMBCurr", true}
 
Gaudi::Property< bool > m_IflMapTrans {this, "EMBEtaTrans", true}
 
Gaudi::Property< bool > m_IflXtalk {this, "EMBXtalk", true}
 
Gaudi::Property< double > m_dstep {this, "EMBdstep", .2*CLHEP::mm}
 
const LArG4BirksLawm_birksLaw {nullptr}
 
Gaudi::Property< bool > m_doHV {this, "EMBHVEnable", false}
 
double m_etaMaxBarrel {0.}
 
double m_zMinBarrel {0.}
 
double m_zMaxBarrel {0.}
 
double m_rMinAccordion {0.}
 
double m_rMaxAccordion {0.}
 
double m_ThickAbs {0.}
 
double m_ThickEle {0.}
 
int m_NCellTot {0}
 
int m_NCellMax {0}
 
bool m_testbeam {false}
 
double m_hv [2][1024][7][2] {}
 

Detailed Description

Definition at line 32 of file LArBarrelCalculator.h.

Constructor & Destructor Documentation

◆ LArBarrelCalculator() [1/2]

LArBarrelCalculator::LArBarrelCalculator ( const std::string &  name,
ISvcLocator *  pSvcLocator 
)

Definition at line 47 of file LArBarrelCalculator.cxx.

48  : LArCalculatorSvcImp(name, pSvcLocator)
49 {
50 }

◆ LArBarrelCalculator() [2/2]

LArBarrelCalculator::LArBarrelCalculator ( const LArBarrelCalculator )
delete

Member Function Documentation

◆ FiducialCuts()

G4bool LArBarrelCalculator::FiducialCuts ( G4double  radloc,
G4double  zloc,
G4double  etaloc 
) const
private

Definition at line 618 of file LArBarrelCalculator.cxx.

619 {
620  if (radloc < m_rMinAccordion ||
621  radloc > m_rMaxAccordion ||
622  etaloc > m_etaMaxBarrel ||
623  zloc < m_zMinBarrel ||
624  zloc > m_zMaxBarrel) return false;
625  else return true;
626 }

◆ finalize()

StatusCode LArBarrelCalculator::finalize ( )
finaloverridevirtual

Definition at line 145 of file LArBarrelCalculator.cxx.

146 {
147  if (m_BirksLaw) delete m_birksLaw;
148  return StatusCode::SUCCESS;
149 }

◆ InitHV()

void LArBarrelCalculator::InitHV ( )
private

Definition at line 632 of file LArBarrelCalculator.cxx.

633 {
634  ATH_MSG_DEBUG(" **** in LArBarrelCalculator::InitHV() ");
635 
636 
637 
638 
639  float defaultHvVal=2000.;
640  ATH_MSG_DEBUG(" defaultHvVal " << defaultHvVal);
641  for (int ipm=0;ipm<2;ipm++) {
642  for (int ielec=0;ielec<1024;ielec++) {
643  for (int ieta=0;ieta<7;ieta++) {
644  for (int iside=0;iside<2;iside++) {
645  m_hv[ipm][ielec][ieta][iside] = defaultHvVal;
646  }
647  }
648  }
649  }
650 
651  if (m_doHV) {
652  // initialize services
653  SmartIF<StoreGateSvc> pDetStore{Gaudi::svcLocator()->service("DetectorStore")};
654  if (!pDetStore) {
655  ATH_MSG_WARNING("LArBarrelCalculator::InitHV() unable to get Detector Store! Use default HV values.");
656  return;
657  }
658 
659  // get EMBHV Manager
660  const LArHVManager *manager = nullptr;
661  if (pDetStore->retrieve(manager)==StatusCode::SUCCESS) {
662  const EMBHVManager& hvManager=manager->getEMBHVManager();
663  const EMBHVManager::EMBHVData hvdata = hvManager.getDataSim();
664  ATH_MSG_DEBUG(" got HV Manager ");
665  // loop over HV modules
666  for (unsigned int iSide=0;iSide<2;iSide++) {
667  for (unsigned int iPhi=0;iPhi<16;iPhi++) {
668  for (unsigned int iSector=0;iSector<2;iSector++) {
669  for (unsigned int iEta=0;iEta<7;iEta++) {
670  const EMBHVModule& hvMod = hvManager.getHVModule(iSide,iEta+1,iPhi,iSector);
671  for (unsigned int ielec=0;ielec<32;ielec++) {
672  const EMBHVElectrode& electrode = hvMod.getElectrode(ielec);
673  unsigned jElec = ielec+32*iSector+64*iPhi;
674  for (unsigned int iGap=0;iGap<2;iGap++) {
675  double hv = hvdata.voltage (electrode, iGap);
676  ATH_MSG_DEBUG(" iSide,jElec,iEta,iGap,hv " << iSide << " " << jElec << " " << iEta << " " << iGap << " " << hv);
677  if (hv>-999.) m_hv[iSide][jElec][iEta][iGap] = hv;
678  }
679  }
680  }
681  }
682  }
683  }
684  }
685  else {
686  ATH_MSG_WARNING(" Unable to find HV Manager ");
687  }
688  }
689 
690  return;
691 }

◆ initialize()

StatusCode LArBarrelCalculator::initialize ( )
finaloverridevirtual

Definition at line 52 of file LArBarrelCalculator.cxx.

53 {
54  ATH_MSG_DEBUG("LArBarrelCalculator: Beginning initialization ");
55  if (m_BirksLaw) {
56  const double Birks_LAr_density = 1.396;
57  m_birksLaw = new LArG4BirksLaw(Birks_LAr_density,m_Birksk);
58  }
59 
60  // Access source of detector parameters.
62 
63  // Get the out-of-time cut from the detector parameters routine.
64  m_OOTcut = parameters->GetValue("LArExpHallOutOfTimeCut"); //FIXME should be done via configurables
65  ATH_MSG_DEBUG("**** OutOfTime cut " << m_OOTcut);
66 
67  // Main Barrel parameters
68  // All the UNITS are implicitly the GEANT4 ONES e.g. mm, rad, etc ...
69  m_etaMaxBarrel = parameters->GetValue("LArEMBMaxEtaAcceptance");
70  m_zMinBarrel = parameters->GetValue("LArEMBfiducialMothZmin");
71  m_zMaxBarrel = parameters->GetValue("LArEMBfiducialMothZmax");
72  m_NCellMax = (int) (parameters->GetValue("LArEMBnoOFPhysPhiCell"));
73  m_rMinAccordion = parameters->GetValue("LArEMBRadiusInnerAccordion");
74  m_rMaxAccordion = parameters->GetValue("LArEMBFiducialRmax");
75 
76  // absorbers and electrodes thickness
77  m_ThickAbs = 0.5*( parameters->GetValue("LArEMBThickAbsGlue")
78  +parameters->GetValue("LArEMBThickAbsIron")
79  +parameters->GetValue("LArEMBThickAbsLead"));
80 
81  G4double check = 0.5*( parameters->GetValue("LArEMBThinAbsGlue")
82  +parameters->GetValue("LArEMBThinAbsIron")
83  +parameters->GetValue("LArEMBThinAbsLead"));
84 
85  if (fabs(check-m_ThickAbs)>0.001)
86  {
87  ATH_MSG_WARNING("Thin and Thick Abs have difference thickness!");
88  }
89 
90  m_ThickEle= 0.5*( parameters->GetValue("LArEMBThickElecCopper")
91  +parameters->GetValue("LArEMBThickElecKapton"));
92 
93  // === GU 11/06/2003 total number of cells in phi
94  // to distinguish 1 module (testbeam case) from full Atlas
95  m_NCellTot = (int) (parameters->GetValue("LArEMBnoOFPhysPhiCell"));
96  if (m_NCellTot != 1024)
97  {
98  m_NCellMax=1024;
99  m_testbeam=true;
100  }
101  // ===
102  // access geometry computation class
103  ATH_CHECK(m_geometry.retrieve());
104 
105  // access current maps if required
106  if (m_IflCur)
107  {
108  ATH_MSG_DEBUG(" LArBarrelCalculator: start reading of current maps");
110  m_etamap1 = std::make_unique<MapEta>(1);
111  m_etamap2 = std::make_unique<MapEta>(2);
112  if (m_IflMapTrans) m_etamap3 = std::make_unique<MapEta>(3);
113  ATH_MSG_DEBUG(" LArBarrelCalculator: end of reading current maps");
114  }
115  // Initialize HV values
116  this->InitHV();
117 
118  ATH_MSG_DEBUG(" LArBarrelCalculator: s_NCellMax " << m_NCellMax);
119  ATH_MSG_DEBUG(" LArBarrelCalculator: s_NCellTot " << m_NCellTot);
120  ATH_MSG_DEBUG(" LArBarrelCalculator: s_rMinAccordion " << m_rMinAccordion);
121  ATH_MSG_DEBUG(" LArBarrelCalculator: s_rMaxAccordion " << m_rMaxAccordion);
122  ATH_MSG_DEBUG(" LArBarrelCalculator: s_zMinBarrel " << m_zMinBarrel);
123  ATH_MSG_DEBUG(" LArBarrelCalculator: s_zMaxBarrel " << m_zMaxBarrel);
124  ATH_MSG_DEBUG(" LArBarrelCalculator: s_etaMaxBarrel " << m_etaMaxBarrel);
125  ATH_MSG_DEBUG(" LArBarrelCalculator: s_ThickAbs " << m_ThickAbs);
126  ATH_MSG_DEBUG(" LArBarrelCalculator: s_ThickEle " << m_ThickEle);
127  if(m_IflCur) ATH_MSG_DEBUG(" LArBarrelCalculator: Deposited Energy dE/dX Corrected ==> CURRENT Option ON");
128  else ATH_MSG_DEBUG(" LArBarrelCalculator: Crude Deposited Energy dE/dX NO CURRENT option");
129  if (m_IflCur && m_IflMapTrans) ATH_MSG_DEBUG(" LArBarrelCalculator: Compute effect of E field around eta=0.8 ");
130  else ATH_MSG_DEBUG(" LArBarrelCalculator: Ignore effect of E field around eta=0.8 ");
131  if(m_BirksLaw)
132  {
133  ATH_MSG_DEBUG(" LArBarrelCalculator: Birks' law ON ");
134  ATH_MSG_DEBUG(" LArBarrelCalculator: parameter k " << m_birksLaw->k());
135  }
136  else
137  {
138  ATH_MSG_DEBUG(" LArBarrelCalculator: Birks' law OFF");
139  }
140 
141  return StatusCode::SUCCESS;
142 }

◆ initializeForSDCreation()

void LArBarrelCalculator::initializeForSDCreation ( )
inlinefinaloverridevirtual

Reimplemented from LArCalculatorSvcImp.

Definition at line 43 of file LArBarrelCalculator.h.

43 { m_geometry->initializeForSDCreation(); };

◆ isInTime()

virtual G4bool LArBarrelCalculator::isInTime ( G4double  hitTime) const
inlinefinaloverridevirtual

Definition at line 48 of file LArBarrelCalculator.h.

49  {
50  return !(std::fabs(hitTime) > m_OOTcut);
51  }

◆ OOTcut()

virtual G4float LArBarrelCalculator::OOTcut ( ) const
inlinefinaloverridevirtual

Definition at line 56 of file LArBarrelCalculator.h.

56 { return m_OOTcut; }

◆ operator=()

LArBarrelCalculator& LArBarrelCalculator::operator= ( const LArBarrelCalculator )
delete

◆ Process()

G4bool LArBarrelCalculator::Process ( const G4Step *  a_step,
std::vector< LArHitData > &  hdata 
) const
finaloverridevirtual

Definition at line 152 of file LArBarrelCalculator.cxx.

153 {
154 
155  hdata.clear();
156 
157  LArG4Identifier identifier2;
158  LArG4Identifier identifier_xt1;
159  LArG4Identifier identifier_xt2;
160 
161  // check the Step content is non trivial
162  G4double thisStepEnergyDeposit = step->GetTotalEnergyDeposit() * step->GetTrack()->GetWeight();
163  G4double thisStepLength = step->GetStepLength() / Units::mm;
164 
165 #ifdef DEBUGSTEP
166  ATH_MSG_DEBUG("****** LArBarrelCalculator: Step energy,length "
167  << thisStepEnergyDeposit << " " << thisStepLength);
168 #endif
169 
170  if(thisStepEnergyDeposit <= 0. || thisStepLength <= 0.)
171  {
172 #ifdef DEBUGSTEP2
173  ATH_MSG_DEBUG(" Invalid hit trivial content");
174 #endif
175  return false;
176  }
177 
178  // Get Step geometrical parameters (first and end)
179  G4StepPoint *thisStepPoint = step->GetPreStepPoint();
180  G4StepPoint *thisStepBackPoint = step->GetPostStepPoint();
181  G4ThreeVector startPoint = thisStepPoint->GetPosition();
182  G4ThreeVector endPoint = thisStepBackPoint->GetPosition();
183 
184 #ifdef DEBUGSTEP
185  ATH_MSG_DEBUG(" Beginning step position "
186  << startPoint.x() << " " << startPoint.y() << " " << startPoint.z());
187 #endif
188 
189  // find zside from volume name
190  G4int zSide = 1;
191  if (!m_testbeam) {
192  const G4NavigationHistory* g4navigation = thisStepPoint->GetTouchable()->GetHistory();
193  G4int ndep = g4navigation->GetDepth();
194  for (G4int ii=0;ii<=ndep;ii++) {
195  G4VPhysicalVolume* v1 = g4navigation->GetVolume(ii);
196  G4String vname = v1->GetName();
197  if ( vname.find("NegPhysical") != std::string::npos) zSide=-1;
198  }
199  }
200 
201 
202  // BACK directly into the LOCAL half_Barrel Z > 0. PART (mother "stac_phys1")
203 
204  const G4AffineTransform transformation =
205  thisStepPoint->GetTouchable()->GetHistory()->GetTopTransform();
206 
207  G4ThreeVector startPointinLocal =
208  transformation.TransformPoint(startPoint);
209  G4ThreeVector endPointinLocal =
210  transformation.TransformPoint (endPoint);
211 
212 #ifdef DEBUGSTEP
213  ATH_MSG_DEBUG(" Local beginning step position "
214  << startPointinLocal.x() << " " << startPointinLocal.y() << " "
215  << startPointinLocal.z());
216  ATH_MSG_DEBUG(" Local end step position "
217  << endPointinLocal.x() << " " << endPointinLocal.y() << " "
218  << endPointinLocal.z());
219 #endif
220 
221  G4double energy = step->GetTotalEnergyDeposit() * step->GetTrack()->GetWeight(); // Despite the name, this is only ionization.
222 
223  if (m_BirksLaw) {
224  const double EField = 10.; // kV/cm, assume constant for now
225  const double wholeStepLengthCm = step->GetStepLength() / CLHEP::cm;
226  energy = (*m_birksLaw)(energy, wholeStepLengthCm, EField);
227  }
228 
229 
230  // compute number of sub steps
231  // = 1 if no charge collection
232  // otherwise 200 microns (dstep) division
233 
234  G4int nsub_step=1;
235  if (m_IflCur) nsub_step=(int) (thisStepLength/m_dstep) + 1;
236  // delta is fraction of step between two sub steps
237  G4double delta=1./((double) nsub_step);
238 #ifdef DEBUGSTEP
239  ATH_MSG_DEBUG(" nsub_step,delta " << nsub_step << " " << delta);
240 #endif
241 
242  energy /= ((float) (nsub_step)); // Divide into substeps.
243 
244  // loop over sub steps
245 
246  for (G4int i=0;i<nsub_step;i++) {
247 
248  bool Xtalk = false;
249 
250  G4double fraction=(((G4double) i)+0.5)*delta;
251  G4ThreeVector subinLocal=startPointinLocal*(1-fraction) + endPointinLocal*fraction;
252  G4double xloc = subinLocal.x();
253  G4double yloc = subinLocal.y();
254  G4double zloc = subinLocal.z();
255  G4double etaloc = subinLocal.pseudoRapidity();
256  G4double philoc = subinLocal.phi();
257  if(philoc<0.) philoc = philoc + 2.*M_PI;
258  G4double radloc = sqrt( xloc*xloc + yloc*yloc );
259 #ifdef DEBUGSTEP
260  ATH_MSG_DEBUG(" local position sub_step "
261  << xloc << " " << yloc << " " << zloc);
262 #endif
263 
264  // apply fiducial cuts
265  if (!FiducialCuts(radloc,zloc,etaloc)) {
266 #ifdef DEBUGSTEP
267  ATH_MSG_DEBUG("LArBarrelCalculator: outside fiducial cuts");
268 #endif
269  continue;
270  }
271  LArG4::Barrel::CalcData currentCellData;
272  m_geometry->findCell(currentCellData,xloc,yloc,zloc,radloc,etaloc,philoc,m_IflCur);
273 
274  if (currentCellData.cellID == 0)
275  {
276 #ifdef DEBUGSTEP
277  ATH_MSG_DEBUG("LArBarrelCalculator: Invalid hit CELLID == 0 ");
278  ATH_MSG_DEBUG("x,y,z local " << xloc << " " << yloc << " " << zloc);
279  ATH_MSG_DEBUG("radius " << radloc <<" etaloc "<< etaloc << " philoc " << philoc);
280 #endif
281  continue;
282  }
283  G4int region = currentCellData.region;
284  G4int etaBin = currentCellData.etaBin;
285  G4int phiBin = currentCellData.phiBin;
286  G4int sampling = currentCellData.sampling;
287 
288  if( zSide == -1 )
289  {
290  // following code for an Y-axis rotation to define Z < 0. Barrel part
291  if( sampling == 1 && region ==0 )
292  {
293  phiBin = 31 - phiBin;
294  if(phiBin < 0 ) phiBin += 64;
295  }
296  if( sampling == 1 && region ==1 )
297  {
298  phiBin = 127 - phiBin;
299  if(phiBin < 0 ) phiBin += 256;
300  }
301  if( sampling >= 2 )
302  {
303  phiBin = 127 - phiBin;
304  if(phiBin < 0 ) phiBin += 256;
305  }
306  }
307 #ifdef DEBUGSTEP
308  ATH_MSG_DEBUG(" region,side,sampling,eta,phi " << region << " " << zSide << " "
309  << sampling << " " << etaBin << " " << phiBin);
310  ATH_MSG_DEBUG(" distance to electrode,abs " << currentCellData.distElec << " "
311  << currentCellData.distAbs);
312  ATH_MSG_DEBUG(" local coordinates " << currentCellData.x0 << " " << currentCellData.y0);
313 #endif
314 
315  if (std::fabs(currentCellData.distElec)< m_ThickEle ||
316  std::fabs(currentCellData.distAbs) < m_ThickAbs) {
317 #ifdef DEBUGSTEP
318  ATH_MSG_DEBUG(" hit in absorber or electrode radius:" << radloc);
319 #endif
320  continue;
321  }
322 
323  identifier2.clear();
324  identifier2 << 4 // LArCalorimeter
325  << 1 // LArEM
326  << zSide
327  << sampling
328  << region
329  << etaBin
330  << phiBin;
331 
332  // time computation is not necessarily correct for test beam
333  // G4double time;
334  // if (m_testbeam)
335  // {
336  // time=0.;
337  // }
338  // else
339  // {
340  // const G4double tof = thisStepPoint->GetGlobalTime() / Units::ns;
341  // time = tof - ( thisStepPoint->GetPosition().mag()/CLHEP::c_light ) / Units::ns;
342  // }
343 
344  const G4double time = (m_testbeam)? 0.0 :
345  ( thisStepPoint->GetGlobalTime() - ( thisStepPoint->GetPosition().mag()/CLHEP::c_light ) ) / Units::ns;
346 
347 #ifdef DEBUGSTEP
348  ATH_MSG_DEBUG(" Energy for sub step " << energy);
349 #endif
350 
351  G4double Current;
352  G4double Current_xt1,Current_xt2;
353  Current_xt1=0;
354  Current_xt2=0;
355  if (!m_IflCur) {
356  // no charge collection Current=E from Geant
357  Current=energy;
358  }
359  else {
360  // full charge collection
361  G4double xmap,ymap;
362  G4int nfold = currentCellData.nfold;
363  G4double x0=currentCellData.x0;
364  G4double y0=currentCellData.y0;
365  if (x0<1500 || x0>1960 || y0>30 || y0<-30) {
366  ATH_MSG_WARNING("weird x0,y0 " << x0 << " " << y0);
367  }
368 #ifdef DEBUGSTEP
369  G4double rr = sqrt(x0*x0+y0*y0);
370  ATH_MSG_DEBUG(" radius,rad0 " << radloc << " " << rr);
371 #endif
372  const CurrMap* mapOfCurrent(nullptr);
373  bool UseFold=false;
374  // Are we close to a fold ? (fold 0 has some pathology)
375  if ((x0 > m_accmap->GetXmin(nfold) || nfold==0) &&
376  x0 < m_accmap->GetXmax(nfold) &&
377  y0 > m_accmap->GetYmin(nfold) &&
378  y0 < m_accmap->GetYmax(nfold) &&
379  (nfold != 13 || x0 < 1957.5) &&
380  (nfold>0 || x0 < 1504.6) ) {
381  xmap=x0;
382  ymap=y0;
383 #ifdef DEBUGSTEP
384  ATH_MSG_DEBUG(" Map for fold xmap,ymap " << nfold << " " << xmap << " " << ymap);
385 #endif
386  UseFold=true;
387  G4int sampMap = currentCellData.sampMap;
388  G4int etaMap = currentCellData.etaMap;
389  mapOfCurrent = m_accmap->GetMap(nfold,region,sampMap,etaMap);
390  // catch problem to find map
391  if (!mapOfCurrent) {
392  ATH_MSG_WARNING(" Problem to access map fold = " << nfold);
393  ATH_MSG_WARNING(" region,sampling,eta,fold " << region << " " << sampMap << " "
394  << etaMap << " " << nfold);
395  return false;
396  }
397  }
398  else {
399  G4int n;
400  G4int nstraight = currentCellData.nstraight;
401  if (nstraight%2==0) n=22;
402  else n=21;
403  mapOfCurrent = m_accmap->GetMap(n,region,sampling,etaBin);
404  // catch problem to find map
405  if (!mapOfCurrent) {
406  ATH_MSG_WARNING(" Problem to access map straight = " << nstraight);
407  return false;
408  }
409  xmap = currentCellData.xl;
410  ymap = currentCellData.distElec;
411  // special case for first straight section, which is shorter
412  if (nstraight==0) xmap = 0.5*(xmap+1.);
413 #ifdef DEBUGSTEP
414  ATH_MSG_DEBUG(" Map for straight xl,delec " << xmap << " " << ymap);
415 #endif
416  } // fold or straight
417  double gap;
418  double current0,current1,current2;
419  // get current for elementary charge
420  mapOfCurrent->GetAll(xmap,ymap,&gap,&current0,&current1,&current2);
421  G4double gap2=std::fabs(currentCellData.distElec)+std::fabs(currentCellData.distAbs)
423 
424  // in which HV cell are we ?
425  int ipm,ielec,ieta,iside;
426  if (zSide==1) ipm=1; // A side
427  else ipm=0; // C side
428  ielec=currentCellData.phiGap;
429  if (zSide==-1) {
430  ielec = 511 - ielec;
431  if(ielec < 0 ) ielec += 1024;
432  }
433  ieta=((int) (etaloc*(1./0.2)));
434  if (ieta>6) ieta=6; //part 1.4 to 1.475 is same HV as 1.2 to 1.4
435  iside=0; // phi lower than electrode 0, 1 for phi higher than electrode
436  if ((currentCellData.distElec>0 && zSide==1)
437  || (currentCellData.distElec<0 && zSide==-1) ) iside=1;
438 
439  // HV extrapolation
440  double current;
441  double hv=m_hv[ipm][ielec][ieta][iside];
442  if (hv>1995.) current=current0;
443  else if (hv>5.) current=ScaleHV(hv,current0,current1,current2);
444  else current=0.;
445 
446  // extrapolation for non nominal gap (allows to include sagging effect)
447  // i ~ (gap/gap2)**1.3
448  // gap = nominal gap from current map
449  // gap2 = real gap from geometry
450  // dgap = gap/gap2 -1
451  // at first order i ~ (1+1.3*dgap)
452  double dgap=0;
453  if (gap>1e-3 && gap2 >1e-3) dgap=gap/gap2-1;
454  current = current * (1. + 1.3*dgap);
455 
456 #ifdef DEBUGSTEP
457  ATH_MSG_DEBUG(" elementary current " << current0);
458  ATH_MSG_DEBUG(" Gap from Map/calculator " << gap << " " << gap2);
459 #endif
460 
461  Current = energy*current;
462 
463  // check if pathology...
464  if (!UseFold && std::fabs(currentCellData.distElec)>2.1 && current0 < 0.1) {
465  ATH_MSG_WARNING(" xl,distEle " << currentCellData.xl << " "
466  << currentCellData.distElec
467  << " str number " << currentCellData.nstraight
468  << " sampling,eta " << sampling << " " << etaBin << " "
469  << " current/E " << current0);
470  }
471 
472  // compute effect around transition at eta=0.8 from electric field
473  // read from map weighting factor ( E**1.3 = E*E**0.3)
474  // in x +-9mm around transition (assume symmetry around 0)
475  // in y distance from electrode in gap
476  bool InTrans=false;
477  //double Current_test=Current;
478  if (m_IflMapTrans) {
479  float etaTrans=0.8;
480  if (fabs(etaloc-etaTrans) < 0.025) {
481  double x=std::fabs(zloc-radloc*sinh(etaTrans))/cosh(etaTrans);
482  double y=std::fabs(currentCellData.distElec);
483  if ( x < m_etamap3->Xmax() ) {
484  double resp;
485  m_etamap3->GetData0(x,y,&resp);
486  Current = Current*resp;
487  InTrans=true;
488  }
489  } // eta = 0.8 +- 0.025
490  } // if m_IflMapTrans
491 
492  // simulate cross-talk effects and transitions in eta between cells
493  // (only in region 0, samplings 1 and 2)
494 
495  if (region==0) {
496 
497  double resp,xt0,xt1,xt2,deta;
498  if (sampling==1) {
499  deta=etaloc-0.003125*((double)etaBin+0.5);
500  m_etamap1->GetData(std::fabs(deta),std::fabs(currentCellData.distElec),
501  &resp,&xt0,&xt1,&xt2);
502 #ifdef DEBUGSTEP
503  ATH_MSG_DEBUG("hit in strip etaloc,etaBin,deta,delec,resp,xt0,xt1,xt2 "
504  << etaloc << " " << etaBin << " " << deta*1000 << " " << currentCellData.distElec
505  << " " << resp << " " << xt0 << " " << xt1 << " " << xt2);
506 #endif
507  if (!InTrans) Current = Current*resp;
508  if (m_IflXtalk && etaBin > 1 && etaBin < 446) {
509  Xtalk=true;
510  if (deta>0) {
511  identifier_xt1.clear();
512  identifier_xt1 << 4 << 1 << zSide << sampling << region << (etaBin+1) << phiBin;
513  Current_xt1 = Current*xt1;
514  identifier_xt2.clear();
515  identifier_xt2 << 4 << 1 << zSide << sampling << region << (etaBin-1) << phiBin;
516  Current_xt2 = Current*xt2;
517  }
518  else {
519  identifier_xt1.clear();
520  identifier_xt1 << 4 << 1 << zSide << sampling << region << (etaBin-1) << phiBin;
521  Current_xt1 = Current*xt1;
522  identifier_xt2.clear();
523  identifier_xt2 << 4 << 1 << zSide << sampling << region << (etaBin+1) << phiBin;
524  Current_xt2 = Current*xt2;
525  }
526  Current = Current*xt0;
527  } // m_IflXtalk = true
528  }
529  else if (sampling==2 && !InTrans) {
530  deta=etaloc-0.025*((double)etaBin+0.5);
531  m_etamap2->GetData0(std::fabs(deta),std::fabs(currentCellData.distElec),
532  &resp);
533 #ifdef DEBUGSTEP
534  ATH_MSG_DEBUG("hit in middle etaloc,etaBin,deta,delec,resp,xt0,xt1,xt2 "
535  << etaloc << " " << etaBin << " " << deta*1000 << " " << currentCellData.distElec
536  << " " << resp);
537 #endif
538  Current = Current*resp;
539  } // sampling =1 or 2
540 
541  } // region=0
542 
543  } // switch for current simulation
544 
545 
546  // check if we have a new hit in a different cell, or if we add this substep
547  // to an already existing hit
548  G4bool found=false;
549  for (unsigned int i=0; i<hdata.size(); i++) {
550  if (hdata[i].id==identifier2) {
551  hdata[i].energy += Current;
552  hdata[i].time += time*Current;
553  found=true;
554  break;
555  }
556  } // loop over hits
557  if (!found) {
558  LArHitData newdata = {identifier2, time*Current, Current};
559  hdata.push_back(newdata);
560  } // hit was not existing before
561 
562  if (Xtalk) {
563  for (unsigned int i=0; i<hdata.size(); i++) {
564  if (hdata[i].id==identifier_xt1) {
565  hdata[i].energy += Current_xt1;
566  hdata[i].time += time*Current_xt1;
567  found=true;
568  break;
569  }
570  } // loop over hits
571  if (!found) {
572  LArHitData newdata = {identifier_xt1, time*Current_xt1, Current_xt1};
573  hdata.push_back(newdata);
574  }
575  found=false;
576  for (unsigned int i=0; i<hdata.size(); i++) {
577  if (hdata[i].id==identifier_xt2) {
578  hdata[i].energy += Current_xt2;
579  hdata[i].time += time*Current_xt2;
580  found=true;
581  break;
582  }
583  } // loop over hits
584  if (!found) {
585  LArHitData newdata = {identifier_xt2, time*Current_xt2, Current_xt2};
586  hdata.push_back(newdata);
587  }
588  } // Xtalk true
589 
590 
591  } // *** End of loop over sub steps
592 
593 #ifdef DEBUGSTEP
594  ATH_MSG_DEBUG("Number of hits for this step " << m_nhits << " "
595  << hdata.size());
596 #endif
597 
598  // finalise time computations
599  for (unsigned int i=0;i<hdata.size();i++) {
600  if (std::fabs(hdata[i].energy)>1e-6) hdata[i].time=hdata[i].time/hdata[i].energy;
601  else hdata[i].time=0.;
602 #ifdef DEBUGSTEP
603  ATH_MSG_DEBUG("Hit Energy/Time "
604  << hdata[i].energy << " " << hdata[i].time);
605 #endif
606  }
607 
608  if (hdata.size()>0) return true;
609  return false;
610 }

◆ ScaleHV()

double LArBarrelCalculator::ScaleHV ( double  hv,
double  curr0,
double  curr1,
double  curr2 
) const
private

Definition at line 700 of file LArBarrelCalculator.cxx.

701 {
702  double b;
703  double resp=0.;
704  if (hv>1000.) {
705  if (std::fabs(curr1)>1e-6) {
706  double x=curr0/curr1;
707  if (x>1e-6) {
708  b=log(x)*(1./(log(2000.)-log(1000.)));
709  resp = curr0*exp(b*(log(hv)-log(2000.)));
710  }
711  }
712  }
713  else if (hv>400.) {
714  if (std::fabs(curr2)>1e-6) {
715  double x=curr1/curr2;
716  if (x>1e-6) {
717  b=log(x)*(1./(log(1000.)-log(400.)));
718  resp = curr1*exp(b*(log(hv)-log(1000.)));
719  }
720  }
721  }
722  else {
723  if (std::fabs(curr2)>1e-6) {
724  b=0.57;
725  resp = curr2*exp(b*(log(hv)-log(400.)));
726  }
727  }
728  return resp;
729 }

Member Data Documentation

◆ m_accmap

const AccMap* LArBarrelCalculator::m_accmap {nullptr}
private

Definition at line 61 of file LArBarrelCalculator.h.

◆ m_Birksk

Gaudi::Property<double> LArCalculatorSvcImp::m_Birksk {this, "Birksk", 0.05832}
protectedinherited

Definition at line 27 of file LArCalculatorSvcImp.h.

◆ m_BirksLaw

Gaudi::Property<bool> LArCalculatorSvcImp::m_BirksLaw {this, "BirksLaw", true}
protectedinherited

Definition at line 23 of file LArCalculatorSvcImp.h.

◆ m_birksLaw

const LArG4BirksLaw* LArBarrelCalculator::m_birksLaw {nullptr}
private

Definition at line 73 of file LArBarrelCalculator.h.

◆ m_doHV

Gaudi::Property<bool> LArBarrelCalculator::m_doHV {this, "EMBHVEnable", false}
private

Definition at line 74 of file LArBarrelCalculator.h.

◆ m_dstep

Gaudi::Property<double> LArBarrelCalculator::m_dstep {this, "EMBdstep", .2*CLHEP::mm}
private

Definition at line 71 of file LArBarrelCalculator.h.

◆ m_etamap1

std::unique_ptr<MapEta> LArBarrelCalculator::m_etamap1 {}
private

Definition at line 62 of file LArBarrelCalculator.h.

◆ m_etamap2

std::unique_ptr<MapEta> LArBarrelCalculator::m_etamap2 {}
private

Definition at line 63 of file LArBarrelCalculator.h.

◆ m_etamap3

std::unique_ptr<MapEta> LArBarrelCalculator::m_etamap3 {}
private

Definition at line 64 of file LArBarrelCalculator.h.

◆ m_etaMaxBarrel

double LArBarrelCalculator::m_etaMaxBarrel {0.}
private

Definition at line 77 of file LArBarrelCalculator.h.

◆ m_geometry

ServiceHandle<ILArBarrelGeometry> LArBarrelCalculator::m_geometry {this, "GeometryCalculator", "LArBarrelGeometry"}
private

Definition at line 60 of file LArBarrelCalculator.h.

◆ m_hv

double LArBarrelCalculator::m_hv[2][1024][7][2] {}
private

Definition at line 100 of file LArBarrelCalculator.h.

◆ m_IflCur

Gaudi::Property<bool> LArBarrelCalculator::m_IflCur {this, "EMBCurr", true}
private

Definition at line 67 of file LArBarrelCalculator.h.

◆ m_IflMapTrans

Gaudi::Property<bool> LArBarrelCalculator::m_IflMapTrans {this, "EMBEtaTrans", true}
private

Definition at line 68 of file LArBarrelCalculator.h.

◆ m_IflXtalk

Gaudi::Property<bool> LArBarrelCalculator::m_IflXtalk {this, "EMBXtalk", true}
private

Definition at line 69 of file LArBarrelCalculator.h.

◆ m_NCellMax

int LArBarrelCalculator::m_NCellMax {0}
private

Definition at line 88 of file LArBarrelCalculator.h.

◆ m_NCellTot

int LArBarrelCalculator::m_NCellTot {0}
private

Definition at line 87 of file LArBarrelCalculator.h.

◆ m_OOTcut

Gaudi::Property<double> LArCalculatorSvcImp::m_OOTcut {this, "OOTcut", 300*CLHEP::ns}
protectedinherited

Definition at line 30 of file LArCalculatorSvcImp.h.

◆ m_rMaxAccordion

double LArBarrelCalculator::m_rMaxAccordion {0.}
private

Definition at line 82 of file LArBarrelCalculator.h.

◆ m_rMinAccordion

double LArBarrelCalculator::m_rMinAccordion {0.}
private

Definition at line 81 of file LArBarrelCalculator.h.

◆ m_testbeam

bool LArBarrelCalculator::m_testbeam {false}
private

Definition at line 91 of file LArBarrelCalculator.h.

◆ m_ThickAbs

double LArBarrelCalculator::m_ThickAbs {0.}
private

Definition at line 84 of file LArBarrelCalculator.h.

◆ m_ThickEle

double LArBarrelCalculator::m_ThickEle {0.}
private

Definition at line 85 of file LArBarrelCalculator.h.

◆ m_zMaxBarrel

double LArBarrelCalculator::m_zMaxBarrel {0.}
private

Definition at line 79 of file LArBarrelCalculator.h.

◆ m_zMinBarrel

double LArBarrelCalculator::m_zMinBarrel {0.}
private

Definition at line 78 of file LArBarrelCalculator.h.


The documentation for this class was generated from the following files:
LArBarrelCalculator::m_IflXtalk
Gaudi::Property< bool > m_IflXtalk
Definition: LArBarrelCalculator.h:69
AllowedVariables::e
e
Definition: AsgElectronSelectorTool.cxx:37
LArBarrelCalculator::m_geometry
ServiceHandle< ILArBarrelGeometry > m_geometry
Definition: LArBarrelCalculator.h:60
fillPileUpNoiseLumi.current
current
Definition: fillPileUpNoiseLumi.py:52
LArBarrelCalculator::m_etamap2
std::unique_ptr< MapEta > m_etamap2
Definition: LArBarrelCalculator.h:63
EMBHVElectrode
Definition: EMBHVElectrode.h:15
LArG4Identifier
Definition: LArG4Identifier.h:121
LArBarrelCalculator::m_ThickEle
double m_ThickEle
Definition: LArBarrelCalculator.h:85
python.SystemOfUnits.mm
float mm
Definition: SystemOfUnits.py:98
LArCalculatorSvcImp::m_Birksk
Gaudi::Property< double > m_Birksk
Definition: LArCalculatorSvcImp.h:27
LArBarrelCalculator::m_NCellTot
int m_NCellTot
Definition: LArBarrelCalculator.h:87
LArGeo::VDetectorParameters
Definition: VDetectorParameters.h:29
LArBarrelCalculator::m_accmap
const AccMap * m_accmap
Definition: LArBarrelCalculator.h:61
LArBarrelCalculator::m_etamap3
std::unique_ptr< MapEta > m_etamap3
Definition: LArBarrelCalculator.h:64
LArBarrelCalculator::m_etamap1
std::unique_ptr< MapEta > m_etamap1
Definition: LArBarrelCalculator.h:62
gridSubmitIDTPM.check
check
Definition: gridSubmitIDTPM.py:45
LArBarrelCalculator::m_IflMapTrans
Gaudi::Property< bool > m_IflMapTrans
Definition: LArBarrelCalculator.h:68
LArBarrelCalculator::m_etaMaxBarrel
double m_etaMaxBarrel
Definition: LArBarrelCalculator.h:77
LArBarrelCalculator::m_zMinBarrel
double m_zMinBarrel
Definition: LArBarrelCalculator.h:78
LArBarrelCalculator::m_hv
double m_hv[2][1024][7][2]
Definition: LArBarrelCalculator.h:100
LArG4Identifier::clear
void clear()
M_PI
#define M_PI
Definition: ActiveFraction.h:11
LArHitData
Definition: ILArCalculatorSvc.h:23
LArBarrelCalculator::m_ThickAbs
double m_ThickAbs
Definition: LArBarrelCalculator.h:84
LArCalculatorSvcImp::LArCalculatorSvcImp
LArCalculatorSvcImp(const std::string &name, ISvcLocator *pSvcLocator)
Definition: LArCalculatorSvcImp.h:15
LArG4::Barrel::CalcData
Definition: ILArBarrelGeometry.h:26
LArBarrelCalculator::m_zMaxBarrel
double m_zMaxBarrel
Definition: LArBarrelCalculator.h:79
LArG4::Barrel::CalcData::phiGap
G4int phiGap
Definition: ILArBarrelGeometry.h:33
drawFromPickle.exp
exp
Definition: drawFromPickle.py:36
LArG4::Barrel::CalcData::x0
G4double x0
Definition: ILArBarrelGeometry.h:39
x
#define x
LArBarrelCalculator::InitHV
void InitHV()
Definition: LArBarrelCalculator.cxx:632
LArG4::Barrel::CalcData::xl
G4double xl
Definition: ILArBarrelGeometry.h:38
CaloSwCorrections.gap
def gap(flags, cells_name, *args, **kw)
Definition: CaloSwCorrections.py:212
LArBarrelCalculator::m_dstep
Gaudi::Property< double > m_dstep
Definition: LArBarrelCalculator.h:71
LArBarrelCalculator::m_NCellMax
int m_NCellMax
Definition: LArBarrelCalculator.h:88
LArBarrelCalculator::m_testbeam
bool m_testbeam
Definition: LArBarrelCalculator.h:91
LArG4BirksLaw::k
double k() const
Definition: LArG4BirksLaw.h:19
EMBHVManager::getHVModule
const EMBHVModule & getHVModule(unsigned int iSide, unsigned int iEta, unsigned int iPhi, unsigned int iSector) const
Definition: EMBHVManager.cxx:194
cm
const double cm
Definition: Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimParametrization/tools/FCAL_ChannelMap.cxx:25
LArG4::Barrel::CalcData::cellID
int cellID
Definition: ILArBarrelGeometry.h:27
LArBarrelCalculator::ScaleHV
double ScaleHV(double, double, double, double) const
Definition: LArBarrelCalculator.cxx:700
xAOD::etaBin
setSAddress setEtaMS setDirPhiMS setDirZMS setBarrelRadius setEndcapAlpha setEndcapRadius setInterceptInner setEtaMap etaBin
Definition: L2StandAloneMuon_v1.cxx:149
ParticleGun_FastCalo_ChargeFlip_Config.energy
energy
Definition: ParticleGun_FastCalo_ChargeFlip_Config.py:78
CurrMap
Definition: CurrMap.h:10
lumiFormat.i
int i
Definition: lumiFormat.py:85
EMBHVManager::EMBHVData::voltage
double voltage(const EMBHVElectrode &electrode, const int &iGap) const
Definition: EMBHVManager.cxx:130
LArBarrelCalculator::m_birksLaw
const LArG4BirksLaw * m_birksLaw
Definition: LArBarrelCalculator.h:73
beamspotman.n
n
Definition: beamspotman.py:729
LArG4::Barrel::CalcData::sampling
G4int sampling
Definition: ILArBarrelGeometry.h:28
ATH_MSG_DEBUG
#define ATH_MSG_DEBUG(x)
Definition: AthMsgStreamMacros.h:29
xAOD::etaMap
setSAddress setEtaMS setDirPhiMS setDirZMS setBarrelRadius setEndcapAlpha setEndcapRadius setInterceptInner etaMap
Definition: L2StandAloneMuon_v1.cxx:145
LArG4::Barrel::CalcData::phiBin
G4int phiBin
Definition: ILArBarrelGeometry.h:31
LArBarrelCalculator::m_rMinAccordion
double m_rMinAccordion
Definition: LArBarrelCalculator.h:81
ATH_CHECK
#define ATH_CHECK
Definition: AthCheckMacros.h:40
LArG4::Barrel::CalcData::nfold
G4int nfold
Definition: ILArBarrelGeometry.h:35
xAOD::double
double
Definition: CompositeParticle_v1.cxx:159
LArG4::Barrel::CalcData::etaMap
G4int etaMap
Definition: ILArBarrelGeometry.h:42
MuonR4::SegmentFit::ParamDefs::x0
@ x0
EMBHVManager::getDataSim
EMBHVData getDataSim() const
Definition: EMBHVManager.cxx:324
LArG4::Barrel::CalcData::sampMap
G4int sampMap
Definition: ILArBarrelGeometry.h:41
LArBarrelCalculator::m_rMaxAccordion
double m_rMaxAccordion
Definition: LArBarrelCalculator.h:82
LArG4::Barrel::CalcData::region
G4int region
Definition: ILArBarrelGeometry.h:29
LArCalculatorSvcImp::m_BirksLaw
Gaudi::Property< bool > m_BirksLaw
Definition: LArCalculatorSvcImp.h:23
MuonR4::SegmentFit::ParamDefs::y0
@ y0
EMBHVManager::EMBHVData
Definition: EMBHVManager.h:39
name
std::string name
Definition: Control/AthContainers/Root/debug.cxx:240
xAOD::phiBin
setSAddress setEtaMS setDirPhiMS setDirZMS setBarrelRadius setEndcapAlpha setEndcapRadius setPhiMap phiBin
Definition: L2StandAloneMuon_v2.cxx:145
plotBeamSpotMon.b
b
Definition: plotBeamSpotMon.py:76
LArBarrelCalculator::FiducialCuts
G4bool FiducialCuts(G4double, G4double, G4double) const
Definition: LArBarrelCalculator.cxx:618
checkTriggerxAOD.found
found
Definition: checkTriggerxAOD.py:328
EMBHVManager
This class provides direct access to information on the HV electrodes within the barrels....
Definition: EMBHVManager.h:36
python.PhysicalConstants.c_light
float c_light
Definition: PhysicalConstants.py:73
AccMap::GetMap
const CurrMap * GetMap(int ifold, int ielecregion) const
Definition: AccMap.cxx:92
Trk::iPhi
@ iPhi
Definition: ParamDefs.h:47
LArG4::Barrel::CalcData::etaBin
G4int etaBin
Definition: ILArBarrelGeometry.h:30
LArG4::Barrel::CalcData::distElec
G4double distElec
Definition: ILArBarrelGeometry.h:36
LArG4::Barrel::CalcData::nstraight
G4int nstraight
Definition: ILArBarrelGeometry.h:34
LArG4::Barrel::CalcData::y0
G4double y0
Definition: ILArBarrelGeometry.h:40
LArG4BirksLaw
Definition: LArG4BirksLaw.h:8
AccMap::GetAccMap
static const AccMap * GetAccMap()
Definition: AccMap.cxx:72
LArHVManager
This class provides access to the High Voltage throughout the LAr. High voltage conditions can also b...
Definition: LArHVManager.h:24
y
#define y
python.CaloAddPedShiftConfig.int
int
Definition: CaloAddPedShiftConfig.py:45
CaloSwCorrections.time
def time(flags, cells_name, *args, **kw)
Definition: CaloSwCorrections.py:242
AccMap::GetYmin
float GetYmin(int ifold) const
Definition: AccMap.h:31
ATH_MSG_WARNING
#define ATH_MSG_WARNING(x)
Definition: AthMsgStreamMacros.h:32
python.CaloCondTools.log
log
Definition: CaloCondTools.py:20
LArCellBinning.step
step
Definition: LArCellBinning.py:158
physics_parameters.parameters
parameters
Definition: physics_parameters.py:144
EMBHVModule::getElectrode
const EMBHVElectrode & getElectrode(unsigned int iElectrode) const
Definition: EMBHVModule.cxx:69
python.Logging.manager
manager
Definition: PhysicsAnalysis/D3PDTools/AnaAlgorithm/python/Logging.py:92
rr
const boost::regex rr(r_r)
hitTime
float hitTime(const AFP_SIDSimHit &hit)
Definition: AFP_SIDSimHit.h:39
EMBHVModule
Describes one HV Module within the EMB.
Definition: EMBHVModule.h:20
LArBarrelCalculator::m_doHV
Gaudi::Property< bool > m_doHV
Definition: LArBarrelCalculator.h:74
LArG4::Barrel::CalcData::distAbs
G4double distAbs
Definition: ILArBarrelGeometry.h:37
xAOD::iEta
setScale setgFexType iEta
Definition: gFexJetRoI_v1.cxx:77
LArGeo::VDetectorParameters::GetInstance
static const VDetectorParameters * GetInstance()
Definition: VDetectorParameters.cxx:29
AccMap::GetXmin
float GetXmin(int ifold) const
Definition: AccMap.h:25
LArCalculatorSvcImp::m_OOTcut
Gaudi::Property< double > m_OOTcut
Definition: LArCalculatorSvcImp.h:30
LArBarrelCalculator::m_IflCur
Gaudi::Property< bool > m_IflCur
Definition: LArBarrelCalculator.h:67
python.SystemOfUnits.ns
float ns
Definition: SystemOfUnits.py:146
python.LArMinBiasAlgConfig.float
float
Definition: LArMinBiasAlgConfig.py:65