Geant4 Cross Reference |
1 // 1 2 // ******************************************* 3 // * License and Disclaimer 4 // * 5 // * The Geant4 software is copyright of th 6 // * the Geant4 Collaboration. It is provided 7 // * conditions of the Geant4 Software License 8 // * LICENSE and available at http://cern.ch/ 9 // * include a list of copyright holders. 10 // * 11 // * Neither the authors of this software syst 12 // * institutes,nor the agencies providing fin 13 // * work make any representation or warran 14 // * regarding this software system or assum 15 // * use. Please see the license in the file 16 // * for the full disclaimer and the limitatio 17 // * 18 // * This code implementation is the result 19 // * technical work of the GEANT4 collaboratio 20 // * By using, copying, modifying or distri 21 // * any work based on the software) you ag 22 // * use in resulting scientific publicati 23 // * acceptance of all terms of the Geant4 Sof 24 // ******************************************* 25 // 26 // 27 /// \file XrayTESdetSteppingAction.cc 28 /// \brief Implementation of the SteppingActio 29 // 30 // Authors: P.Dondero (paolo.dondero@cern.ch), 31 // 32 //....oooOO0OOooo........oooOO0OOooo........oo 33 //....oooOO0OOooo........oooOO0OOooo........oo 34 35 36 #include "XrayTESdetSteppingAction.hh" 37 #include "G4AnalysisManager.hh" 38 39 #include "G4SteppingManager.hh" 40 #include "G4RunManager.hh" 41 #include "G4SystemOfUnits.hh" 42 #include "G4PhysicalConstants.hh" 43 44 //....oooOO0OOooo........oooOO0OOooo........oo 45 46 G4String clean_name (const G4String& vol_name) 47 { 48 G4String cleaned_name = ""; 49 50 // Search if the vol_name contains one of th 51 if (((std::string)vol_name).find("Bipxl") != 52 cleaned_name = "Bipxl"; 53 if (((std::string)vol_name).find("membrane") 54 cleaned_name = "membrane"; 55 if (((std::string)vol_name).find("gridpiece" 56 cleaned_name = "gridpiece"; 57 if (((std::string)vol_name).find("trapezoid" 58 cleaned_name = "Mesh"; 59 if (cleaned_name == "") cleaned_name = vol_n 60 return cleaned_name; 61 } 62 63 //....oooOO0OOooo........oooOO0OOooo........oo 64 65 void XrayTESdetSteppingAction::UserSteppingAct 66 { 67 // Define quantities to use in the 68 G4RunManager *rm = G4RunManager::GetRunManag 69 auto analysisManager = G4AnalysisManager::In 70 71 G4int eventID = rm->GetCurrentEvent()->GetEv 72 G4Track* track = step->GetTrack(); 73 G4StepPoint* pre_step_point = step->GetPreSt 74 G4int parentID = track->GetParentID(); 75 G4int trackID = track->GetTrackID(); 76 const G4String& vol_name = track->GetVolume( 77 G4String mother_name = ""; 78 G4double init_kinetic_energy = 0; 79 G4double kinetic_energy = track->GetKineticE 80 G4double pre_kinetic_energy = pre_step_point 81 const G4String& particle_name = track->GetPa 82 const G4VProcess* pre_step_proc = pre_step_p 83 const G4VProcess* post_step_proc = pre_step_ 84 85 bool proceed = false; 86 G4String pre_step_name_s = "Undefined"; 87 G4String post_step_name_s = "Undefined"; 88 G4String creator_process_name = "Undefined"; 89 G4String init_creator_process = "Undefined"; 90 const G4String& pre_step_name = pre_step_nam 91 const G4String& post_step_name = post_step_n 92 93 if (track->GetCreatorProcess() != nullptr) 94 { 95 creator_process_name = track->GetCreatorPr 96 } 97 98 if (pre_step_proc != nullptr) 99 { 100 pre_step_name_s = pre_step_proc->GetProces 101 } 102 103 if (post_step_proc != nullptr) 104 { 105 post_step_name_s = post_step_proc->GetProc 106 } 107 108 // Saves the first energy of the secondaries 109 if (eventID != fPrev_eventID) 110 { 111 fInit_energy.clear(); 112 fCreator_proc.clear(); 113 } 114 115 G4int current_step_number = -1; 116 current_step_number = track->GetCurrentStepN 117 if (current_step_number != -1) 118 { 119 if (current_step_number == 1 && parentID = 120 { 121 proceed = true; 122 fInit_energy.insert(std::pair<G4int, G4d 123 fCreator_proc.insert(std::pair<G4int, G4 124 pre_step_name_s = "InitStep"; 125 } 126 } 127 128 // Fill the map for the secondaries 129 if (trackID > 1) 130 { 131 // First value for each trackIDs; this ser 132 if (fInit_energy.count((int)trackID) == 0) 133 { 134 fInit_energy.insert(std::pair<G4int, G4d 135 fCreator_proc.insert(std::pair<G4int, G4 136 } 137 } 138 139 // Exclude the world volume, which has no mo 140 if (((std::string)vol_name).find("expHall") 141 { 142 if (track->GetVolume() != nullptr) 143 { 144 if (track->GetVolume()->GetMotherLogical 145 { 146 mother_name = track->GetVolume()->GetM 147 } 148 } 149 if (((std::string)mother_name).find("Detec 150 { 151 proceed = true; 152 } 153 } 154 155 init_kinetic_energy = fInit_energy[(int)trac 156 init_creator_process = fCreator_proc[trackID 157 158 G4String cleaned_name = ""; 159 cleaned_name = clean_name(vol_name); 160 161 // Saves on tuple 162 if (proceed) 163 { 164 G4double x = pre_step_point->GetPosition() 165 G4double y = pre_step_point->GetPosition() 166 G4double z = pre_step_point->GetPosition() 167 G4ThreeVector direction = pre_step_point-> 168 G4double theta = direction.theta(); 169 G4double phi = direction.phi(); 170 G4int pixel_number = track->GetVolume()->G 171 G4int def_replica_number = -999; 172 G4int step_number = track->GetCurrentStepN 173 G4double step_energy_dep = step->GetTotalE 174 175 // Fill tuple 176 analysisManager->FillNtupleIColumn(0, even 177 analysisManager->FillNtupleSColumn(1, clea 178 analysisManager->FillNtupleIColumn(2, trac 179 analysisManager->FillNtupleDColumn(3, x); 180 analysisManager->FillNtupleDColumn(4, y); 181 analysisManager->FillNtupleDColumn(5, z); 182 analysisManager->FillNtupleDColumn(6, thet 183 analysisManager->FillNtupleDColumn(7, phi) 184 analysisManager->FillNtupleIColumn(8, pare 185 if (((std::string)vol_name).find("Bipxl") 186 { 187 analysisManager->FillNtupleIColumn(9, pi 188 } 189 else 190 { 191 analysisManager->FillNtupleIColumn(9, de 192 } 193 analysisManager->FillNtupleDColumn(10, ste 194 analysisManager->FillNtupleIColumn(11, ste 195 analysisManager->FillNtupleDColumn(12, ini 196 analysisManager->FillNtupleDColumn(13, kin 197 analysisManager->FillNtupleSColumn(14, par 198 analysisManager->FillNtupleSColumn(15, pre 199 analysisManager->FillNtupleSColumn(16, pos 200 analysisManager->FillNtupleSColumn(17, ini 201 analysisManager->AddNtupleRow(0); 202 } 203 fPrev_eventID = eventID; 204 fPrev_trackID = trackID; 205 } 206 207 //....oooOO0OOooo........oooOO0OOooo........oo 208