31 #include "G4AffineTransform.hh"
32 #include "G4NavigationHistory.hh"
34 #include "G4ThreeVector.hh"
35 #include "G4StepPoint.hh"
47 namespace BarrelPresampler {
50 : base_class(
name, pSvcLocator)
62 m_end_module[0]=(
m_mod[0][0]*
m_cmm+2*eps)+m_zminPS+eps;
63 for (
int i=1;
i<8;
i++) m_end_module[
i]=m_end_module[
i-1]+(
m_mod[
i][0]*
m_cmm+2*eps)+eps;
65 for (
int i=0;
i<8;
i++) {
72 m_first_cathod[0]=m_zminPS+
m_mod[0][5]*
m_cmm+m_cat_th/2.+2*eps;
73 for (
int i=1;
i<8;
i++) m_first_cathod[
i]=m_end_module[
i-1]+
m_mod[
i][5]*
m_cmm+m_cat_th/2.+2*eps;
80 for (
int i=0;
i<7;
i++) m_ncell_module[
i]=8;
87 for (
int i=0;
i<7;
i++) m_ngap_cell[
i]=(
int)((
m_mod[
i][1]+0.1)/m_ncell_module[
i]);
96 return StatusCode::SUCCESS;
111 LArG4Identifier Geometry::CalculateIdentifier(
const G4Step* a_step)
const
113 const static G4String fullPSName = (m_detectorName.empty()) ?
"LAr::Barrel::Presampler" : G4String(m_detectorName+
"::LAr::Barrel::Presampler");
114 const static G4String fullCryoName = (m_detectorName.empty()) ?
"LAr::TBBarrel::Cryostat::LAr" : G4String(m_detectorName+
"::LAr::TBBarrel::Cryostat::LAr");
115 const static G4String fullModuleName = (m_detectorName.empty()) ?
"LAr::Barrel::Presampler::Module" : G4String(m_detectorName+
"::LAr::Barrel::Presampler::Module");
118 const G4NavigationHistory* g4navigation = a_step->GetPreStepPoint()->GetTouchable()->GetHistory();
119 const G4int ndep = g4navigation->GetDepth();
124 for (G4int ii=0;ii<=ndep;ii++) {
126 const G4String& vname = g4navigation->GetVolume(ii)->GetName();
127 if (idep<0 && vname==fullPSName) idep=ii;
128 else if (!iactive && vname==fullModuleName) iactive=
true;
131 if (idep < 0) std::abort();
134 return this->CalculatePSActiveIdentifier( a_step , idep );
136 return this->CalculatePS_DMIdentifier( a_step , idep );
145 LArG4Identifier Geometry::CalculatePSActiveIdentifier(
const G4Step* a_step,
const G4int ind)
const
149 const G4ThreeVector
p = (a_step->GetPostStepPoint()->GetPosition() + a_step->GetPreStepPoint()->GetPosition()) * 0.5;
152 ATH_MSG_VERBOSE(
"Position of the step in the ATLAS frame (x,y,z) --> " <<
p.x() <<
" " <<
p.y() <<
" " <<
p.z());
153 ATH_MSG_VERBOSE(
"Eta and Phi in the atlas frame --> " <<
p.eta() <<
" " <<
p.phi());
158 const G4NavigationHistory* g4navigation = a_step->GetPreStepPoint()->GetTouchable()->GetHistory();
159 const G4ThreeVector ploc = g4navigation->GetTransform(ind).TransformPoint(
p);
161 const G4int zSide(this->determineZSide(
p.z()));
164 (void)this->findCell(currentCellData,ploc.x(),ploc.y(),ploc.z());
170 if(currentCellData.
phiBin < 0 ) currentCellData.
phiBin += 64;
180 << currentCellData.
phiBin;
183 ATH_MSG_VERBOSE(
"Here the identifier for the presampler ACTIVE ----> ");
196 LArG4Identifier Geometry::CalculatePS_DMIdentifier(
const G4Step* a_step,
const G4int ind)
const
220 const G4ThreeVector
p = (a_step->GetPostStepPoint()->GetPosition() + a_step->GetPreStepPoint()->GetPosition()) * 0.5;
223 ATH_MSG_VERBOSE(
"Position of the step in the ATLAS frame (x,y,z) --> " <<
p.x() <<
" " <<
p.y() <<
" " <<
p.z());
224 ATH_MSG_VERBOSE(
"Eta and Phi in the atlas frame --> " <<
p.eta() <<
" " <<
p.phi());
230 const G4NavigationHistory* g4navigation = a_step->GetPreStepPoint()->GetTouchable()->GetHistory();
231 const G4ThreeVector ploc = g4navigation->GetTransform(ind).TransformPoint(
p);
232 const G4double
radius=sqrt(ploc.x()*ploc.x() + ploc.y()*ploc.y());
235 const G4ThreeVector ploc2(ploc.x(),ploc.y(),ploc.z()+m_zpres+m_zminPS);
238 ATH_MSG_VERBOSE(
"Position of the step after traslation (x,y,z) --> " << ploc2.x() <<
" " << ploc2.y() <<
" " << ploc2.z());
239 ATH_MSG_VERBOSE(
"Eta and Phi after translation --> " << ploc2.eta() <<
" " << ploc2.phi());
243 const G4int zSide(this->determineZSide(
p.z()));
246 const G4double
phi = (ploc2.phi() < 0.) ? ploc2.phi()+2.*
M_PI : ploc2.phi();
247 const G4double
eta = ploc2.eta();
258 const G4int numberPhiMod = 32;
259 const G4double dphi = ( 2.*
M_PI ) / numberPhiMod;
260 const G4double inv_dphi = 1. / dphi;
261 const G4double PSModuleRmean = 1420 ;
263 const G4double Rcheck = PSModuleRmean /
cos(
phicheck);
266 currentCellData.
region = 3;
268 currentCellData.
region = 2;
271 const G4double detaDM = 0.1 ;
272 const G4double dphiDM = ( 2 *
M_PI ) / 64. ;
274 currentCellData.
phiBin = G4int(
phi * (1. / dphiDM) );
275 currentCellData.
etaBin = G4int(
eta * (1. / detaDM) );
280 if(currentCellData.
phiBin < 0 ) currentCellData.
phiBin += 64;
284 if ( currentCellData.
phiBin == 64 ) currentCellData.
phiBin = 0;
294 << currentCellData.
phiBin;
297 ATH_MSG_VERBOSE(
"Here the identifier for the presampler DEAD materials ---->");
333 bool Geometry::findCell(CalcData & currentCellData, G4double xloc,G4double yloc,G4double zloc)
const
336 currentCellData.sampling = 0;
337 currentCellData.region = 0;
340 G4double
phi = atan2(yloc,xloc);
342 const G4double z2=fabs(zloc+m_zpres+m_zminPS);
345 const G4int numberPhiBins = 64;
346 const G4double dphi = ( 2.*
M_PI ) / numberPhiBins;
347 const G4double inv_dphi = 1. / dphi;
349 currentCellData.phiBin = G4int(
phi * inv_dphi );
350 if (currentCellData.phiBin >63) currentCellData.phiBin=63;
351 if (currentCellData.phiBin <0) currentCellData.phiBin=0;
356 if (z2 < m_zminPS ) {
357 currentCellData.etaBin=0;
360 if (z2 > m_end_module[7]) {
361 currentCellData.etaBin=60;
366 currentCellData.module=0;
367 for (
int i=1;
i<8;
i++) {
368 if (m_first_cathod[
i]>=z2)
break;
369 currentCellData.module++;
371 if (currentCellData.module <0 || currentCellData.module > 7)
373 G4cerr <<
"LArBarrelPresampler: invalid module/hit " << currentCellData.module <<
" " << z2 << G4endl;
374 if (currentCellData.module<0) currentCellData.etaBin=0;
375 if (currentCellData.module>7) currentCellData.etaBin=60;
386 ATH_MSG_VERBOSE(
"sector number, dist along height " << isect <<
" " << currentCellData.dist);
387 ATH_MSG_VERBOSE(
"z2,module number,m_first_cathod " << z2 <<
" " << currentCellData.module <<
" "
388 << m_first_cathod[currentCellData.module]);
393 G4double deltaz=z2-(m_first_cathod[currentCellData.module]+currentCellData.dist*
tan(m_tilt[currentCellData.module]));
395 if (currentCellData.module>0) {
396 currentCellData.module=currentCellData.module-1;
397 deltaz=z2-(m_first_cathod[currentCellData.module]+currentCellData.dist*
tan(m_tilt[currentCellData.module]));
403 currentCellData.gap = ((
int)(deltaz/m_pitch[currentCellData.module]));
406 ATH_MSG_VERBOSE(
"deltaz from first cathode,gap number " << deltaz <<
" " << currentCellData.gap);
410 currentCellData.etaBin= currentCellData.gap/m_ngap_cell[currentCellData.module];
414 if (currentCellData.etaBin >= m_ncell_module[currentCellData.module]) currentCellData.etaBin=m_ncell_module[currentCellData.module]-1;
416 for (
int i=0;
i<currentCellData.module;
i++) currentCellData.etaBin=currentCellData.etaBin+m_ncell_module[
i];
421 if (currentCellData.etaBin < 0 || currentCellData.etaBin > 60) {
422 ATH_MSG_WARNING(
"LArBarrelPresamplerGeometry::findCell etaBin outside range " << currentCellData.etaBin);
426 const G4double zmiddle=m_first_cathod[currentCellData.module]+((
double)(currentCellData.gap+0.5))*m_pitch[currentCellData.module];
431 currentCellData.xElec=currentCellData.dist*
cos(m_tilt[currentCellData.module])+(z2-zmiddle)*
sin(m_tilt[currentCellData.module]);
432 currentCellData.distElec=(z2-zmiddle)*
cos(m_tilt[currentCellData.module]) - currentCellData.dist*
sin(m_tilt[currentCellData.module]);
434 ATH_MSG_VERBOSE(
"zmiddle,xloc,yloc " << zmiddle <<
" " << currentCellData.distElec <<
" " << currentCellData.xElec);
438 if (fabs(currentCellData.dist)>m_halfThickLAr) {