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 RE05/src/RE05StackingAction.cc 28 /// \brief Implementation of the RE05StackingA 29 // 30 31 #include "RE05StackingAction.hh" 32 33 #include "RE05StackingActionMessenger.hh" 34 35 #include "G4Event.hh" 36 #include "G4HCofThisEvent.hh" 37 #include "G4ParticleDefinition.hh" 38 #include "G4ParticleTypes.hh" 39 #include "G4RunManager.hh" 40 #include "G4SDManager.hh" 41 #include "G4SystemOfUnits.hh" 42 #include "G4Track.hh" 43 #include "G4TrackStatus.hh" 44 #include "G4ios.hh" 45 46 //....oooOO0OOooo........oooOO0OOooo........oo 47 48 RE05StackingAction::RE05StackingAction() 49 : G4UserStackingAction(), 50 fTrkHits(0), 51 fMuonHits(0), 52 fMessenger(0), 53 fStage(0), 54 fReqMuon(2), 55 fReqIso(10), 56 fAngRoI(30.0 * deg) 57 { 58 fMessenger = new RE05StackingActionMessenger 59 } 60 61 //....oooOO0OOooo........oooOO0OOooo........oo 62 63 RE05StackingAction::~RE05StackingAction() 64 { 65 delete fMessenger; 66 } 67 68 //....oooOO0OOooo........oooOO0OOooo........oo 69 70 G4ClassificationOfNewTrack RE05StackingAction: 71 { 72 G4ClassificationOfNewTrack classification = 73 switch (fStage) { 74 case 0: // Stage 0 : Primary muons only 75 if (aTrack->GetParentID() == 0) { 76 G4ParticleDefinition* particleType = a 77 if ((particleType == G4MuonPlus::MuonP 78 || (particleType == G4MuonMinus::M 79 { 80 classification = fUrgent; 81 } 82 } 83 break; 84 85 case 1: // Stage 1 : Charged primaries on 86 // Suspended tracks wil 87 if (aTrack->GetParentID() != 0) { 88 break; 89 } 90 if (aTrack->GetTrackStatus() == fSuspend 91 break; 92 } 93 if (aTrack->GetDefinition()->GetPDGCharg 94 break; 95 } 96 classification = fUrgent; 97 break; 98 99 default: // Stage 2 : Accept all primarie 100 // Accept all secondar 101 // Kill secondaries ou 102 if (aTrack->GetParentID() == 0) { 103 classification = fUrgent; 104 break; 105 } 106 if ((fAngRoI < 0.) || InsideRoI(aTrack, 107 classification = fUrgent; 108 break; 109 } 110 classification = fKill; 111 } 112 return classification; 113 } 114 115 //....oooOO0OOooo........oooOO0OOooo........oo 116 117 G4bool RE05StackingAction::InsideRoI(const G4T 118 { 119 if (!fMuonHits) { 120 fMuonHits = (RE05MuonHitsCollection*)GetCo 121 } 122 if (!fMuonHits) { 123 G4cerr << "muonCollection NOT FOUND" << G4 124 return true; 125 } 126 127 G4int nhits = fMuonHits->entries(); 128 129 const G4ThreeVector trPos = aTrack->GetPosit 130 for (G4int i = 0; i < nhits; i++) { 131 G4ThreeVector muHitPos = (*fMuonHits)[i]-> 132 G4double angl = muHitPos.angle(trPos); 133 if (angl < ang) { 134 return true; 135 } 136 } 137 138 return false; 139 } 140 141 //....oooOO0OOooo........oooOO0OOooo........oo 142 143 G4VHitsCollection* RE05StackingAction::GetColl 144 { 145 G4SDManager* SDMan = G4SDManager::GetSDMpoin 146 G4RunManager* runMan = G4RunManager::GetRunM 147 int colID = SDMan->GetCollectionID(colName); 148 if (colID >= 0) { 149 const G4Event* currentEvent = runMan->GetC 150 G4HCofThisEvent* HCE = currentEvent->GetHC 151 return HCE->GetHC(colID); 152 } 153 return 0; 154 } 155 156 //....oooOO0OOooo........oooOO0OOooo........oo 157 158 void RE05StackingAction::NewStage() 159 { 160 fStage++; 161 G4int nhits; 162 if (fStage == 1) { 163 // Stage 0->1 : check if at least "fReqMuo 164 // otherwise abort current ev 165 if (!fMuonHits) { 166 fMuonHits = (RE05MuonHitsCollection*)Get 167 } 168 if (!fMuonHits) { 169 G4cerr << "muonCollection NOT FOUND" << 170 return; 171 } 172 nhits = fMuonHits->entries(); 173 //// G4cout << "Stage 0->1 : " << nhits 174 //// << G4endl; 175 if (nhits < fReqMuon) { 176 stackManager->clear(); 177 //// G4cout << "++++++++ event abor 178 return; 179 } 180 stackManager->ReClassify(); 181 return; 182 } 183 184 else if (fStage == 2) { 185 // Stage 1->2 : check the isolation of muo 186 // at least "fReqIsoMuon" iso 187 // otherwise abort current ev 188 // Isolation requires "fReqIs 189 // (including own hits) in th 190 // in the tracker layers. 191 nhits = fMuonHits->entries(); 192 if (!fTrkHits) { 193 fTrkHits = (RE05TrackerHitsCollection*)G 194 } 195 if (!fTrkHits) { 196 G4cerr << "trackerCollection NOT FOUND" 197 return; 198 } 199 G4int nTrkhits = fTrkHits->entries(); 200 G4int isoMuon = 0; 201 for (G4int j = 0; j < nhits; j++) { 202 G4ThreeVector hitPos = (*fMuonHits)[j]-> 203 G4int nhitIn = 0; 204 for (G4int jj = 0; (jj < nTrkhits) && (n 205 G4ThreeVector trkhitPos = (*fTrkHits)[ 206 if (trkhitPos.angle(hitPos) < fAngRoI) 207 } 208 if (nhitIn <= fReqIso) isoMuon++; 209 } 210 //// G4cout << "Stage 1->2 : " << isoMu 211 if (isoMuon < fReqIsoMuon) { 212 stackManager->clear(); 213 //// G4cout << "++++++++ event abor 214 return; 215 } 216 stackManager->ReClassify(); 217 return; 218 } 219 220 else { 221 // Other fStage change : just re-classify 222 stackManager->ReClassify(); 223 } 224 } 225 226 //....oooOO0OOooo........oooOO0OOooo........oo 227 228 void RE05StackingAction::PrepareNewEvent() 229 { 230 fStage = 0; 231 fTrkHits = 0; 232 fMuonHits = 0; 233 } 234 235 //....oooOO0OOooo........oooOO0OOooo........oo 236