File indexing completed on 2024-04-06 12:32:15
0001 #include "Validation/Geometry/interface/MaterialBudgetData.h"
0002
0003 #include "G4Step.hh"
0004 #include "G4Material.hh"
0005 #include "G4EventManager.hh"
0006 #include "G4Event.hh"
0007 #include "DD4hep/Filter.h"
0008
0009 MaterialBudgetData::MaterialBudgetData() {
0010
0011 myMaterialBudgetCategorizer = nullptr;
0012 allStepsToTree = false;
0013 isHGCal = false;
0014 isMtd = false;
0015 densityConvertionFactor = 6.24E18;
0016 }
0017
0018 MaterialBudgetData::~MaterialBudgetData() = default;
0019
0020 void MaterialBudgetData::SetAllStepsToTree() { allStepsToTree = true; }
0021
0022 void MaterialBudgetData::dataStartTrack(const G4Track* aTrack) {
0023 const G4ThreeVector& dir = aTrack->GetMomentum();
0024
0025 if (myMaterialBudgetCategorizer == nullptr) {
0026 if (isHGCal) {
0027 myMaterialBudgetCategorizer = std::make_unique<MaterialBudgetCategorizer>("HGCal");
0028 } else if (isMtd) {
0029 myMaterialBudgetCategorizer = std::make_unique<MaterialBudgetCategorizer>("Mtd");
0030 } else {
0031 myMaterialBudgetCategorizer = std::make_unique<MaterialBudgetCategorizer>("Tracker");
0032 }
0033 }
0034
0035 theStepN = 0;
0036 theTotalMB = 0;
0037 theTotalIL = 0;
0038 theEta = 0;
0039 thePhi = 0;
0040 theID = 0;
0041 thePt = 0;
0042 theEnergy = 0;
0043 theMass = 0;
0044
0045 theSupportMB = 0.f;
0046 theSensitiveMB = 0.f;
0047 theCoolingMB = 0.f;
0048 theElectronicsMB = 0.f;
0049 theOtherMB = 0.f;
0050
0051
0052 theAirMB = 0.f;
0053 theCablesMB = 0.f;
0054 theCopperMB = 0.f;
0055 theH_ScintillatorMB = 0.f;
0056 theLeadMB = 0.f;
0057 theEpoxyMB = 0.f;
0058 theKaptonMB = 0.f;
0059 theAluminiumMB = 0.f;
0060 theHGC_G10_FR4MB = 0.f;
0061 theSiliconMB = 0.f;
0062 theStainlessSteelMB = 0.f;
0063 theWCuMB = 0.f;
0064 thePolystyreneMB = 0.f;
0065 theHGC_EEConnectorMB = 0.f;
0066 theHGC_HEConnectorMB = 0.f;
0067
0068 theSupportIL = 0.f;
0069 theSensitiveIL = 0.f;
0070 theCoolingIL = 0.f;
0071 theElectronicsIL = 0.f;
0072 theOtherIL = 0.f;
0073
0074
0075 theAirIL = 0.f;
0076 theCablesIL = 0.f;
0077 theCopperIL = 0.f;
0078 theH_ScintillatorIL = 0.f;
0079 theLeadIL = 0.f;
0080 theEpoxyIL = 0.f;
0081 theKaptonIL = 0.f;
0082 theAluminiumIL = 0.f;
0083 theHGC_G10_FR4IL = 0.f;
0084 theSiliconIL = 0.f;
0085 theStainlessSteelIL = 0.f;
0086 theWCuIL = 0.f;
0087 thePolystyreneIL = 0.f;
0088 theHGC_EEConnectorIL = 0.f;
0089 theHGC_HEConnectorIL = 0.f;
0090
0091 theSupportFractionMB = 0.f;
0092 theSensitiveFractionMB = 0.f;
0093 theCoolingFractionMB = 0.f;
0094 theElectronicsFractionMB = 0.f;
0095 theOtherFractionMB = 0.f;
0096
0097 theAirFractionMB = 0.f;
0098 theCablesFractionMB = 0.f;
0099 theCopperFractionMB = 0.f;
0100 theH_ScintillatorFractionMB = 0.f;
0101 theLeadFractionMB = 0.f;
0102 theEpoxyFractionMB = 0.f;
0103 theKaptonFractionMB = 0.f;
0104 theAluminiumFractionMB = 0.f;
0105 theHGC_G10_FR4FractionMB = 0.f;
0106 theSiliconFractionMB = 0.f;
0107 theStainlessSteelFractionMB = 0.f;
0108 theWCuFractionMB = 0.f;
0109 thePolystyreneFractionMB = 0.f;
0110 theHGC_EEConnectorFractionMB = 0.f;
0111 theHGC_HEConnectorFractionMB = 0.f;
0112
0113 theSupportFractionIL = 0.f;
0114 theSensitiveFractionIL = 0.f;
0115 theCoolingFractionIL = 0.f;
0116 theElectronicsFractionIL = 0.f;
0117 theOtherFractionIL = 0.f;
0118
0119 theAirFractionIL = 0.f;
0120 theCablesFractionIL = 0.f;
0121 theCopperFractionIL = 0.f;
0122 theH_ScintillatorFractionIL = 0.f;
0123 theLeadFractionIL = 0.f;
0124 theEpoxyFractionIL = 0.f;
0125 theKaptonFractionIL = 0.f;
0126 theAluminiumFractionIL = 0.f;
0127 theHGC_G10_FR4FractionIL = 0.f;
0128 theSiliconFractionIL = 0.f;
0129 theStainlessSteelFractionIL = 0.f;
0130 theWCuFractionIL = 0.f;
0131 thePolystyreneFractionIL = 0.f;
0132 theHGC_EEConnectorFractionIL = 0.f;
0133 theHGC_HEConnectorFractionIL = 0.f;
0134
0135 theID = (int)(aTrack->GetDefinition()->GetPDGEncoding());
0136 thePt = dir.perp();
0137 if (dir.theta() != 0) {
0138 theEta = dir.eta();
0139 } else {
0140 theEta = -99;
0141 }
0142 thePhi = dir.phi();
0143 theEnergy = aTrack->GetTotalEnergy();
0144 theMass = aTrack->GetDefinition()->GetPDGMass();
0145 }
0146
0147 void MaterialBudgetData::dataEndTrack(const G4Track* aTrack) {
0148 LogDebug("MaterialBudget") << "MaterialBudgetData: [OVAL] MaterialBudget "
0149 << G4EventManager::GetEventManager()->GetConstCurrentEvent()->GetEventID()
0150 << " Eta:" << theEta << " Phi:" << thePhi << " TotalMB" << theTotalMB;
0151
0152 LogDebug("MaterialBudget") << "MaterialBudgetData:" << theStepN << "Recorded steps ";
0153
0154 if (!isHGCal) {
0155 LogDebug("Material Budget") << "MaterialBudgetData: Radiation Length "
0156 << G4EventManager::GetEventManager()->GetConstCurrentEvent()->GetEventID() << " Eta"
0157 << theEta << " Phi" << thePhi << " TotalMB" << theTotalMB << " SUP " << theSupportMB
0158 << " SEN " << theSensitiveMB << " CAB " << theCablesMB << " COL " << theCoolingMB
0159 << " ELE " << theElectronicsMB << " other " << theOtherMB << " Air " << theAirMB;
0160
0161 LogDebug("Material Budget") << "MaterialBudgetData: Interaction Length "
0162 << G4EventManager::GetEventManager()->GetConstCurrentEvent()->GetEventID() << " Eta "
0163 << theEta << " Phi " << thePhi << " TotalIL " << theTotalIL << " SUP " << theSupportIL
0164 << " SEN " << theSensitiveIL << " CAB " << theCablesIL << " COL " << theCoolingIL
0165 << " ELE " << theElectronicsIL << " other " << theOtherIL << " Air " << theAirIL
0166 << std::endl;
0167
0168 } else {
0169 LogDebug("MaterialBudget") << "MaterialBudgetData: HGCal Material Budget: Radiation Length "
0170 << G4EventManager::GetEventManager()->GetConstCurrentEvent()->GetEventID() << " Eta "
0171 << theEta << " Phi " << thePhi << " TotaLMB" << theTotalMB << " theCopperMB "
0172 << theCopperMB << " theH_ScintillatorMB " << theH_ScintillatorMB << " CAB "
0173 << theCablesMB << " theLeadMB " << theLeadMB << " theEpoxyMB " << theEpoxyMB
0174 << " theKaptonMB " << theKaptonMB << " theAluminiumMB " << theAluminiumMB
0175 << " theHGC_G10_FR4MB " << theHGC_G10_FR4MB << " theSiliconMB " << theSiliconMB
0176 << " theAirMB " << theAirMB << " theStainlessSteelMB " << theStainlessSteelMB
0177 << " theWCuMB " << theWCuMB << " thePolystyreneMB " << thePolystyreneMB
0178 << " theHGC_EEConnectorMB " << theHGC_EEConnectorMB << " theHGC_HEConnectorMB "
0179 << theHGC_HEConnectorMB;
0180
0181 LogDebug("MaterialBudget") << "MaterialBudgetData: HGCal Material Budget: Interaction Length "
0182 << G4EventManager::GetEventManager()->GetConstCurrentEvent()->GetEventID() << " Eta "
0183 << theEta << " Phi " << thePhi << " TotalIL " << theTotalIL << " theCopperIL "
0184 << theCopperIL << " theH_ScintillatorIL " << theH_ScintillatorIL << " CAB "
0185 << theCablesIL << " theLeadIL " << theLeadIL << " theEpoxyIL " << theEpoxyIL
0186 << " theKaptonIL " << theKaptonIL << " theAluminiumIL " << theAluminiumIL
0187 << " theHGC_G10_FR4IL " << theHGC_G10_FR4IL << " theSiliconIL " << theSiliconIL
0188 << " Air " << theAirIL << " theStainlessSteelIL " << theStainlessSteelIL << " theWCuIL "
0189 << theWCuIL << " thePolystyreneIL " << thePolystyreneIL << " theHGC_EEConnectorIL "
0190 << theHGC_EEConnectorIL << " theHGC_HEConnectorIL " << theHGC_HEConnectorIL << std::endl;
0191 }
0192 }
0193
0194 void MaterialBudgetData::dataPerStep(const G4Step* aStep) {
0195 assert(aStep);
0196 G4StepPoint* prePoint = aStep->GetPreStepPoint();
0197 G4StepPoint* postPoint = aStep->GetPostStepPoint();
0198 assert(prePoint);
0199 assert(postPoint);
0200 G4Material* theMaterialPre = prePoint->GetMaterial();
0201 assert(theMaterialPre);
0202
0203 CLHEP::Hep3Vector prePos = prePoint->GetPosition();
0204 CLHEP::Hep3Vector postPos = postPoint->GetPosition();
0205
0206 G4double steplen = aStep->GetStepLength();
0207
0208 G4double radlen = theMaterialPre->GetRadlen();
0209 G4double intlen = theMaterialPre->GetNuclearInterLength();
0210 G4double density = theMaterialPre->GetDensity() / densityConvertionFactor;
0211
0212 G4String materialName = static_cast<std::string>(dd4hep::dd::noNamespace(theMaterialPre->GetName()));
0213
0214 LogDebug("MaterialBudget") << "MaterialBudgetData: Material " << materialName << " steplen " << steplen << " radlen "
0215 << radlen << " mb " << steplen / radlen;
0216
0217 G4String volumeName = static_cast<std::string>(
0218 dd4hep::dd::noNamespace(aStep->GetPreStepPoint()->GetTouchable()->GetVolume(0)->GetLogicalVolume()->GetName()));
0219
0220 LogDebug("MaterialBudget") << "MaterialBudgetData: Volume " << volumeName << " Material " << materialName;
0221
0222
0223 assert(myMaterialBudgetCategorizer);
0224 int volumeID = myMaterialBudgetCategorizer->volume(volumeName);
0225 int materialID = myMaterialBudgetCategorizer->material(materialName);
0226
0227 LogDebug("MaterialBudget") << "MaterialBudgetData: Volume ID " << volumeID << " Material ID " << materialID;
0228
0229
0230
0231
0232 if (!isHGCal) {
0233 bool isCtgOk = !myMaterialBudgetCategorizer->x0fraction(materialName).empty() &&
0234 !myMaterialBudgetCategorizer->l0fraction(materialName).empty() &&
0235 (myMaterialBudgetCategorizer->x0fraction(materialName).size() == 7)
0236 && (myMaterialBudgetCategorizer->l0fraction(materialName).size() == 7);
0237
0238 if (!isCtgOk) {
0239 if (materialName.compare("Air") == 0) {
0240 theAirFractionMB = 1.f;
0241 theAirFractionIL = 1.f;
0242 } else {
0243 theOtherFractionMB = 1.f;
0244 theOtherFractionIL = 1.f;
0245 edm::LogVerbatim("MaterialBudget") << "MaterialBudgetData: Material forced to 'Other': " << materialName
0246 << " in volume " << volumeName << ". Check Categorization.";
0247 }
0248 } else {
0249 theSupportFractionMB = myMaterialBudgetCategorizer->x0fraction(materialName)[0];
0250 theSensitiveFractionMB = myMaterialBudgetCategorizer->x0fraction(materialName)[1];
0251 theCablesFractionMB = myMaterialBudgetCategorizer->x0fraction(materialName)[2];
0252 theCoolingFractionMB = myMaterialBudgetCategorizer->x0fraction(materialName)[3];
0253 theElectronicsFractionMB = myMaterialBudgetCategorizer->x0fraction(materialName)[4];
0254 theOtherFractionMB = myMaterialBudgetCategorizer->x0fraction(materialName)[5];
0255 theAirFractionMB = myMaterialBudgetCategorizer->x0fraction(materialName)[6];
0256
0257 if (theOtherFractionMB > 0.f)
0258 edm::LogVerbatim("MaterialBudget")
0259 << "MaterialBudgetData: Material found with no category: " << materialName << " in volume " << volumeName;
0260
0261 theSupportFractionIL = myMaterialBudgetCategorizer->l0fraction(materialName)[0];
0262 theSensitiveFractionIL = myMaterialBudgetCategorizer->l0fraction(materialName)[1];
0263 theCablesFractionIL = myMaterialBudgetCategorizer->l0fraction(materialName)[2];
0264 theCoolingFractionIL = myMaterialBudgetCategorizer->l0fraction(materialName)[3];
0265 theElectronicsFractionIL = myMaterialBudgetCategorizer->l0fraction(materialName)[4];
0266 theOtherFractionIL = myMaterialBudgetCategorizer->l0fraction(materialName)[5];
0267 theAirFractionIL = myMaterialBudgetCategorizer->l0fraction(materialName)[6];
0268
0269 if (theOtherFractionIL > 0.f)
0270 edm::LogVerbatim("MaterialBudget")
0271 << "MaterialBudgetData: Material found with no category: " << materialName << " in volume " << volumeName;
0272 }
0273 } else {
0274
0275 bool isHGCalx0fractionEmpty = myMaterialBudgetCategorizer->HGCalx0fraction(materialName).empty();
0276 bool isHGCall0fractionEmpty = myMaterialBudgetCategorizer->HGCall0fraction(materialName).empty();
0277
0278 if (isHGCalx0fractionEmpty && isHGCall0fractionEmpty) {
0279 theOtherFractionMB = 1.f;
0280 theOtherFractionIL = 1.f;
0281
0282 edm::LogVerbatim("MaterialBudget") << "MaterialBudgetData: Material forced to 'Other': " << materialName
0283 << " in volume " << volumeName;
0284 } else {
0285 theAirFractionMB = myMaterialBudgetCategorizer->HGCalx0fraction(materialName)[0];
0286 theCablesFractionMB = myMaterialBudgetCategorizer->HGCalx0fraction(materialName)[1];
0287 theCopperFractionMB = myMaterialBudgetCategorizer->HGCalx0fraction(materialName)[2];
0288 theH_ScintillatorFractionMB = myMaterialBudgetCategorizer->HGCalx0fraction(materialName)[3];
0289 theLeadFractionMB = myMaterialBudgetCategorizer->HGCalx0fraction(materialName)[4];
0290 theHGC_G10_FR4FractionMB = myMaterialBudgetCategorizer->HGCalx0fraction(materialName)[5];
0291 theSiliconFractionMB = myMaterialBudgetCategorizer->HGCalx0fraction(materialName)[6];
0292 theStainlessSteelFractionMB = myMaterialBudgetCategorizer->HGCalx0fraction(materialName)[7];
0293 theWCuFractionMB = myMaterialBudgetCategorizer->HGCalx0fraction(materialName)[8];
0294 theOtherFractionMB = myMaterialBudgetCategorizer->HGCalx0fraction(materialName)[9];
0295 theEpoxyFractionMB = myMaterialBudgetCategorizer->HGCalx0fraction(materialName)[10];
0296 theKaptonFractionMB = myMaterialBudgetCategorizer->HGCalx0fraction(materialName)[11];
0297 theAluminiumFractionMB = myMaterialBudgetCategorizer->HGCalx0fraction(materialName)[12];
0298 thePolystyreneFractionMB = myMaterialBudgetCategorizer->HGCalx0fraction(materialName)[13];
0299 theHGC_EEConnectorFractionMB = myMaterialBudgetCategorizer->HGCalx0fraction(materialName)[14];
0300 theHGC_HEConnectorFractionMB = myMaterialBudgetCategorizer->HGCalx0fraction(materialName)[15];
0301
0302 if (theOtherFractionMB != 0)
0303
0304
0305 std::cout << "MaterialBudgetData: Material found with no category: " << materialName << " in volume "
0306 << volumeName << std::endl;
0307
0308 theAirFractionIL = myMaterialBudgetCategorizer->HGCall0fraction(materialName)[0];
0309 theCablesFractionIL = myMaterialBudgetCategorizer->HGCall0fraction(materialName)[1];
0310 theCopperFractionIL = myMaterialBudgetCategorizer->HGCall0fraction(materialName)[2];
0311 theH_ScintillatorFractionIL = myMaterialBudgetCategorizer->HGCall0fraction(materialName)[3];
0312 theLeadFractionIL = myMaterialBudgetCategorizer->HGCall0fraction(materialName)[4];
0313 theHGC_G10_FR4FractionIL = myMaterialBudgetCategorizer->HGCall0fraction(materialName)[5];
0314 theSiliconFractionIL = myMaterialBudgetCategorizer->HGCall0fraction(materialName)[6];
0315 theStainlessSteelFractionIL = myMaterialBudgetCategorizer->HGCall0fraction(materialName)[7];
0316 theWCuFractionIL = myMaterialBudgetCategorizer->HGCall0fraction(materialName)[8];
0317 theOtherFractionIL = myMaterialBudgetCategorizer->HGCall0fraction(materialName)[9];
0318 theEpoxyFractionIL = myMaterialBudgetCategorizer->HGCall0fraction(materialName)[10];
0319 theKaptonFractionIL = myMaterialBudgetCategorizer->HGCall0fraction(materialName)[11];
0320 theAluminiumFractionIL = myMaterialBudgetCategorizer->HGCall0fraction(materialName)[12];
0321 thePolystyreneFractionIL = myMaterialBudgetCategorizer->HGCall0fraction(materialName)[13];
0322 theHGC_EEConnectorFractionIL = myMaterialBudgetCategorizer->HGCall0fraction(materialName)[14];
0323 theHGC_HEConnectorFractionIL = myMaterialBudgetCategorizer->HGCall0fraction(materialName)[15];
0324
0325 if (theOtherFractionIL != 0)
0326 edm::LogVerbatim("MaterialBudget") << "MaterialBudgetData: Material found with no category " << materialName
0327 << " in volume " << volumeName << std::endl;
0328 }
0329 }
0330
0331 float dmb = steplen / radlen;
0332 float dil = steplen / intlen;
0333
0334 G4VPhysicalVolume* pv = aStep->GetPreStepPoint()->GetPhysicalVolume();
0335 const G4VTouchable* t = aStep->GetPreStepPoint()->GetTouchable();
0336 const G4ThreeVector& objectTranslation = t->GetTranslation();
0337 const G4RotationMatrix* objectRotation = t->GetRotation();
0338 const G4VProcess* interactionPre = prePoint->GetProcessDefinedStep();
0339 const G4VProcess* interactionPost = postPoint->GetProcessDefinedStep();
0340
0341 G4Track* track = aStep->GetTrack();
0342 if (theStepN == 0)
0343 LogDebug("MaterialBudget") << "MaterialBudgetData: Simulated Particle " << theID << "\tMass " << theMass
0344 << " MeV/c2"
0345 << "\tPt = " << thePt << " MeV/c"
0346 << "\tEta = " << theEta << "\tPhi = " << thePhi << "\tEnergy = " << theEnergy << " MeV";
0347
0348
0349 if (allStepsToTree) {
0350 assert(theStepN < MAXNUMBERSTEPS);
0351 if (theStepN > MAXNUMBERSTEPS)
0352 theStepN = MAXNUMBERSTEPS - 1;
0353 theDmb[theStepN] = dmb;
0354 theDil[theStepN] = dil;
0355 theSupportDmb[theStepN] = (dmb * theSupportFractionMB);
0356 theSensitiveDmb[theStepN] = (dmb * theSensitiveFractionMB);
0357 theCoolingDmb[theStepN] = (dmb * theCoolingFractionMB);
0358 theElectronicsDmb[theStepN] = (dmb * theElectronicsFractionMB);
0359 theOtherDmb[theStepN] = (dmb * theOtherFractionMB);
0360
0361 theAirDmb[theStepN] = (dmb * theAirFractionMB);
0362 theCablesDmb[theStepN] = (dmb * theCablesFractionMB);
0363 theCopperDmb[theStepN] = (dmb * theCopperFractionMB);
0364 theH_ScintillatorDmb[theStepN] = (dmb * theH_ScintillatorFractionMB);
0365 theLeadDmb[theStepN] = (dmb * theLeadFractionMB);
0366 theEpoxyDmb[theStepN] = (dmb * theEpoxyFractionMB);
0367 theKaptonDmb[theStepN] = (dmb * theKaptonFractionMB);
0368 theAluminiumDmb[theStepN] = (dmb * theAluminiumFractionMB);
0369 theHGC_G10_FR4Dmb[theStepN] = (dmb * theHGC_G10_FR4FractionMB);
0370 theSiliconDmb[theStepN] = (dmb * theSiliconFractionMB);
0371 theStainlessSteelDmb[theStepN] = (dmb * theStainlessSteelFractionMB);
0372 theWCuDmb[theStepN] = (dmb * theWCuFractionMB);
0373 thePolystyreneDmb[theStepN] = (dmb * thePolystyreneFractionMB);
0374 theHGC_EEConnectorDmb[theStepN] = (dmb * theHGC_EEConnectorFractionMB);
0375 theHGC_HEConnectorDmb[theStepN] = (dmb * theHGC_HEConnectorFractionMB);
0376
0377 theSupportDil[theStepN] = (dil * theSupportFractionIL);
0378 theSensitiveDil[theStepN] = (dil * theSensitiveFractionIL);
0379 theCoolingDil[theStepN] = (dil * theCoolingFractionIL);
0380 theElectronicsDil[theStepN] = (dil * theElectronicsFractionIL);
0381 theOtherDil[theStepN] = (dil * theOtherFractionIL);
0382
0383 theAirDil[theStepN] = (dil * theAirFractionIL);
0384 theCablesDil[theStepN] = (dil * theCablesFractionIL);
0385 theCopperDil[theStepN] = (dil * theCopperFractionIL);
0386 theH_ScintillatorDil[theStepN] = (dil * theH_ScintillatorFractionIL);
0387 theLeadDil[theStepN] = (dil * theLeadFractionIL);
0388 theEpoxyDil[theStepN] = (dil * theEpoxyFractionIL);
0389 theKaptonDil[theStepN] = (dil * theKaptonFractionIL);
0390 theAluminiumDil[theStepN] = (dil * theAluminiumFractionIL);
0391 theHGC_G10_FR4Dil[theStepN] = (dil * theHGC_G10_FR4FractionIL);
0392 theSiliconDil[theStepN] = (dil * theSiliconFractionIL);
0393 theStainlessSteelDil[theStepN] = (dil * theStainlessSteelFractionIL);
0394 theWCuDil[theStepN] = (dil * theWCuFractionIL);
0395 thePolystyreneDil[theStepN] = (dil * thePolystyreneFractionIL);
0396 theHGC_EEConnectorDil[theStepN] = (dil * theHGC_EEConnectorFractionIL);
0397 theHGC_HEConnectorDil[theStepN] = (dil * theHGC_HEConnectorFractionIL);
0398
0399 theInitialX[theStepN] = prePos.x();
0400 theInitialY[theStepN] = prePos.y();
0401 theInitialZ[theStepN] = prePos.z();
0402 theFinalX[theStepN] = postPos.x();
0403 theFinalY[theStepN] = postPos.y();
0404 theFinalZ[theStepN] = postPos.z();
0405 theVolumeID[theStepN] = volumeID;
0406 theVolumeName[theStepN] = volumeName;
0407 theVolumeCopy[theStepN] = pv->GetCopyNo();
0408 theVolumeX[theStepN] = objectTranslation.x();
0409 theVolumeY[theStepN] = objectTranslation.y();
0410 theVolumeZ[theStepN] = objectTranslation.z();
0411 theVolumeXaxis1[theStepN] = objectRotation->xx();
0412 theVolumeXaxis2[theStepN] = objectRotation->xy();
0413 theVolumeXaxis3[theStepN] = objectRotation->xz();
0414 theVolumeYaxis1[theStepN] = objectRotation->yx();
0415 theVolumeYaxis2[theStepN] = objectRotation->yy();
0416 theVolumeYaxis3[theStepN] = objectRotation->yz();
0417 theVolumeZaxis1[theStepN] = objectRotation->zx();
0418 theVolumeZaxis2[theStepN] = objectRotation->zy();
0419 theVolumeZaxis3[theStepN] = objectRotation->zz();
0420 theMaterialID[theStepN] = materialID;
0421 theMaterialName[theStepN] = materialName;
0422 theMaterialX0[theStepN] = radlen;
0423 theMaterialLambda0[theStepN] = intlen;
0424 theMaterialDensity[theStepN] = density;
0425 theStepID[theStepN] = track->GetDefinition()->GetPDGEncoding();
0426 theStepInitialPt[theStepN] = prePoint->GetMomentum().perp();
0427 theStepInitialEta[theStepN] = prePoint->GetMomentum().eta();
0428 theStepInitialPhi[theStepN] = prePoint->GetMomentum().phi();
0429 theStepInitialEnergy[theStepN] = prePoint->GetTotalEnergy();
0430 theStepInitialPx[theStepN] = prePoint->GetMomentum().x();
0431 theStepInitialPy[theStepN] = prePoint->GetMomentum().y();
0432 theStepInitialPz[theStepN] = prePoint->GetMomentum().z();
0433 theStepInitialBeta[theStepN] = prePoint->GetBeta();
0434 theStepInitialGamma[theStepN] = prePoint->GetGamma();
0435 theStepInitialMass[theStepN] = prePoint->GetMass();
0436 theStepFinalPt[theStepN] = postPoint->GetMomentum().perp();
0437 theStepFinalEta[theStepN] = postPoint->GetMomentum().eta();
0438 theStepFinalPhi[theStepN] = postPoint->GetMomentum().phi();
0439 theStepFinalEnergy[theStepN] = postPoint->GetTotalEnergy();
0440 theStepFinalPx[theStepN] = postPoint->GetMomentum().x();
0441 theStepFinalPy[theStepN] = postPoint->GetMomentum().y();
0442 theStepFinalPz[theStepN] = postPoint->GetMomentum().z();
0443 theStepFinalBeta[theStepN] = postPoint->GetBeta();
0444 theStepFinalGamma[theStepN] = postPoint->GetGamma();
0445 theStepFinalMass[theStepN] = postPoint->GetMass();
0446 int preProcType = -99;
0447 int postProcType = -99;
0448 if (interactionPre)
0449 preProcType = interactionPre->GetProcessType();
0450 theStepPreProcess[theStepN] = preProcType;
0451 if (interactionPost)
0452 postProcType = interactionPost->GetProcessType();
0453 theStepPostProcess[theStepN] = postProcType;
0454
0455 LogDebug("MaterialBudget") << "MaterialBudgetData: Step " << theStepN << "\tDelta MB = " << theDmb[theStepN]
0456 << std::endl
0457 << " Support " << theSupportDmb[theStepN] << " Sensitive " << theSensitiveDmb[theStepN]
0458 << " Cables " << theCablesDmb[theStepN] << " Cooling " << theCoolingDmb[theStepN]
0459 << " Electronics " << theElectronicsDmb[theStepN] << " Other " << theOtherDmb[theStepN]
0460 << " Air " << theAirDmb[theStepN] << std::endl
0461 << "\tDelta IL = " << theDil[theStepN] << std::endl
0462 << " Support " << theSupportDil[theStepN] << " Sensitive " << theSensitiveDil[theStepN]
0463 << " Cables " << theCablesDil[theStepN] << " Cooling " << theCoolingDil[theStepN]
0464 << " Electronics " << theElectronicsDil[theStepN] << " Other " << theOtherDil[theStepN]
0465 << " Air " << theAirDil[theStepN];
0466
0467 if (interactionPre)
0468 LogDebug("MaterialBudget") << "MaterialBudgetData: Process Pre " << interactionPre->GetProcessName() << " type "
0469 << theStepPreProcess[theStepN] << " name "
0470 << interactionPre->GetProcessTypeName(G4ProcessType(theStepPreProcess[theStepN]));
0471 if (interactionPost)
0472 LogDebug("MaterialBudget")
0473 << "MaterialBudgetData: Process Post " << interactionPost->GetProcessName() << " type "
0474 << theStepPostProcess[theStepN] << " name "
0475 << interactionPost->GetProcessTypeName(G4ProcessType(theStepPostProcess[theStepN]))
0476 << " Pre x = " << theInitialX[theStepN] << " y = " << theInitialY[theStepN]
0477 << " z = " << theInitialZ[theStepN]
0478 << " Polar Radius = " << sqrt(prePos.x() * prePos.x() + prePos.y() * prePos.y())
0479 << " Pt = " << theStepInitialPt[theStepN] << " Energy = " << theStepInitialEnergy[theStepN] << " Final: "
0480 << " Post x = " << theFinalX[theStepN] << " y = " << theFinalY[theStepN] << " z = " << theFinalZ[theStepN]
0481 << " Polar Radius = " << sqrt(postPos.x() * postPos.x() + postPos.y() * postPos.y())
0482 << " Pt = " << theStepFinalPt[theStepN] << " Energy = " << theStepFinalEnergy[theStepN] << std::endl
0483 << " Volume " << volumeID << " name " << theVolumeName[theStepN] << " copy number " << theVolumeCopy[theStepN]
0484 << " material " << theMaterialID[theStepN] << " " << theMaterialName[theStepN]
0485 << " Density = " << theMaterialDensity[theStepN] << " g/cm3"
0486 << " X0 = " << theMaterialX0[theStepN] << " mm"
0487 << " Lambda0 = " << theMaterialLambda0[theStepN] << " mm" << std::endl
0488 << " Particle " << theStepID[theStepN] << " Initial Pt = " << theStepInitialPt[theStepN] << " MeV/c"
0489 << " eta = " << theStepInitialEta[theStepN] << " phi = " << theStepInitialPhi[theStepN]
0490 << " Energy = " << theStepInitialEnergy[theStepN] << " MeV"
0491 << " Mass = " << theStepInitialMass[theStepN] << " MeV/c2"
0492 << " Beta = " << theStepInitialBeta[theStepN] << " Gamma = " << theStepInitialGamma[theStepN] << std::endl
0493 << " Particle " << theStepID[theStepN] << " Final Pt = " << theStepFinalPt[theStepN] << " MeV/c"
0494 << " eta = " << theStepFinalEta[theStepN] << " phi = " << theStepFinalPhi[theStepN]
0495 << " Energy = " << theStepFinalEnergy[theStepN] << " MeV"
0496 << " Mass = " << theStepFinalMass[theStepN] << " MeV/c2"
0497 << " Beta = " << theStepFinalBeta[theStepN] << " Gamma = " << theStepFinalGamma[theStepN] << std::endl
0498 << " Volume Centre x = " << theVolumeX[theStepN] << " y = " << theVolumeY[theStepN]
0499 << " z = " << theVolumeZ[theStepN] << "Polar Radius = "
0500 << sqrt(theVolumeX[theStepN] * theVolumeX[theStepN] + theVolumeY[theStepN] * theVolumeY[theStepN])
0501 << std::endl
0502 << " x axis = (" << theVolumeXaxis1[theStepN] << "," << theVolumeXaxis2[theStepN] << ","
0503 << theVolumeXaxis3[theStepN] << ")" << std::endl
0504 << " y axis = (" << theVolumeYaxis1[theStepN] << "," << theVolumeYaxis2[theStepN] << ","
0505 << theVolumeYaxis3[theStepN] << ")" << std::endl
0506 << " z axis = (" << theVolumeZaxis1[theStepN] << "," << theVolumeZaxis2[theStepN] << ","
0507 << theVolumeZaxis3[theStepN] << ")" << std::endl
0508 << " Secondaries" << std::endl;
0509
0510 for (G4TrackVector::const_iterator iSec = aStep->GetSecondary()->begin(); iSec != aStep->GetSecondary()->end();
0511 iSec++) {
0512 LogDebug("MaterialBudget") << "MaterialBudgetData: tid " << (*iSec)->GetDefinition()->GetPDGEncoding()
0513 << " created through process type " << (*iSec)->GetCreatorProcess()->GetProcessType()
0514 << (*iSec)->GetCreatorProcess()->GetProcessTypeName(
0515 G4ProcessType((*iSec)->GetCreatorProcess()->GetProcessType()));
0516 }
0517 }
0518
0519 theTrkLen = aStep->GetTrack()->GetTrackLength();
0520 thePVname = static_cast<std::string>(dd4hep::dd::noNamespace(pv->GetName()));
0521 thePVcopyNo = pv->GetCopyNo();
0522 theRadLen = radlen;
0523 theIntLen = intlen;
0524 theTotalMB += dmb;
0525 theTotalIL += dil;
0526
0527 theSupportMB += (dmb * theSupportFractionMB);
0528 theSensitiveMB += (dmb * theSensitiveFractionMB);
0529 theCoolingMB += (dmb * theCoolingFractionMB);
0530 theElectronicsMB += (dmb * theElectronicsFractionMB);
0531 theOtherMB += (dmb * theOtherFractionMB);
0532
0533
0534 theAirMB += (dmb * theAirFractionMB);
0535 theCablesMB += (dmb * theCablesFractionMB);
0536 theCopperMB += (dmb * theCopperFractionMB);
0537 theH_ScintillatorMB += (dmb * theH_ScintillatorFractionMB);
0538 theLeadMB += (dmb * theLeadFractionMB);
0539 theEpoxyMB += (dmb * theEpoxyFractionMB);
0540 theKaptonMB += (dmb * theKaptonFractionMB);
0541 theAluminiumMB += (dmb * theAluminiumFractionMB);
0542 theHGC_G10_FR4MB += (dmb * theHGC_G10_FR4FractionMB);
0543 theSiliconMB += (dmb * theSiliconFractionMB);
0544 theStainlessSteelMB += (dmb * theStainlessSteelFractionMB);
0545 theWCuMB += (dmb * theWCuFractionMB);
0546 thePolystyreneMB += (dmb * thePolystyreneFractionMB);
0547 theHGC_EEConnectorMB += (dmb * theHGC_EEConnectorFractionMB);
0548 theHGC_HEConnectorMB += (dmb * theHGC_HEConnectorFractionMB);
0549
0550 theSupportIL += (dil * theSupportFractionIL);
0551 theSensitiveIL += (dil * theSensitiveFractionIL);
0552 theCoolingIL += (dil * theCoolingFractionIL);
0553 theElectronicsIL += (dil * theElectronicsFractionIL);
0554 theOtherIL += (dil * theOtherFractionIL);
0555
0556 theAirIL += (dil * theAirFractionIL);
0557 theCablesIL += (dil * theCablesFractionIL);
0558 theCopperIL += (dil * theCopperFractionIL);
0559 theH_ScintillatorIL += (dil * theH_ScintillatorFractionIL);
0560 theLeadIL += (dil * theLeadFractionIL);
0561 theEpoxyIL += (dil * theEpoxyFractionIL);
0562 theKaptonIL += (dil * theKaptonFractionIL);
0563 theAluminiumIL += (dil * theAluminiumFractionIL);
0564 theHGC_G10_FR4IL += (dil * theHGC_G10_FR4FractionIL);
0565 theSiliconIL += (dil * theSiliconFractionIL);
0566 theStainlessSteelIL += (dil * theStainlessSteelFractionIL);
0567 theWCuIL += (dil * theWCuFractionIL);
0568 thePolystyreneIL += (dil * thePolystyreneFractionIL);
0569 theHGC_EEConnectorIL += (dil * theHGC_EEConnectorFractionIL);
0570 theHGC_HEConnectorIL += (dil * theHGC_HEConnectorFractionIL);
0571
0572
0573
0574 theStepN++;
0575 }