ATLAS Offline Software
Loading...
Searching...
No Matches
LArG4::Barrel::Geometry Class Reference

#include <LArBarrelGeometry.h>

Inheritance diagram for LArG4::Barrel::Geometry:
Collaboration diagram for LArG4::Barrel::Geometry:

Public Member Functions

 Geometry (const std::string &name, ISvcLocator *pSvcLocator)
virtual ~Geometry ()=default
virtual StatusCode initialize () override final
virtual StatusCode finalize () override final
virtual void initializeForSDCreation () override final
virtual LArG4Identifier CalculateIdentifier (const G4Step *) const override final
virtual LArG4Identifier CalculateSuperResolutionIdentifier (const G4Step *a_step) const override final
virtual void findCell (CalcData &currentCellData, const double &x, const double &y, const double &z, const double &r, const double &eta, const double &phi, const bool detail) const override final

Protected Member Functions

 Geometry ()

Private Member Functions

LArG4Identifier CalculateECAMIdentifier (const G4Step *, const G4int indEcam, const bool inSTAC=true, int zside=1) const
bool CheckLArIdentifier (int sampling, int region, int eta, int phi) const
bool CheckDMIdentifier (int type, int sampling, int region, int eta, int phi) const
double Distance_Ele (const double &x, const double &y, const int &PhiC, int &Num_Straight, const int &Num_Coude, double &xl) const
double Distance_Abs (const double &x, const double &y, const int &nabs, const int &Num_Straight, const int &Num_Coude) const
G4int SampSeg (G4double, G4double, G4double, G4int &, G4int &, G4int &, G4int &, G4int &) const
void GetRphi ()
 phi vs r of first absorber in nominal geometry
G4double Phi0 (G4double) const
G4int PhiGap (const double &, const double &, const double &) const

Private Attributes

Gaudi::Property< std::string > m_detectorName {this, "DetectorName", "LArMgr"}
Gaudi::Property< bool > m_testbeam {this, "TestBeam", false}
G4String m_ecamName
double m_rMinAccordion {0.}
double m_rMaxAccordion {0.}
double m_zMinBarrel {0.}
double m_zMaxBarrel {0.}
double m_zMaxBarrelDMMargin {0.}
double m_etaMaxBarrel {0.}
int m_NCellTot {0}
int m_NCellMax {0}
int m_Nbrt {0}
int m_Nbrt1 {0}
double m_gam0 {0.}
double m_rint_eleFib {0.}
double * m_rc {nullptr}
double * m_phic {nullptr}
double * m_xc {nullptr}
double * m_yc {nullptr}
double * m_delta {nullptr}
int m_parity {0}
const LArCoudeElectrodesm_coudeelec {nullptr}
const LArCoudeAbsorbersm_coudeabs {nullptr}
const LArStraightElectrodesm_electrode {nullptr}
const LArStraightAbsorbersm_absorber {nullptr}
bool m_iflSAG {false}
G4int m_NRphi {0}
G4double m_Rmin {0.}
G4double m_Rmax {0.}
G4double m_Rphi [5000] = {0}
G4double m_dR {0.}

Detailed Description

Definition at line 36 of file LArBarrelGeometry.h.

Constructor & Destructor Documentation

◆ Geometry() [1/2]

LArG4::Barrel::Geometry::Geometry ( const std::string & name,
ISvcLocator * pSvcLocator )

Definition at line 43 of file LArBarrelGeometry.cxx.

44 : base_class(name, pSvcLocator)
45 {
46 }

◆ ~Geometry()

virtual LArG4::Barrel::Geometry::~Geometry ( )
virtualdefault

◆ Geometry() [2/2]

LArG4::Barrel::Geometry::Geometry ( )
protected

Member Function Documentation

◆ CalculateECAMIdentifier()

LArG4Identifier LArG4::Barrel::Geometry::CalculateECAMIdentifier ( const G4Step * a_step,
const G4int indEcam,
const bool inSTAC = true,
int zside = 1 ) const
private

Definition at line 1120 of file LArBarrelGeometry.cxx.

1121 {
1122
1123 LArG4Identifier result = LArG4Identifier();
1124
1125 // Get all the information about the step
1126
1127 const G4StepPoint *thisStepPoint = a_step->GetPreStepPoint();
1128 const G4StepPoint *thisStepBackPoint = a_step->GetPostStepPoint();
1129 const G4ThreeVector startPoint = thisStepPoint->GetPosition();
1130 const G4ThreeVector endPoint = thisStepBackPoint->GetPosition();
1131 const G4ThreeVector p = (thisStepPoint->GetPosition() + thisStepBackPoint->GetPosition()) * 0.5;
1132
1133#ifdef DEBUGHITS
1134 ATH_MSG_VERBOSE("Position of the step in the ATLAS frame (x,y,z) --> " << p.x() << " " << p.y() << " " << p.z());
1135 ATH_MSG_VERBOSE("Eta and Phi in the ATLAS frame --> " << p.eta() << " " << p.phi());
1136#endif
1137
1138 // BACK directly into the LOCAL half_Barrel. All the variables in this LOCAL framework get the SUFFIX Zpos
1139
1140 const G4NavigationHistory* g4navigation = thisStepPoint->GetTouchable()->GetHistory();
1141 const G4AffineTransform transformation = g4navigation->GetTransform(indECAM);
1142 const G4ThreeVector startPointinLocal = transformation.TransformPoint(startPoint);
1143 const G4ThreeVector endPointinLocal = transformation.TransformPoint (endPoint);
1144 const G4ThreeVector midinLocal = (startPointinLocal+endPointinLocal)*0.5;
1145
1146#ifdef DEBUGHITS
1147 ATH_MSG_VERBOSE("Position of the step in the LOCAL frame (x,y,z) --> " << midinLocal.x() << " " << midinLocal.y() << " " << midinLocal.z());
1148 ATH_MSG_VERBOSE("Eta and Phi of the step in LOCAL frame --> " << midinLocal.eta() << " " << midinLocal.phi());
1149#endif
1150
1151 // coordinates in the local frame
1152
1153 const G4double xZpos = midinLocal.x();
1154 const G4double yZpos = midinLocal.y();
1155 const G4double zZpos = midinLocal.z();
1156 const G4double etaZpos = midinLocal.pseudoRapidity();
1157 const G4double phiZpos = (midinLocal.phi()<0.) ? midinLocal.phi() + 2.*M_PI : midinLocal.phi();
1158 const G4double radius2Zpos = xZpos*xZpos + yZpos*yZpos;
1159 const G4double radiusZpos = sqrt(radius2Zpos);
1160
1161 CalcData currentCellData;
1162 if (m_testbeam) {
1163 currentCellData.zSide = 1;
1164 }
1165 else {
1166 currentCellData.zSide = zside;
1167 }
1168
1169 // Check if the hit is in the fiducial range and in the STAC volume
1170 // if yes this is active or inactive material
1171
1172 if (inSTAC && radiusZpos >=m_rMinAccordion && radiusZpos <= m_rMaxAccordion &&
1173 zZpos <= m_zMaxBarrel && zZpos >= m_zMinBarrel && etaZpos <= m_etaMaxBarrel) {
1174
1175#ifdef DEBUGHITS
1176 ATH_MSG_VERBOSE("This hit is in the STAC volume !!!!! ");
1177#endif
1178
1179 // DETERMINATION of currentCellData.cellID,
1180 // currentCellData.zSide, currentCellData.sampling,
1181 // currentCellData.phiBin, currentCellData.etaBin,
1182 // m_stackNumID
1183 const bool MapDetail(false);
1184 this->findCell( currentCellData, xZpos, yZpos, zZpos, radiusZpos, etaZpos, phiZpos, MapDetail );
1185
1186 // adjust phi in the negative half barrel frame
1187
1188 if( currentCellData.zSide == -1 )
1189 {
1190 if( currentCellData.sampling == 1 && currentCellData.region ==0 )
1191 {
1192 currentCellData.phiBin = 31 - currentCellData.phiBin;
1193 if(currentCellData.phiBin < 0 ) currentCellData.phiBin += 64;
1194 }
1195 if( currentCellData.sampling == 1 && currentCellData.region ==1 )
1196 {
1197 currentCellData.phiBin = 127 - currentCellData.phiBin;
1198 if(currentCellData.phiBin < 0 ) currentCellData.phiBin += 256;
1199 }
1200 if( currentCellData.sampling >= 2 )
1201 {
1202 currentCellData.phiBin = 127 - currentCellData.phiBin;
1203 if(currentCellData.phiBin < 0 ) currentCellData.phiBin += 256;
1204 }
1205 }
1206
1207 // there are few hundred microns between the 4mm nominal beginning of the active region
1208 // and the real beginning of the first strip at eta=0.025/8
1209 // to avoid inactive energy with strip=0 assign this to strip=1
1210 if (currentCellData.sampling==1 && currentCellData.region==0 && currentCellData.etaBin==0) {
1211 currentCellData.etaBin=1;
1212 }
1213
1214 result << 4 // LArCalorimeter
1215 << 1 // LArEM
1216 << currentCellData.zSide
1217 << currentCellData.sampling
1218 << currentCellData.region
1219 << currentCellData.etaBin
1220 << currentCellData.phiBin;
1221
1222#ifdef DEBUGHITS
1223 ATH_MSG_VERBOSE("Here the identifier for the barrel ACTIVE ----> ");
1224 ATH_MSG_VERBOSE("eta in local frame --> " << etaZpos);
1225 ATH_MSG_VERBOSE("currentCellData.zSide ----> " << currentCellData.zSide);
1226 ATH_MSG_VERBOSE("currentCellData.sampling ----> " << currentCellData.sampling);
1227 ATH_MSG_VERBOSE("currentCellData.region ----> " << currentCellData.region);
1228 ATH_MSG_VERBOSE("currentCellData.etaBin ----> " << currentCellData.etaBin);
1229 ATH_MSG_VERBOSE("currentCellData.phiBin ----> " << currentCellData.phiBin);
1230 ATH_MSG_VERBOSE("And also etafirst ----> " << thisStepPoint->GetPosition().pseudoRapidity());
1231#endif
1232
1233 // if (!Geometry::CheckLArIdentifier(currentCellData.sampling,currentCellData.region,currentCellData.etaBin,currentCellData.phiBin)) {
1234 // ATH_MSG_ERROR(" ** Bad LAr identifier " << currentCellData.sampling << " " << currentCellData.region << " "
1235 // << currentCellData.etaBin << " " << currentCellData.phiBin);
1236 // ATH_MSG_ERROR(" x,y,z,eta,phi " << xZpos << " " << yZpos << " " << zZpos
1237 // << " " << radiusZpos << " " << etaZpos << " " << phiZpos);
1238 // }
1239
1240
1241 }
1242 // hits in dead material part
1243 else {
1244
1245 G4int sampling=0;
1246 G4int region=0;
1247 const G4int numDeadPhiBins = 64;
1248 double abs_eta = std::fabs(etaZpos);
1249 const double DM1EtaWidth = 0.1 ;
1250 const double DM1PhiWidth = 2.*M_PI / numDeadPhiBins ;
1251 currentCellData.etaBin = (G4int) ( abs_eta * (1./DM1EtaWidth) ) ;
1252 currentCellData.phiBin = (G4int) (phiZpos/ DM1PhiWidth );
1253 G4int type=1;
1254 // protect against rounding error for phi ~2pi
1255 if (currentCellData.phiBin==numDeadPhiBins) currentCellData.phiBin=currentCellData.phiBin-1;
1256
1257 // adjust phi for negative half barrel
1258
1259 if ( currentCellData.zSide == -1 ) {
1260 currentCellData.phiBin = 31 - currentCellData.phiBin;
1261 if (currentCellData.phiBin < 0 ) currentCellData.phiBin +=64 ;
1262 }
1263
1264 // material in front of the active accordion
1265 if ( radiusZpos < m_rMinAccordion ) {
1266 sampling =1 ;
1267 region = 3 ;
1268 if (currentCellData.etaBin > 14) currentCellData.etaBin=14;
1269
1270#ifdef DEBUGHITS
1271 ATH_MSG_VERBOSE("This hit is in the ECAM volume in front of the accordion (DEAD MATERIAL) !!!!! ");
1272#endif
1273
1274 } else if (radiusZpos >= m_rMaxAccordion){ // material behind the active accordion
1275 sampling = 2;
1276
1277 if (abs_eta < 1.0 ) {
1278 region = 0 ;
1279#ifdef DEBUGHITS
1280 ATH_MSG_VERBOSE("This hit is in the ECAM volume behind accordion (DEAD MATERIAL 0) !!!!! ");
1281#endif
1282 } else if ( abs_eta >= 1.0 && abs_eta < 1.5) {
1283 region = 2;
1284 currentCellData.etaBin = currentCellData.etaBin - 10; // to have etabin between 0 and 4
1285#ifdef DEBUGHITS
1286 ATH_MSG_VERBOSE("This hit is in the ECAM volume behind accordion (DEAD MATERIAL 2) !!!!! ");
1287#endif
1288 } else {
1289 ATH_MSG_ERROR(" LArBarrelGeometry: hit behind accordion at eta>1.5 !!! ");
1290 region = 2;
1291 currentCellData.etaBin = 4;
1292 }
1293
1294 } else if (zZpos <= m_zMinBarrel) { // inactive matter between two EMB halves
1295 type=2;
1296 region=0;
1297 const G4int phisave=currentCellData.phiBin;
1298 const G4bool MapDetail(false);
1299 this->findCell( currentCellData, xZpos, yZpos, zZpos, radiusZpos, etaZpos, phiZpos, MapDetail );
1300 sampling = currentCellData.sampling; // sampling as in normal definition
1301 currentCellData.etaBin=0;
1302 currentCellData.phiBin=phisave;
1303
1304 } else if (zZpos >= m_zMaxBarrel - m_zMaxBarrelDMMargin || abs_eta >= 1.40) { // inactive matter between EMB and scintillator including some margin for error
1305 if (abs_eta >1.0 && abs_eta < 1.5) {
1306 sampling=2;
1307 region=2;
1308 currentCellData.etaBin = currentCellData.etaBin - 10; // to have etabin between 0 and 4
1309 } else if (abs_eta < 1.6) {
1310 sampling=1;
1311 region=4;
1312 currentCellData.etaBin=0;
1313 } else {
1314 ATH_MSG_ERROR(" LArBarrelGeometry: hit at eta>1.6 !!! ");
1315 sampling=1;
1316 region=4;
1317 currentCellData.etaBin=0;
1318 }
1319 } else {
1320 if (!m_testbeam) {
1321 const G4double thisStepEnergyDeposit = a_step->GetTotalEnergyDeposit() * a_step->GetTrack()->GetWeight();
1322 std::ostringstream dmLog;
1323 dmLog << "LArBarrelGeometry: cannot find region for DM hit..." << std::endl;
1324 dmLog << "LArBarrelGeometry: cannot find region for DM hit..." << std::endl;
1325 dmLog << "m_zMinBarrel: " << m_zMinBarrel << std::endl;
1326 dmLog << "m_zMaxBarrel: " << m_zMaxBarrel << std::endl;
1327 dmLog << "m_rMinAccordion: " << m_rMinAccordion << std::endl;
1328 dmLog << "m_rMaxAccordion: " << m_rMaxAccordion << std::endl;
1329 dmLog << "r,z,eta,phi " << radiusZpos << " " << zZpos << " " << etaZpos << " " << phiZpos << std::endl;
1330 dmLog << "x,y,z (Atlas) " << p.x() << " " << p.y() << " " << p.z() << std::endl;
1331 dmLog << " inSTAC " << inSTAC << std::endl;
1332 dmLog << " eDeposited " << thisStepEnergyDeposit << std::endl;
1333 const G4VPhysicalVolume* vol = thisStepPoint->GetPhysicalVolume();
1334 const G4String volName = vol->GetName();
1335 dmLog << " volName " << volName << std::endl;
1336 const G4int ndep = g4navigation->GetDepth();
1337 for (G4int ii=0;ii<=ndep;ii++) {
1338 const G4VPhysicalVolume* v1 = g4navigation->GetVolume(ii);
1339 const G4String vname = v1->GetName();
1340 dmLog << "vname " << vname << std::endl;
1341 }
1342 if (thisStepEnergyDeposit > 1.*CLHEP::MeV) {
1343 ATH_MSG_ERROR(dmLog.str());
1344 } else {
1345 ATH_MSG_WARNING(dmLog.str());
1346 }
1347 }
1348 // in test beam case, we can get there for leakage on the side (in phi) of the module
1349 // in this case, we attribute an identifier like inactive material
1350 else
1351 {
1352 G4bool MapDetail=false;
1353 this->findCell( currentCellData, xZpos, yZpos, zZpos, radiusZpos, etaZpos, phiZpos, MapDetail );
1354 // ATH_MSG_ERROR(" Lateral lakage r,eta,phi " << radiusZpos << " " << etaZpos << " "
1355 // << phiZpos << " sampling/region/eta/phi " << currentCellData.sampling << " " <<
1356 // currentCellData.region << " " << currentCellData.etaBin << " " << currentCellData.phiBin);
1357 // protect against small space between z=4m and real beginning of ieta=1 in strips
1358 if (currentCellData.sampling==1 && currentCellData.region==0 && currentCellData.etaBin==0) {
1359 currentCellData.etaBin=1;
1360 // ATH_MSG_ERROR("S1R0 etabin 0 found r,z,phi local " << radiusZpos << " "
1361 // << " " << zZpos << " " << phiZpos);
1362 }
1363 result << 4 // LArCalorimeter
1364 << 1 // LArEM
1365 << currentCellData.zSide
1366 << currentCellData.sampling
1367 << currentCellData.region
1368 << currentCellData.etaBin
1369 << currentCellData.phiBin;
1370 return result;
1371 }
1372 }
1373
1374 result << 10 // LArCalorimeter
1375 << currentCellData.zSide * 4 // LArBarrel
1376 << type
1377 << sampling
1378 << region
1379 << currentCellData.etaBin
1380 << currentCellData.phiBin;
1381
1382#ifdef DEBUGHITS
1383 ATH_MSG_VERBOSE("Here the identifier for the barrel DEAD materials ---->");
1384 ATH_MSG_VERBOSE("Type ----> " << type);
1385 ATH_MSG_VERBOSE("Sampling ----> " << sampling);
1386 ATH_MSG_VERBOSE("Region ----> " << region);
1387 ATH_MSG_VERBOSE("zSide ----> " << currentCellData.zSide*4);
1388 ATH_MSG_VERBOSE("etaBin ----> " << currentCellData.etaBin);
1389 ATH_MSG_VERBOSE("phiBin ----> " << currentCellData.phiBin);
1390#endif
1391
1392 // if (!Geometry::CheckDMIdentifier(type,sampling,region,currentCellData.etaBin,currentCellData.phiBin)) {
1393 // ATH_MSG_ERROR(" ** Bad DM identifier " << type << " " << sampling << " " << region << " "
1394 // << currentCellData.etaBin << " " << currentCellData.phiBin);
1395 // ATH_MSG_ERROR("x,y,z,r,eta,phi" << xZpos << " " << yZpos << " " << zZpos <<
1396 // " " << radiusZpos << " " << etaZpos << " " << phiZpos);
1397 // }
1398
1399 }
1400
1401 return result;
1402
1403 }
#define M_PI
#define ATH_MSG_ERROR(x)
#define ATH_MSG_VERBOSE(x)
#define ATH_MSG_WARNING(x)
Gaudi::Property< bool > m_testbeam
virtual void findCell(CalcData &currentCellData, const double &x, const double &y, const double &z, const double &r, const double &eta, const double &phi, const bool detail) const override final

◆ CalculateIdentifier()

LArG4Identifier LArG4::Barrel::Geometry::CalculateIdentifier ( const G4Step * a_step) const
finaloverridevirtual

Definition at line 923 of file LArBarrelGeometry.cxx.

924 {
925
926 // The default result is a blank identifier.
927 LArG4Identifier result = LArG4Identifier();
928
929 // Get all the required information from the current step
930 const G4NavigationHistory* g4navigation = a_step->GetPreStepPoint()->GetTouchable()->GetHistory();
931 const G4int ndep = g4navigation->GetDepth();
932 bool inSTAC = false;
933 int zside=1;
934 G4int indECAM = -999;
935
936 // Now navigate through the volumes hierarchy
937 for (G4int ii=0;ii<=ndep;ii++) {
938 const G4String& vname = g4navigation->GetVolume(ii)->GetName();
939 // FIXME Need to find a way to avoid these string-comparisons
940 if ( indECAM<0 && vname == m_ecamName ) indECAM=ii;
941 if ( !inSTAC && vname.find("STAC") !=std::string::npos) inSTAC=true;
942 if ( vname.find("NegPhysical") != std::string::npos) zside=-1;
943 }
944 if (indECAM>=0)
945 result = this->CalculateECAMIdentifier( a_step , indECAM, inSTAC, zside) ;
946 else
947 ATH_MSG_ERROR("LArBarrel::Geometry::CalculateIdentifier ECAM volume not found in hierarchy");
948
949 return result;
950 }
LArG4Identifier CalculateECAMIdentifier(const G4Step *, const G4int indEcam, const bool inSTAC=true, int zside=1) const

◆ CalculateSuperResolutionIdentifier()

LArG4Identifier LArG4::Barrel::Geometry::CalculateSuperResolutionIdentifier ( const G4Step * a_step) const
finaloverridevirtual

Definition at line 952 of file LArBarrelGeometry.cxx.

952 {
953
954 LArG4Identifier result;
955
956 // Convert step to local half-barrel coordinates
957 const G4StepPoint *preStep = a_step->GetPreStepPoint();
958 const G4StepPoint *postStep = a_step->GetPostStepPoint();
959 G4ThreeVector midGlobal = (preStep->GetPosition() + postStep->GetPosition()) * 0.5;
960
961 const G4NavigationHistory* g4nav = preStep->GetTouchable()->GetHistory();
962
963 G4int indECAM = -1;
964 const G4int ndep = g4nav->GetDepth();
965
966 for (G4int ii = 0; ii <= ndep; ++ii) {
967 if (g4nav->GetVolume(ii)->GetName() == m_ecamName) {
968 indECAM = ii;
969 break;
970 }
971 }
972
973 if (indECAM == -1) {
974 ATH_MSG_ERROR("ECAM volume not found!");
975 return result;
976 }
977
978 G4AffineTransform localTrans = g4nav->GetTransform(indECAM);
979 G4ThreeVector midLocal = localTrans.TransformPoint(midGlobal);
980
981 double x = midLocal.x(), y = midLocal.y(), z = midLocal.z();
982 double radius = sqrt(x*x + y*y);
983 double eta = midLocal.pseudoRapidity();
984 double phi = midLocal.phi() < 0 ? midLocal.phi() + 2*M_PI : midLocal.phi();
985
986 // --- Find the standard (low-resolution) cell ---
987 CalcData currentCellData;
988 bool MapDetail = false;
989 findCell(currentCellData, x, y, z, radius, eta, phi, MapDetail);
990
991 int etaBin = currentCellData.etaBin;
992 int phiBin = currentCellData.phiBin;
993 int sampling = currentCellData.sampling;
994 int region = currentCellData.region;
995
996 // --- Compute minimum boundaries (low resolution cell) ---
997 double eta_min = 0.0;
998 double delta_eta = 0.025; // default
999 if (sampling == 1) {
1000 if (region == 0) {
1001 delta_eta = 0.003125;
1002 eta_min = etaBin * delta_eta;
1003 } else {
1004 delta_eta = 0.025;
1005 eta_min = 1.4 + etaBin * delta_eta;
1006 }
1007 } else if (sampling == 2) {
1008 delta_eta = 0.025;
1009 eta_min = etaBin * delta_eta;
1010 } else if (sampling == 3) {
1011 delta_eta = 0.05;
1012 eta_min = etaBin * delta_eta;
1013 }
1014
1015 double phiGranularity = (sampling == 1 && region == 0) ? (2*M_PI/64) : (2*M_PI/256);
1016 double phi_min = phiBin * phiGranularity;
1017
1018 double radius_min = 0.0;
1019 double radius_max = 0.0;
1020
1021 if (sampling == 1) {
1022 radius_min = (etaBin == 0) ? m_rMinAccordion : Rmax1[etaBin - 1];
1023 radius_max = Rmax1[etaBin];
1024 } else if (sampling == 2) {
1025 radius_min = Rmax1[etaBin];
1026 radius_max = Rmax2[etaBin];
1027 } else if (sampling == 3) {
1028 radius_min = Rmax2[etaBin];
1029 radius_max = m_rMaxAccordion;
1030 }
1031
1032 // --- Super-Resolution granularity (bins within cell) ---
1033 const int nBinsEtaSR = 16; // 3 bits
1034 const int nBinsPhiSR = 16; // 2 bits
1035 const int nBinsRadiusSR = 5; // 3 bits
1036
1037 double delta_eta_SR = delta_eta / nBinsEtaSR;
1038 double delta_phi_SR = phiGranularity / nBinsPhiSR;
1039 double delta_r_SR = (radius_max - radius_min) / nBinsRadiusSR;
1040
1041 auto computeBin = [](double val, double min_lr, double delta_sr, int maxBins) {
1042 const int overflow = maxBins - 1;
1043 if (delta_sr == 0.0) return overflow;
1044 int bin = int((val - min_lr) / delta_sr);
1045 return std::clamp(bin, 0, overflow);;
1046 };
1047
1048 // --- Compute SR bin numbers (local within the cell) ---
1049 int etaBinSR = computeBin(eta, eta_min, delta_eta_SR, nBinsEtaSR);
1050 int phiBinSR = computeBin(phi, phi_min, delta_phi_SR, nBinsPhiSR);
1051 int rBinSR = computeBin(radius, radius_min, delta_r_SR, nBinsRadiusSR);
1052
1053 // Packing (16 bits):
1054 // Bit 15: Marker bit (always 1)
1055 // Bits 9-14 (6 bits): etaBinSR (shifted by 1)
1056 // Bits 3-8 (6 bits): phiBinSR (shifted by 1)
1057 // Bits 0-2 (3 bits): rBinSR (reduced by 1 bit)
1058
1059 unsigned short int packedID = ( 0x8000 ) // Marker bit always set (bit 15)
1060 | ( (etaBinSR & 0x3F) << 9 ) // 6 bits for eta (bits 9-14)
1061 | ( (phiBinSR & 0x3F) << 3 ) // 6 bits for phi (bits 3-8)
1062 | ( rBinSR & 0x07 ); // 3 bits for r (bits 0-2)
1063
1064 // --- Store packedID into LArG4Identifier ---
1065 result << packedID;
1066
1067 return result;
1068 }
Scalar eta() const
pseudorapidity method
Scalar phi() const
phi method
#define y
#define x
#define z
static const G4double Rmax2[56]
maximum values of active S2 region from electrode drawing
static const G4double Rmax1[448]
maximum values of active S1 region from electrode drawing (first eta value checked by D....
setSAddress setEtaMS setDirPhiMS setDirZMS setBarrelRadius setEndcapAlpha setEndcapRadius setPhiMap phiBin
setSAddress setEtaMS setDirPhiMS setDirZMS setBarrelRadius setEndcapAlpha setEndcapRadius setInterceptInner setEtaMap etaBin

◆ CheckDMIdentifier()

bool LArG4::Barrel::Geometry::CheckDMIdentifier ( int type,
int sampling,
int region,
int eta,
int phi ) const
private

Definition at line 1405 of file LArBarrelGeometry.cxx.

1406 {
1407
1408 if (type <1 || type > 2) return false;
1409 if (type==1) {
1410 if (sampling<1 || sampling>2) return false;
1411 if (sampling==1) {
1412 if (region!=3 && region !=4) return false;
1413 if (phi<0 || phi>63) return false;
1414 if (region==3) {
1415 if (eta<0 || eta>14) return false;
1416 }
1417 if (region==4) {
1418 if (eta !=0) return false;
1419 }
1420 }
1421 if (sampling==2) {
1422 if (region !=0 && region !=2) return false;
1423 if (phi<0 || phi>63) return false;
1424 if (region==0){
1425 if (eta<0 || eta>9) return false;
1426 }
1427 if (region==2) {
1428 if (eta<0 || eta>4) return false;
1429 }
1430 }
1431 }
1432 if (type==2) {
1433 if (sampling<1 || sampling >3) return false;
1434 if (region !=0) return false;
1435 if (eta!=0) return false;
1436 if (phi<0 || phi>63) return false;
1437 }
1438 return true;
1439 }

◆ CheckLArIdentifier()

bool LArG4::Barrel::Geometry::CheckLArIdentifier ( int sampling,
int region,
int eta,
int phi ) const
private

Definition at line 1442 of file LArBarrelGeometry.cxx.

1443 {
1444 if (sampling<0 || sampling >3) return false;
1445 if (sampling==0) {
1446 if (region!=0) return false;
1447 if (eta<0 || eta>60) return false;
1448 if (phi<0 || phi>63) return false;
1449 }
1450 if (sampling==1) {
1451 if (region<0 || region >1) return false;
1452 if (region==0) {
1453 if (eta<1 || eta>447) return false;
1454 if (phi<0 || phi>63) return false;
1455 }
1456 if (region==1) {
1457 if (eta<0 || eta>2) return false;
1458 if (phi<0 || phi>255) return false;
1459 }
1460 }
1461 if (sampling==2) {
1462 if (region<0 || region >1) return false;
1463 if (region==0) {
1464 if (eta<0 || eta>55) return false;
1465 if (phi<0 || phi>255) return false;
1466 }
1467 if (region==1) {
1468 if (eta!=0) return false;
1469 if (phi<0 || phi>255) return false;
1470 }
1471 }
1472 if (sampling==3) {
1473 if (region !=0) return false;
1474 if (eta<0 || eta>26) return false;
1475 if (phi<0 || phi>255) return false;
1476 }
1477 return true;
1478 }

◆ Distance_Abs()

double LArG4::Barrel::Geometry::Distance_Abs ( const double & x,
const double & y,
const int & nabs,
const int & Num_Straight,
const int & Num_Coude ) const
private

Definition at line 227 of file LArBarrelGeometry.cxx.

230 {
231 //
232 // FrameWork is consistent with the one used to PhiCell determination
233 // e.g. it assumes HERE to be the LOCAL one of "stac_phys1",
234 // (mother of ACCordion volumes) from which Z> 0. and Z < 0. half_barrel
235 // parts are then defined.
236 //
237 // One needs POINTERS to Electrode neutral fibers
238 // either for straight parts or for folds
239 //
240 // u unit 2D_Vector along straight part of the electrode neutral fiber
241 const G4double u[2] = { m_absorber->Cosu(Num_Straight, PhiCell), m_absorber->Sinu(Num_Straight, PhiCell) };
242 // Middle m_coordinates of this straight part of the electrode neutral fiber
243 const G4double Xm[2] = { m_absorber->XCentAbs(Num_Straight, PhiCell), m_absorber->YCentAbs(Num_Straight, PhiCell) };
244 // m_Hit Vector components
245 double dx = xhit - Xm[0]; double dy = yhit - Xm[1];
246
247 // First compute algebric distance hit (2D) the 2D_projection of the
248 // m_Hit Vector on this electrode neutral fiber.
249 const double hit = dx*u[0] + dy*u[1];
250
251 //
252 // Flat of Fold Region ?
253 //
254 if(std::fabs(hit) < m_absorber->HalfLength(Num_Straight,PhiCell)) {
255 // Flat Region
256 return u[0]*dy - u[1]*dx;
257 }
258 else {
259 // Fold Center c_coordinates
260 const G4double Xc[2] = { m_coudeabs->XCentCoude(Num_Coude, PhiCell), m_coudeabs->YCentCoude(Num_Coude, PhiCell) };
261 // c_Hit Vector components and its length
262 dx = xhit - Xc[0];
263 dy = yhit - Xc[1];
264 const double dr = sqrt( dx*dx + dy*dy);
265 return (Num_Coude%2 == m_parity) ? (m_rint_eleFib-dr) : (dr - m_rint_eleFib);
266 } // end of Fold Regions
267 } // end of the function Distance_Abs
const LArCoudeAbsorbers * m_coudeabs
const LArStraightAbsorbers * m_absorber
@ u
Enums for curvilinear frames.
Definition ParamDefs.h:77

◆ Distance_Ele()

double LArG4::Barrel::Geometry::Distance_Ele ( const double & x,
const double & y,
const int & PhiC,
int & Num_Straight,
const int & Num_Coude,
double & xl ) const
private

Definition at line 153 of file LArBarrelGeometry.cxx.

156 {
157 //
158 // FrameWork is consistent with the one used to PhiCell determination
159 // e.g. it assumes HERE to be the LOCAL one of "stac_phys1",
160 // (mother of ACCordion volumes) from which Z> 0. and Z < 0. half_barrel
161 // parts are then defined.
162 //
163 // One needs POINTERS to Electrode neutral fibers
164 // either for straight parts or for folds
165 //
166 // Fold Center ccoordinates
167 const G4double Xc[2] = { m_coudeelec->XCentCoude(Num_Coude, PhiCell), m_coudeelec->YCentCoude(Num_Coude, PhiCell) };
168 const G4double radfold = sqrt(Xc[0]*Xc[0]+Xc[1]*Xc[1]);
169 const G4double radhit = sqrt(xhit*xhit+yhit*yhit);
170
171 // check if the assumed straight number is the correct one
172 // (this can be wrong when we are close to a fold and there is sagging)
173 if (Num_Coude == Num_Straight && radhit <radfold) {
174 if (Num_Straight>0) { Num_Straight = Num_Straight-1; }
175 // ATH_MSG_VERBOSE("radhit,radfold " << radhit << " " << radfold << " change straight by +1");
176 }
177 if (Num_Coude == (Num_Straight+1) && radhit > radfold) {
178 if (Num_Straight<12) { Num_Straight = Num_Straight+1; }
179 // ATH_MSG_VERBOSE("radhit,radfold " << radhit << " " << radfold << " change straight by -1");
180 }
181
182 // u unit 2D_Vector along straight part of the electrode neutral fiber
183 const double u[2] = { m_electrode->Cosu(Num_Straight, PhiCell), m_electrode->Sinu(Num_Straight, PhiCell) };
184 // Middle m_coordinates of this straight part of the electrode neutral fiber
185 const double Xm[2] = { m_electrode->XCentEle(Num_Straight, PhiCell), m_electrode->YCentEle(Num_Straight, PhiCell) };
186 // m_Hit Vector components
187 double dx = xhit - Xm[0];
188 double dy = yhit - Xm[1];
189
190 // First compute algebric distance m_hit (2D) the 2D_projection of the
191 // m_Hit Vector on this electrode neutral fiber.
192 const double hit = dx*u[0] + dy*u[1];
193
194 //
195 // Flat of Fold Region ?
196 //
197 const G4double Half_Elec(m_electrode->HalfLength(Num_Straight,PhiCell));
198
199 if(std::fabs(hit) < Half_Elec) {
200 // Flat Region
201 xl=hit/Half_Elec;
202 return u[0]*dy - u[1]*dx;
203 }
204 else {
205 // Fold region
206 // c_Hit Vector components and its length
207 dx = xhit - Xc[0];
208 dy = yhit - Xc[1];
209 const double dr = sqrt( dx*dx + dy*dy);
210 if (Num_Coude==Num_Straight) { xl=-1.; }
211 else xl=+1;
212 return (Num_Coude%2 == m_parity) ? (m_rint_eleFib-dr) : (dr - m_rint_eleFib);
213 } // end of Fold Regions
214 } // end of the function Distance_Ele
const LArCoudeElectrodes * m_coudeelec
const LArStraightElectrodes * m_electrode

◆ finalize()

StatusCode LArG4::Barrel::Geometry::finalize ( )
finaloverridevirtual

Definition at line 124 of file LArBarrelGeometry.cxx.

125 {
126 if (m_rc) delete [] m_rc;
127 if (m_phic) delete [] m_phic;
128 if (m_delta) delete [] m_delta;
129 if (m_xc) delete [] m_xc;
130 if (m_yc) delete [] m_yc;
131
132 return StatusCode::SUCCESS;
133 }

◆ findCell()

void LArG4::Barrel::Geometry::findCell ( CalcData & currentCellData,
const double & x,
const double & y,
const double & z,
const double & r,
const double & eta,
const double & phi,
const bool detail ) const
finaloverridevirtual

Definition at line 580 of file LArBarrelGeometry.cxx.

588 {
589
590 currentCellData.cellID = 0;
591
592 if (aRadius < m_rc[0] || aRadius >= m_rc[m_Nbrt1-1]) {
593#ifdef DEBUGHITS
594 ATH_MSG_VERBOSE(" Outside Accordion " << aRadius << " " << m_rc[0] << " " << m_rc[m_Nbrt1-1]);
595#endif
596 return; // outside accordion
597 }
598
599 // set the straight section number
600 currentCellData.nstraight=0;
601 for (int i=1;i<m_Nbrt1;i++) {
602 if (m_rc[i] > aRadius) { break; }
603 currentCellData.nstraight++;
604 }
605 if (currentCellData.nstraight <0 || currentCellData.nstraight >= m_Nbrt) {
606 ATH_MSG_ERROR("Invalid straight number " << currentCellData.nstraight << " " << aRadius);
607 return;
608 }
609
610 // get the closest fold number
611 currentCellData.nfold=currentCellData.nstraight;
612 if (std::fabs(aRadius-m_rc[currentCellData.nfold]) > std::fabs(aRadius-m_rc[currentCellData.nfold+1]) ) {
613 currentCellData.nfold +=1;
614 }
615 if (currentCellData.nfold <0 || currentCellData.nfold >= m_Nbrt1) {
616 ATH_MSG_ERROR("Invalid fold number " << currentCellData.nfold);
617 return;
618 }
619
620
621#ifdef DEBUGHITS
622 ATH_MSG_VERBOSE(" BarrelGeometry: radius,eta,phi " << aRadius << " " << anEta << " ");
623 ATH_MSG_VERBOSE(" Straight/Fold numbers " << currentCellData.nstraight << " " << currentCellData.nfold);
624#endif
625
626 // eta and longitudinal segmentations
627 G4int ireg,isamp,ieta,isamp2,ieta2;
628 currentCellData.cellID = this->SampSeg(anEta,aRadius,zPosition,ireg,isamp,ieta,isamp2,ieta2);
629
630 currentCellData.etaBin = ieta;
631 currentCellData.sampling = isamp;
632 currentCellData.region = ireg;
633 currentCellData.etaMap = ieta2;
634 currentCellData.sampMap = isamp2;
635
636 // compute electrode number in phi
637 int phicell = this->PhiGap(aRadius,xPosition,yPosition);
638 if (phicell<0) phicell=0;
639 // for test beam, some protection
640 if (m_NCellTot !=1024) {
641 if (phicell>=m_NCellTot) {
642 if (phicell<512) { phicell=m_NCellTot-1; }
643 else { phicell=0; }
644 currentCellData.cellID=0;
645 }
646 }
647
648#ifdef DEBUGHITS
649 ATH_MSG_VERBOSE(" phigap " << phicell);
650#endif
651
652 // compute readout cell number
653 int sampling_phi_nGaps=4;
654 if (currentCellData.region==0 && currentCellData.sampling==1) { sampling_phi_nGaps=16; }
655
656 if (currentCellData.cellID==0) {
657 currentCellData.phiBin = (G4int) ( phicell/sampling_phi_nGaps );
658 currentCellData.distElec=9999.;
659 return;
660 }
661
662 // compute distance to electrode
663 G4double xl;
664 G4int nstr = currentCellData.nstraight;
665 const G4double distElec = this->Distance_Ele(xPosition,yPosition,phicell,nstr,currentCellData.nfold,xl);
666
667#ifdef DEBUGHITS
668 ATH_MSG_VERBOSE(" distElec " << distElec);
669#endif
670
671 // if distance is < 2.5mm we are sure to be in the correct gap
672 if (std::fabs(distElec) > 2.5) {
673 // try +-2 electrode in phi to get minimum distance
674 double dElecMin=distElec;
675 double xlmin=xl;
676 int phicellmin=phicell;
677 for (int ii=-2;ii<3;ii++) {
678 if (ii==0) { continue; }
679 int phicellnew = phicell+ii;
680 // for test beam no phi wrapping
681 if (m_NCellTot != 1024 && ( phicellnew<0 || phicellnew >= m_NCellTot)) { continue; }
682 if (phicellnew < 0) { phicellnew += m_NCellTot; }
683 if (phicellnew >= m_NCellTot) { phicellnew -= m_NCellTot; }
684 double xln;
685 int nstr2=currentCellData.nstraight;
686 double dElec = Distance_Ele(xPosition,yPosition,phicellnew,nstr2,currentCellData.nfold,xln);
687 if (std::fabs(dElec)<std::fabs(dElecMin)) {
688 phicellmin=phicellnew;
689 xlmin=xln;
690 dElecMin = dElec;
691 nstr=nstr2;
692 }
693 }
694 currentCellData.phiGap = phicellmin;
695 currentCellData.distElec = dElecMin;
696 currentCellData.xl = xlmin;
697 currentCellData.nstraight = nstr;
698 } // end distance >2.5mm
699 else {
700 currentCellData.phiGap=phicell;
701 currentCellData.distElec=distElec;
702 currentCellData.xl=xl;
703 currentCellData.nstraight=nstr;
704 }
705
706#ifdef DEBUGHITS
707 ATH_MSG_VERBOSE(" final phiGap,distElec,xl " << currentCellData.phiGap << " " << currentCellData.distElec << " "
708 << currentCellData.xl);
709#endif
710
711 // compute distance to absorber
712
713 G4int nabs;
714 if (currentCellData.distElec<0) nabs=currentCellData.phiGap;
715 else nabs=currentCellData.phiGap+1;
716 if (nabs >= m_NCellMax) nabs -= m_NCellMax;
717 currentCellData.distAbs = Distance_Abs(xPosition,yPosition,nabs,currentCellData.nstraight,currentCellData.nfold);
718#ifdef DEBUGHITS
719 ATH_MSG_VERBOSE(" nabs,distAbs " << nabs << " " << currentCellData.distAbs);
720#endif
721
722 // in some rare cases near fold, the closest distance could give the wrong gap
723 // in such case, the signs of distAbs and distElec are not opposite as they should
724 if ((currentCellData.distAbs>0. && currentCellData.distElec>0) ||
725 (currentCellData.distAbs<0. && currentCellData.distElec<0) ) {
726 // ATH_MSG_VERBOSE("distElec and distAbs same sign " << currentCellData.distElec << " " << currentCellData.distAbs);
727 // ATH_MSG_VERBOSE(" currentCellData.phiGap " << currentCellData.phiGap);
728 if (std::fabs(currentCellData.distElec)>std::fabs(currentCellData.distAbs)) {
729 if (currentCellData.distAbs>0) { currentCellData.phiGap += 1; }
730 if (currentCellData.distAbs<0) { currentCellData.phiGap -= 1; }
731 if (m_NCellTot != 1024) {
732 if (currentCellData.phiGap <0) { currentCellData.phiGap=0; }
733 if (currentCellData.phiGap >= m_NCellTot) { currentCellData.phiGap = m_NCellTot-1; }
734 }
735 else {
736 if (currentCellData.phiGap < 0) { currentCellData.phiGap += m_NCellTot; }
737 if (currentCellData.phiGap >= m_NCellTot) { currentCellData.phiGap -= m_NCellTot; }
738 }
739 currentCellData.distElec = Distance_Ele(xPosition,yPosition,currentCellData.phiGap,currentCellData.nstraight,currentCellData.nfold,currentCellData.xl);
740 // ATH_MSG_VERBOSE(" new phiGap,distElec " << currentCellData.phiGap << " " << currentCellData.distElec);
741 }
742 }
743
744 currentCellData.phiBin = (G4int) ( currentCellData.phiGap/sampling_phi_nGaps );
745
746 if (MapDetail) {
747 // compute x0,y0 coordinates in local electrode frame, using closest fold
748 // as reference
749 const G4double alpha = m_coudeelec->PhiRot(currentCellData.nfold,currentCellData.phiGap);
750 const G4double dx=xPosition-m_coudeelec->XCentCoude(currentCellData.nfold,currentCellData.phiGap);
751 const G4double dy=yPosition-m_coudeelec->YCentCoude(currentCellData.nfold,currentCellData.phiGap);
752 const G4double dx1=dx*cos(alpha)-dy*sin(alpha);
753 const G4double dy1=dx*sin(alpha)+dy*cos(alpha);
754 currentCellData.x0 = dx1 + m_xc[currentCellData.nfold];
755 currentCellData.y0 = dy1 + m_yc[currentCellData.nfold];
756 if (m_parity==1) { currentCellData.y0 = -1*currentCellData.y0; }
757 }
758
759
760 } // end of findCell method
double Distance_Ele(const double &x, const double &y, const int &PhiC, int &Num_Straight, const int &Num_Coude, double &xl) const
G4int PhiGap(const double &, const double &, const double &) const
G4int SampSeg(G4double, G4double, G4double, G4int &, G4int &, G4int &, G4int &, G4int &) const
double Distance_Abs(const double &x, const double &y, const int &nabs, const int &Num_Straight, const int &Num_Coude) const

◆ GetRphi()

void LArG4::Barrel::Geometry::GetRphi ( )
private

phi vs r of first absorber in nominal geometry

Initialize r-phi reference map (called from constructor)

Definition at line 764 of file LArBarrelGeometry.cxx.

765 {
766 const G4double dl=0.001;
767 const G4double inv_dl = 1. / dl;
768 G4double cenx[15],ceny[15];
769 //G4double xl,xl2;
770 std::vector<G4double> sum1(5000);
771 std::vector<G4double> sumx(5000);
772 //xl=0;
773 //xl2=0.;
774 m_NRphi=5000;
775 m_Rmin=1500.;
776 m_dR=0.10;
777 m_Rmax=0.;
778
779 const G4double rint= m_rint_eleFib;
780 const G4double inv_rint = 1. / rint;
781 const G4double dt=dl * inv_rint;
782 const G4double inv_dt = 1. / dt;
783
784 for (G4int i=0;i<m_NRphi;i++) {
785 sum1[i]=0.;
786 sumx[i]=0.;
787 }
788 for (G4int i=0;i<15;i++) {
789 cenx[i]=m_rc[i]*cos(m_phic[i]);
790 ceny[i]=m_rc[i]*sin(m_phic[i]);
791 }
792
793 for (G4int i=0; i<15; i++) {
794
795 // fold
796 G4double phi0,phi1;
797 if (i==0) {
798 // first fold goes up
799 if (m_parity==0) {
800 phi0=-CLHEP::pi/2.;
801 phi1=-m_delta[0];
802 }
803 // first fold goes down
804 else {
805 phi0=m_delta[0];
806 phi1=CLHEP::pi/2;
807 }
808 }
809 else if (i==14) {
810 if (m_parity==0) {
811 phi0=-CLHEP::pi+m_delta[13];
812 phi1=-CLHEP::pi/2.;
813 }
814 else {
815 phi0=CLHEP::pi/2;
816 phi1=CLHEP::pi - m_delta[13];
817 }
818 }
819 else {
820 if (i%2==(1-m_parity)) {
821 phi0=m_delta[i];
822 phi1=CLHEP::pi-m_delta[i-1];
823 }
824 else {
825 phi0=-CLHEP::pi+m_delta[i-1];
826 phi1=-m_delta[i];
827 }
828 }
829 //xl2+=rint*std::fabs(phi1-phi0);
830 const G4int nstep=int((phi1-phi0)*inv_dt)+1;
831 for (int ii=0;ii<nstep;ii++) {
832 //xl+=dl;
833 const G4double phi=phi0+dt*((G4double)ii);
834 const G4double x=cenx[i]+rint*cos(phi);
835 const G4double y=ceny[i]+rint*sin(phi);
836 const G4double radius=sqrt(x*x+y*y);
837 if (radius>m_Rmax) { m_Rmax=radius; }
838 const G4double phid=atan(y/x);
839 const G4int ir=((int) ((radius-m_Rmin)/m_dR) );
840 if (ir>=0 && ir < m_NRphi) {
841 sum1[ir]+=1.;
842 sumx[ir]+=phid;
843 }
844 }
845
846 // straight section
847 if (i<14) {
848 const G4double dx=cenx[i+1]-cenx[i];
849 const G4double dy=ceny[i+1]-ceny[i];
850 const G4double along=std::sqrt(dx*dx+dy*dy-4.*rint*rint);
851 const G4double x0=0.5*(cenx[i+1]+cenx[i]);
852 const G4double y0=0.5*(ceny[i+1]+ceny[i]);
853 const G4double phi = (i%2==m_parity) ? CLHEP::pi/2-m_delta[i] : -CLHEP::pi/2.+m_delta[i];
854 const G4double x1=x0-0.5*along*cos(phi);
855 const G4double y1=y0-0.5*along*sin(phi);
856 //xl2+=along;
857 const int nstep=int(along*inv_dl)+1;
858 for (int ii=0;ii<nstep;ii++) {
859 //xl+=dl;
860 const G4double x=x1+dl*((G4double)ii)*cos(phi);
861 const G4double y=y1+dl*((G4double)ii)*sin(phi);
862 const G4double radius=sqrt(x*x+y*y);
863 if (radius>m_Rmax) { m_Rmax=radius; }
864 const G4double phid=atan(y/x);
865 const G4int ir=((int) ((radius-m_Rmin)/m_dR) );
866 if (ir>=0 && ir < m_NRphi) {
867 sum1[ir]+=1.;
868 sumx[ir]+=phid;
869 }
870 }
871 }
872 }
873 // ATH_MSG_VERBOSE("total electrode length " << xl << " " << xl2);
874 // ATH_MSG_VERBOSE("rmax in accordion " << m_Rmax);
875 for (int i=0; i<m_NRphi; i++) {
876 if (sum1[i]>0) {
877 m_Rphi[i]=sumx[i]/sum1[i];
878 // Not used:
879 //G4double radius = m_Rmin + ((G4double(i))+0.5)*m_dR;
880 //ATH_MSG_VERBOSE(" GUTEST r,phi0 " << radius << " " << m_Rphi[i]);
881 }
882 else { m_Rphi[i]=0.; }
883 }
884 }
#define pi
int ir
counter of the current depth
Definition fastadd.cxx:49

◆ initialize()

StatusCode LArG4::Barrel::Geometry::initialize ( )
finaloverridevirtual

Definition at line 50 of file LArBarrelGeometry.cxx.

51 {
52 // initialize the geometry.
53 // Access source of detector parameters.
54
56
57 // number of straight sections (should be 14)
58 m_Nbrt = (int) (parameters->GetValue("LArEMBnoOFAccZigs"));
59 // Number of ZIGs + 1 i.e. 15 = number of folds
60 m_Nbrt1 = m_Nbrt + 1;
61 // phi of first absorber
62 m_gam0 = parameters->GetValue("LArEMBAbsPhiFirst");
63 // radius of curvature of neutral fiber in the folds
64 m_rint_eleFib = parameters->GetValue("LArEMBNeutFiberRadius");
65
66 m_rc = new double[m_Nbrt1];
67 m_phic = new double[m_Nbrt1];
68 m_delta = new double[m_Nbrt1];
69 m_xc = new double[m_Nbrt1];
70 m_yc = new double[m_Nbrt1];
71
72 // r,phi positions of the centre of the folds (nominal geometry)
73 for (G4int idat = 0; idat < m_Nbrt1 ; idat++) {
74 m_rc[idat] = (double) parameters->GetValue("LArEMBRadiusAtCurvature",idat);
75 m_phic[idat] = (double) parameters->GetValue("LArEMBPhiAtCurvature",idat);
76 m_delta[idat] = (double) parameters->GetValue("LArEMBDeltaZigAngle",idat);
77 m_xc[idat] = m_rc[idat]*cos(m_phic[idat]);
78 m_yc[idat] = m_rc[idat]*sin(m_phic[idat]);
79 }
80 // define parity of accordion waves: =0 if first wave goes up, 1 if first wave goes down in the local frame
81 m_parity=0;
82 if (m_phic[0]<0.) { m_parity=1; }
83 //
84 m_rMinAccordion = parameters->GetValue("LArEMBRadiusInnerAccordion");
85 m_rMaxAccordion = parameters->GetValue("LArEMBFiducialRmax");
86 m_etaMaxBarrel = parameters->GetValue("LArEMBMaxEtaAcceptance");
87 m_zMinBarrel = parameters->GetValue("LArEMBfiducialMothZmin");
88 m_zMaxBarrel = parameters->GetValue("LArEMBfiducialMothZmax");
89 m_zMaxBarrelDMMargin = 10.0; // 10 mm margin
90 // === GU 11/06/2003 total number of cells in phi
91 // to distinguish 1 module (testbeam case) from full Atlas
92 m_NCellTot = (int) (parameters->GetValue("LArEMBnoOFPhysPhiCell"));
93 // total number of cells in phi to distinguish 1 module (testbeam case) from full Atlas
94 m_testbeam=false;
95 if (m_NCellTot != 1024) {
96 m_testbeam=true;
97 }
98 m_NCellMax=1024;
99 // ===
100
101 // Initialize r-phi reference map
102 this->GetRphi();
103
104 if (m_detectorName.empty()) m_ecamName = "LAr::EMB::ECAM";
105 else m_ecamName = m_detectorName+"::LAr::EMB::ECAM";
106
107
108 return StatusCode::SUCCESS;
109 }
LArGeo::VDetectorParameters LArVG4DetectorParameters
Gaudi::Property< std::string > m_detectorName
void GetRphi()
phi vs r of first absorber in nominal geometry
static const VDetectorParameters * GetInstance()

◆ initializeForSDCreation()

void LArG4::Barrel::Geometry::initializeForSDCreation ( )
finaloverridevirtual

Definition at line 113 of file LArBarrelGeometry.cxx.

114 {
115 // get pointers to access G4 geometry
120 }
static const LArCoudeAbsorbers * GetInstance(const std::string &strDetector="")
static const LArCoudeElectrodes * GetInstance(const std::string &strDetector="")
static const LArStraightAbsorbers * GetInstance(const std::string &strDetector="")
static const LArStraightElectrodes * GetInstance(const std::string &strDetector="")

◆ Phi0()

G4double LArG4::Barrel::Geometry::Phi0 ( G4double radius) const
private

Definition at line 889 of file LArBarrelGeometry.cxx.

890 {
891 // TODO This function could be simplified.
892 G4int ir;
893 if (radius < m_Rmin) { ir=0; }
894 else {
895 if (radius > m_Rmax) radius=m_Rmax-0.0001;
896 ir=((int) ((radius-m_Rmin)/m_dR) );
897 }
898 return m_Rphi[ir];
899 }

◆ PhiGap()

G4int LArG4::Barrel::Geometry::PhiGap ( const double & radius,
const double & xhit,
const double & yhit ) const
private

Definition at line 904 of file LArBarrelGeometry.cxx.

905 {
906 const G4double phi0=Phi0(radius)+m_gam0; // from -pi to pi
907 const G4double phi_hit=atan2(yhit,xhit); // from -pi to pi
908 G4double dphi=phi_hit-phi0;
909 // bring back to 0-2pi
910 if (dphi<0) dphi=dphi+2*M_PI;
911 if (dphi>=2*M_PI) dphi=dphi-2*M_PI;
912 dphi=dphi/(2*M_PI)*1024;
913 const G4int ngap=((int) dphi);
914#ifdef DEBUGHITS
915 ATH_MSG_VERBOSE(" phi0 " << phi0 << " dphi, ngap " << dphi << " " << ngap);
916#endif
917 return ngap;
918 }
G4double Phi0(G4double) const

◆ SampSeg()

G4int LArG4::Barrel::Geometry::SampSeg ( G4double eta,
G4double radius,
G4double z,
G4int & iregion,
G4int & isampling,
G4int & ieta,
G4int & isamp2,
G4int & ieta2 ) const
private

Definition at line 294 of file LArBarrelGeometry.cxx.

297 {
298 // Helper struct to hold pre-calculated geometry information
299 struct Geo {
300 G4double Eta_max,R_max_acc,Z_max_acc,R_min_acc,Z_max_lowr,dzdr;
301 G4double zmax1,zmax2,zmax3,zmax4,zmax5,zmax6,zmax7,rmax1,rmax2,rmax3,rmax4;
302 };
303
304 // Fill it once
305 static const Geo g = []() {
307 Geo g{};
308
309 // maximum eta barrel 1.475 (at r=1500.024)
310 g.Eta_max = parameters->GetValue("LArEMBMaxEtaAcceptance");
311 // minimum active radius 1500.024
312 g.R_min_acc= parameters->GetValue("LArEMBRadiusInnerAccordion");
313 // maximum active radius 1960.
314 g.R_max_acc = parameters->GetValue("LArEMBFiducialRmax");
315 // maximum active z (before subtracting edge for signal readout)
316 // currently 3150, should be changed in database to become 3164
317 g.Z_max_acc = parameters->GetValue("LArEMBfiducialMothZmax");
318 // minimum radius at z max for active region
319 const G4double R_min_highz=1548.; //FIXME should be taken from database
320
321 // compute z edge taken by readout strips on the edge
322 const G4double deltaz=7.; // this include copper 6mm + 2*0.5mm empty space on each side
323 g.zmax1=g.Z_max_acc-deltaz;
324 g.zmax2=g.Z_max_acc-2.*deltaz;
325 g.zmax3=g.Z_max_acc-3.*deltaz;
326 g.zmax4=g.Z_max_acc-4.*deltaz;
327 g.zmax5=g.Z_max_acc-5.*deltaz;
328 g.zmax6=g.Z_max_acc-6.*deltaz;
329 g.zmax7=g.Z_max_acc-7.*deltaz;
330 g.rmax1=g.R_max_acc-deltaz;
331 g.rmax2=g.R_max_acc-2.*deltaz;
332 g.rmax3=g.R_max_acc-3.*deltaz;
333 g.rmax4=g.R_max_acc-4.*deltaz;
334
335 // maximum z at r=1500.024
336 g.Z_max_lowr = sinh(g.Eta_max)*g.R_min_acc;
337 // slope of z vs r edge (which is not projective in eta...)
338 g.dzdr = (g.Z_max_acc-g.Z_max_lowr)/(R_min_highz-g.R_min_acc);
339
340 // ATH_MSG_VERBOSE("Initialization of SampSet ");
341 // ATH_MSG_VERBOSE(" Rmin/Rmax " << g.R_min_acc << " " << g.R_max_acc);
342 // ATH_MSG_VERBOSE(" Zmax/Zmax_lowR " << g.Z_max_acc << " " << g.Z_max_lowr);
343 // ATH_MSG_VERBOSE(" Rmin_highz " << R_min_highz);
344 // ATH_MSG_VERBOSE(" dzdr " << g.dzdr);
345 return g;
346 }();
347
348 // inactive thickness between S1 and S2 FIXME should be taken from database
349 const G4double Dr_s12=1.1;
350 const G4double Eta_max_s1=1.4; // maximum eta region 0
351 const G4double Eta_max_s3=1.325; // maximum eta for S3 in region 0
352 const G4double deta=0.025; // basic granularity
353
354
355 // iactive = 1 if active region, 0 if region considered as inactive
356 G4int iactive = 1;
357
358 const G4double aeta=std::fabs(eta);
359
360 G4double r12=-1.;
361 G4double r23=-1.;
362
363 // Not used: G4double rmax=Z_max_acc/sinh(aeta);
364
365 G4int istrip,imid;
366
367 // eta < 1.4
368
369 if (aeta<Eta_max_s1) {
370
371 // get radius for end of strips
372 istrip=(int) (aeta/deta*8.);
373 if (istrip<0 || istrip >=448) {
374 ATH_MSG_ERROR(" Problem aeta,istrip " << aeta << " " << istrip);
375 return 0;
376 }
377 r12=Rmax1[istrip];
378
379 // get radius for end of middle
380 imid = (int) (aeta/deta);
381 if (imid <0 || imid >=56) {
382 ATH_MSG_ERROR(" Problem aeta,imid " << aeta << " " << imid);
383 return 0;
384 }
385 r23=Rmax2[imid];
386
387 iregion=0;
388
389 // strips
390 if (radius <= r12) {
391 isampling=1;
392 ieta=istrip;
393 if (ieta==0) iactive=0;
394 isamp2=1;
395 ieta2=istrip;
396 }
397
398 // region between strips and middle => not active, same identifier as strips
399 else if (radius < (r12+Dr_s12)) {
400 isampling=1;
401 ieta=istrip;
402 iactive=0;
403 isamp2=1;
404 ieta2=istrip;
405 }
406
407 else {
408
409 // eta<1.325, we can be in the back
410 if (aeta<Eta_max_s3) {
411 // radius<r23 we are in the middle
412 if (radius <= r23) {
413 isampling=2;
414 ieta=imid;
415 isamp2=2;
416 ieta2=imid;
417 }
418 // for radius >r23 we have to take care of the readout strips at high z
419 // and attribute some of the energy to other cells
420 else { // radius>r23
421 if (z>g.zmax1) {
422 isampling=2;
423 ieta=55;
424 }
425 else if (z>g.zmax2) {
426 isampling=2;
427 ieta=54;
428 }
429 else if (z>g.zmax3) {
430 isampling=2;
431 ieta=53;
432 }
433 else if (z>g.zmax4) {
434 isampling=3;
435 ieta=26;
436 }
437 else if (aeta<1.3 && z>g.zmax5) {
438 isampling=2;
439 ieta=52;
440 }
441 else if (aeta<1.3 && z>g.zmax6) {
442 isampling=2;
443 ieta=51;
444 }
445 else if (radius>g.rmax4 && z<g.zmax5 && aeta>1.2) {
446 if (radius>g.rmax1) {
447 isampling=2;
448 ieta=51;
449 }
450 else if(radius>g.rmax2) {
451 isampling=3;
452 ieta=25;
453 }
454 else if (radius>g.rmax3) {
455 if (z<g.zmax7) {
456 isampling=2;
457 ieta=50;
458 }
459 else {
460 isampling=3;
461 ieta=25;
462 }
463 }
464 else {
465 if (aeta<1.25) {
466 isampling=2;
467 ieta=49;
468 }
469 else {
470 isampling=3;
471 ieta=25;
472 }
473 }
474 }
475 // normal back cell
476 else {
477 isampling=3;
478 ieta=imid/2;
479 isamp2=3;
480 ieta2=ieta;
481 }
482 isamp2=3;
483 ieta2=imid/2;
484 } // end radius>r23
485 // put into middle energy deposited along readout strips across the back
486 if (isampling==3 && z<g.zmax4 && (radius<g.rmax4 || aeta<1.2) ) {
487 const double etastr = (imid%2==0) ? 0.025*imid : 0.025*(imid+1);
488 const double delta=radius*(sinh(etastr)-sinh(aeta))/cosh(etastr);
489 double deltastr;
490 if (aeta<0.475) { deltastr=1.5;}
491 else if (aeta<0.80) { deltastr=2.75;}
492 else if (aeta<0.85) { deltastr=1.5;}
493 else if (aeta<1.1) { deltastr=2.75;}
494 else { deltastr=3.25;}
495
496 if (std::fabs(delta)<deltastr) {
497 isampling=2;
498 ieta=imid;
499 }
500 } // end if sampling==3
501 } // end if eta<1.325
502 else {
503 isampling=2;
504 ieta=imid;
505 if (z>g.zmax1) {
506 ieta=55;
507 }
508 else if (z>g.zmax2 && aeta<1.375) {
509 ieta=54;
510 }
511 isamp2=2;
512 ieta2=imid;
513 } // end eta>1.352
514 } // end radius selection
515 } // end eta1.4
516
517 // eta between 1.4 and 1.475
518
519 if (aeta>=Eta_max_s1 && aeta<g.Eta_max) {
520 r12 = Rmax1[447]; // radius for end of sampling 1
521 r23=g.Z_max_acc/sinh(aeta); // radius of end of sampling 2, bounded by high z end
522
523 const double zmax = g.Z_max_lowr + g.dzdr*(radius-g.R_min_acc);
524
525 iregion=1;
526 if (radius <=r12) {
527 isampling=1;
528 ieta=int((aeta-Eta_max_s1)/deta);
529 if (z>zmax) { iactive=0; }
530 }
531 else if (radius < (r12+Dr_s12)) {
532 isampling=1;
533 ieta=int((aeta-Eta_max_s1)/deta);
534 iactive=0;
535 }
536 else if (radius <= r23) {
537 isampling=2;
538 ieta=0;
539 if (z>zmax) { iactive=0; }
540 }
541 else {
542 isampling=2;
543 ieta=0;
544 iactive=0;
545 }
546 isamp2=isampling;
547 ieta2=ieta;
548 }
549 // eta above 1.475, not fiducial region, but still returns something
550 // for calibration hits
551 if (aeta>g.Eta_max) {
552 iregion=1;
553 r12 = Rmax1[447];
554 if (radius <=r12) {
555 isampling=1;
556 ieta=2;
557 }
558 else {
559 isampling=2;
560 ieta=0;
561 }
562 isamp2=isampling;
563 ieta2=ieta;
564 iactive=0;
565 }
566
567 // cross-check of active region
568 if (z>g.Z_max_acc || radius>g.R_max_acc || radius<g.R_min_acc || aeta > g.Eta_max) iactive=0;
569
570 return iactive;
571 }

Member Data Documentation

◆ m_absorber

const LArStraightAbsorbers* LArG4::Barrel::Geometry::m_absorber {nullptr}
private

Definition at line 106 of file LArBarrelGeometry.h.

106{nullptr};

◆ m_coudeabs

const LArCoudeAbsorbers* LArG4::Barrel::Geometry::m_coudeabs {nullptr}
private

Definition at line 104 of file LArBarrelGeometry.h.

104{nullptr};

◆ m_coudeelec

const LArCoudeElectrodes* LArG4::Barrel::Geometry::m_coudeelec {nullptr}
private

Definition at line 103 of file LArBarrelGeometry.h.

103{nullptr};

◆ m_delta

double* LArG4::Barrel::Geometry::m_delta {nullptr}
private

Definition at line 99 of file LArBarrelGeometry.h.

99{nullptr};

◆ m_detectorName

Gaudi::Property<std::string> LArG4::Barrel::Geometry::m_detectorName {this, "DetectorName", "LArMgr"}
private

Definition at line 68 of file LArBarrelGeometry.h.

68{this, "DetectorName", "LArMgr"};

◆ m_dR

G4double LArG4::Barrel::Geometry::m_dR {0.}
private

Definition at line 115 of file LArBarrelGeometry.h.

115{0.};

◆ m_ecamName

G4String LArG4::Barrel::Geometry::m_ecamName
private

Definition at line 73 of file LArBarrelGeometry.h.

◆ m_electrode

const LArStraightElectrodes* LArG4::Barrel::Geometry::m_electrode {nullptr}
private

Definition at line 105 of file LArBarrelGeometry.h.

105{nullptr};

◆ m_etaMaxBarrel

double LArG4::Barrel::Geometry::m_etaMaxBarrel {0.}
private

Definition at line 81 of file LArBarrelGeometry.h.

81{0.};

◆ m_gam0

double LArG4::Barrel::Geometry::m_gam0 {0.}
private

Definition at line 92 of file LArBarrelGeometry.h.

92{0.}; //phi position for the first absorber neutral fiber

◆ m_iflSAG

bool LArG4::Barrel::Geometry::m_iflSAG {false}
private

Definition at line 108 of file LArBarrelGeometry.h.

108{false};

◆ m_Nbrt

int LArG4::Barrel::Geometry::m_Nbrt {0}
private

Definition at line 88 of file LArBarrelGeometry.h.

88{0}; // number of straight sections (=14)

◆ m_Nbrt1

int LArG4::Barrel::Geometry::m_Nbrt1 {0}
private

Definition at line 89 of file LArBarrelGeometry.h.

89{0}; // number of folds (=15)

◆ m_NCellMax

int LArG4::Barrel::Geometry::m_NCellMax {0}
private

Definition at line 85 of file LArBarrelGeometry.h.

85{0}; // 1024

◆ m_NCellTot

int LArG4::Barrel::Geometry::m_NCellTot {0}
private

Definition at line 84 of file LArBarrelGeometry.h.

84{0}; // either 64 or 1024 for TestBeam or Atlas

◆ m_NRphi

G4int LArG4::Barrel::Geometry::m_NRphi {0}
private

Definition at line 111 of file LArBarrelGeometry.h.

111{0};

◆ m_parity

int LArG4::Barrel::Geometry::m_parity {0}
private

Definition at line 100 of file LArBarrelGeometry.h.

100{0};

◆ m_phic

double* LArG4::Barrel::Geometry::m_phic {nullptr}
private

Definition at line 96 of file LArBarrelGeometry.h.

96{nullptr};

◆ m_rc

double* LArG4::Barrel::Geometry::m_rc {nullptr}
private

Definition at line 95 of file LArBarrelGeometry.h.

95{nullptr};

◆ m_rint_eleFib

double LArG4::Barrel::Geometry::m_rint_eleFib {0.}
private

Definition at line 93 of file LArBarrelGeometry.h.

93{0.}; //2.78

◆ m_Rmax

G4double LArG4::Barrel::Geometry::m_Rmax {0.}
private

Definition at line 113 of file LArBarrelGeometry.h.

113{0.};

◆ m_rMaxAccordion

double LArG4::Barrel::Geometry::m_rMaxAccordion {0.}
private

Definition at line 77 of file LArBarrelGeometry.h.

77{0.};

◆ m_Rmin

G4double LArG4::Barrel::Geometry::m_Rmin {0.}
private

Definition at line 112 of file LArBarrelGeometry.h.

112{0.};

◆ m_rMinAccordion

double LArG4::Barrel::Geometry::m_rMinAccordion {0.}
private

Definition at line 76 of file LArBarrelGeometry.h.

76{0.};

◆ m_Rphi

G4double LArG4::Barrel::Geometry::m_Rphi[5000] = {0}
private

Definition at line 114 of file LArBarrelGeometry.h.

114{0};

◆ m_testbeam

Gaudi::Property<bool> LArG4::Barrel::Geometry::m_testbeam {this, "TestBeam", false}
private

Definition at line 71 of file LArBarrelGeometry.h.

71{this, "TestBeam", false};

◆ m_xc

double* LArG4::Barrel::Geometry::m_xc {nullptr}
private

Definition at line 97 of file LArBarrelGeometry.h.

97{nullptr};

◆ m_yc

double* LArG4::Barrel::Geometry::m_yc {nullptr}
private

Definition at line 98 of file LArBarrelGeometry.h.

98{nullptr};

◆ m_zMaxBarrel

double LArG4::Barrel::Geometry::m_zMaxBarrel {0.}
private

Definition at line 79 of file LArBarrelGeometry.h.

79{0.};

◆ m_zMaxBarrelDMMargin

double LArG4::Barrel::Geometry::m_zMaxBarrelDMMargin {0.}
private

Definition at line 80 of file LArBarrelGeometry.h.

80{0.};

◆ m_zMinBarrel

double LArG4::Barrel::Geometry::m_zMinBarrel {0.}
private

Definition at line 78 of file LArBarrelGeometry.h.

78{0.};

The documentation for this class was generated from the following files: