Geant4 Cross Reference |
1 // 1 // 2 // ******************************************* 2 // ******************************************************************** 3 // * License and Disclaimer 3 // * License and Disclaimer * 4 // * 4 // * * 5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of * 6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and * 7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file * 8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These * 9 // * include a list of copyright holders. 9 // * include a list of copyright holders. * 10 // * 10 // * * 11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing * 12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this * 13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, * 14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its * 15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above * 16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. * 17 // * 17 // * * 18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and * 19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. * 20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or * 21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its * 22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your * 23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. * 24 // ******************************************* 24 // ******************************************************************** 25 // 25 // 26 // 26 // >> 27 // $Id$ >> 28 // GEANT4 tag $Name: geant4-09-04 $ 27 // 29 // 28 // G4PSPassageCellFlux 30 // G4PSPassageCellFlux 29 #include "G4PSPassageCellFlux.hh" 31 #include "G4PSPassageCellFlux.hh" 30 32 31 #include "G4SystemOfUnits.hh" 33 #include "G4SystemOfUnits.hh" 32 #include "G4StepStatus.hh" 34 #include "G4StepStatus.hh" 33 #include "G4Track.hh" 35 #include "G4Track.hh" 34 #include "G4VSolid.hh" 36 #include "G4VSolid.hh" 35 #include "G4VPhysicalVolume.hh" 37 #include "G4VPhysicalVolume.hh" 36 #include "G4VPVParameterisation.hh" 38 #include "G4VPVParameterisation.hh" 37 #include "G4UnitsTable.hh" 39 #include "G4UnitsTable.hh" 38 #include "G4VScoreHistFiller.hh" << 39 << 40 ////////////////////////////////////////////// 40 //////////////////////////////////////////////////////////////////////////////// 41 // (Description) 41 // (Description) 42 // This is a primitive scorer class for scor 42 // This is a primitive scorer class for scoring cell flux. 43 // The Cell Flux is defined by a track leng 43 // The Cell Flux is defined by a track length divided by a geometry 44 // volume, where only tracks passing through << 44 // volume, where only tracks passing through the geometry are taken 45 // into account. e.g. the unit of Cell Flux i 45 // into account. e.g. the unit of Cell Flux is mm/mm3. 46 // 46 // 47 // If you want to score all tracks in the ge 47 // If you want to score all tracks in the geometry volume, 48 // please use G4PSCellFlux. 48 // please use G4PSCellFlux. 49 // 49 // 50 // Created: 2005-11-14 Tsukasa ASO, Akinori K 50 // Created: 2005-11-14 Tsukasa ASO, Akinori Kimura. 51 // 2010-07-22 Introduce Unit specification. 51 // 2010-07-22 Introduce Unit specification. 52 // 2010-07-22 Add weighted option 52 // 2010-07-22 Add weighted option 53 // 2020-10-06 Use G4VPrimitivePlotter and fi << 53 // 54 // vs. cell flux * track weight ( << 55 // << 56 ////////////////////////////////////////////// 54 /////////////////////////////////////////////////////////////////////////////// 57 55 58 G4PSPassageCellFlux::G4PSPassageCellFlux(const << 56 G4PSPassageCellFlux::G4PSPassageCellFlux(G4String name, G4int depth) 59 : G4PSPassageCellFlux(name, "percm2", depth) << 57 : G4VPrimitiveScorer(name,depth),HCID(-1),fCurrentTrkID(-1),fCellFlux(0), 60 {} << 58 weighted(true) 61 << 59 { 62 G4PSPassageCellFlux::G4PSPassageCellFlux(const << 60 DefineUnitAndCategory(); 63 G4int << 61 SetUnit("percm2"); 64 : G4VPrimitivePlotter(name, depth) << 62 } 65 , HCID(-1) << 63 66 , fCurrentTrkID(-1) << 64 G4PSPassageCellFlux::G4PSPassageCellFlux(G4String name, const G4String& unit, 67 , fCellFlux(0) << 65 G4int depth) 68 , EvtMap(nullptr) << 66 : G4VPrimitiveScorer(name,depth),HCID(-1),fCurrentTrkID(-1),fCellFlux(0), 69 , weighted(true) << 67 weighted(true) 70 { << 68 { 71 DefineUnitAndCategory(); << 69 DefineUnitAndCategory(); 72 SetUnit(unit); << 70 SetUnit(unit); 73 } << 71 } 74 << 72 75 G4bool G4PSPassageCellFlux::ProcessHits(G4Step << 73 G4PSPassageCellFlux::~G4PSPassageCellFlux() 76 { << 74 {;} 77 if(IsPassed(aStep)) << 75 78 { << 76 G4bool G4PSPassageCellFlux::ProcessHits(G4Step* aStep,G4TouchableHistory*) 79 G4int idx = << 77 { 80 ((G4TouchableHistory*) (aStep->GetPreSte << 78 81 ->GetReplicaNumber(indexDepth); << 79 if ( IsPassed(aStep) ) { >> 80 G4int idx = ((G4TouchableHistory*) >> 81 (aStep->GetPreStepPoint()->GetTouchable())) >> 82 ->GetReplicaNumber(indexDepth); 82 G4double cubicVolume = ComputeVolume(aStep 83 G4double cubicVolume = ComputeVolume(aStep, idx); 83 84 84 fCellFlux /= cubicVolume; 85 fCellFlux /= cubicVolume; 85 G4int index = GetIndex(aStep); 86 G4int index = GetIndex(aStep); 86 EvtMap->add(index, fCellFlux); << 87 EvtMap->add(index,fCellFlux); 87 << 88 if(!hitIDMap.empty() && hitIDMap.find(inde << 89 { << 90 auto filler = G4VScoreHistFiller::Instan << 91 if(filler == nullptr) << 92 { << 93 G4Exception( << 94 "G4PSPassageCellFlux::ProcessHits", << 95 "G4TScoreHistFiller is not instantia << 96 } << 97 else << 98 { << 99 filler->FillH1(hitIDMap[index], << 100 aStep->GetPreStepPoint( << 101 } << 102 } << 103 } 88 } 104 89 105 return true; << 90 return TRUE; 106 } 91 } 107 92 108 G4bool G4PSPassageCellFlux::IsPassed(G4Step* a << 93 G4bool G4PSPassageCellFlux::IsPassed(G4Step* aStep){ 109 { << 94 G4bool Passed = FALSE; 110 G4bool Passed = false; << 111 95 112 G4bool IsEnter = aStep->GetPreStepPoint()->G 96 G4bool IsEnter = aStep->GetPreStepPoint()->GetStepStatus() == fGeomBoundary; 113 G4bool IsExit = aStep->GetPostStepPoint()-> 97 G4bool IsExit = aStep->GetPostStepPoint()->GetStepStatus() == fGeomBoundary; 114 98 115 G4int trkid = aStep->GetTrack()->GetT << 99 G4int trkid = aStep->GetTrack()->GetTrackID(); 116 G4double trklength = aStep->GetStepLength(); << 100 G4double trklength = aStep->GetStepLength(); 117 if(weighted) << 101 if ( weighted ) trklength *= aStep->GetPreStepPoint()->GetWeight(); 118 trklength *= aStep->GetPreStepPoint()->Get << 102 119 << 103 if ( IsEnter &&IsExit ){ // Passed at one step 120 if(IsEnter && IsExit) << 104 fCellFlux = trklength; // Track length is absolutely given. 121 { // Passed at one s << 105 Passed = TRUE; 122 fCellFlux = trklength; // Track length is << 106 }else if ( IsEnter ){ // Enter a new geometry 123 Passed = true; << 107 fCurrentTrkID = trkid; // Resetting the current track. 124 } << 108 fCellFlux = trklength; 125 else if(IsEnter) << 109 }else if ( IsExit ){ // Exit a current geometry 126 { // Enter a new geo << 110 if ( fCurrentTrkID == trkid ) {// if the track is same as entered, 127 fCurrentTrkID = trkid; // Resetting the c << 111 fCellFlux += trklength; // add the track length to current one. 128 fCellFlux = trklength; << 112 Passed = TRUE; 129 } << 130 else if(IsExit) << 131 { // Exit a current geometry << 132 if(fCurrentTrkID == trkid) << 133 { // if the track << 134 fCellFlux += trklength; // add the trac << 135 Passed = true; << 136 } 113 } 137 } << 114 }else{ // Inside geometry 138 else << 115 if ( fCurrentTrkID == trkid ){ // if the track is same as entered, 139 { // Inside geometry << 116 fCellFlux += trklength; // adding the track length to current one. 140 if(fCurrentTrkID == trkid) << 141 { // if the track << 142 fCellFlux += trklength; // adding the t << 143 } 117 } 144 } 118 } 145 119 146 return Passed; 120 return Passed; 147 } 121 } 148 122 149 void G4PSPassageCellFlux::Initialize(G4HCofThi 123 void G4PSPassageCellFlux::Initialize(G4HCofThisEvent* HCE) 150 { 124 { 151 fCurrentTrkID = -1; 125 fCurrentTrkID = -1; 152 126 153 EvtMap = new G4THitsMap<G4double>(detector-> << 127 EvtMap = new G4THitsMap<G4double>(detector->GetName(), 154 if(HCID < 0) << 128 GetName()); 155 HCID = GetCollectionID(0); << 129 if ( HCID < 0 ) HCID = GetCollectionID(0); 156 HCE->AddHitsCollection(HCID, EvtMap); << 130 HCE->AddHitsCollection(HCID,EvtMap); >> 131 157 } 132 } 158 133 159 void G4PSPassageCellFlux::clear() { EvtMap->cl << 134 void G4PSPassageCellFlux::EndOfEvent(G4HCofThisEvent*) >> 135 {;} >> 136 >> 137 void G4PSPassageCellFlux::clear(){ >> 138 EvtMap->clear(); >> 139 } >> 140 >> 141 void G4PSPassageCellFlux::DrawAll() >> 142 {;} 160 143 161 void G4PSPassageCellFlux::PrintAll() 144 void G4PSPassageCellFlux::PrintAll() 162 { 145 { 163 G4cout << " MultiFunctionalDet " << detecto 146 G4cout << " MultiFunctionalDet " << detector->GetName() << G4endl; 164 G4cout << " PrimitiveScorer " << GetName() < << 147 G4cout << " PrimitiveScorer " << GetName() <<G4endl; 165 G4cout << " Number of entries " << EvtMap->e 148 G4cout << " Number of entries " << EvtMap->entries() << G4endl; 166 for(const auto& [copy, flux] : *(EvtMap->Get << 149 std::map<G4int,G4double*>::iterator itr = EvtMap->GetMap()->begin(); 167 { << 150 for(; itr != EvtMap->GetMap()->end(); itr++) { 168 G4cout << " copy no.: " << copy << 151 G4cout << " copy no.: " << itr->first 169 << " cell flux : " << *(flux) / Ge << 152 << " cell flux : " << *(itr->second)/GetUnitValue() 170 << GetUnit() << G4endl; << 153 << " [" << GetUnit() >> 154 << G4endl; 171 } 155 } 172 } 156 } 173 157 174 void G4PSPassageCellFlux::SetUnit(const G4Stri 158 void G4PSPassageCellFlux::SetUnit(const G4String& unit) 175 { 159 { 176 CheckAndSetUnit(unit, "Per Unit Surface"); << 160 CheckAndSetUnit(unit,"Per Unit Surface"); 177 } 161 } 178 162 179 void G4PSPassageCellFlux::DefineUnitAndCategor << 163 void G4PSPassageCellFlux::DefineUnitAndCategory(){ 180 { << 164 // Per Unit Surface 181 // Per Unit Surface << 165 new G4UnitDefinition("percentimeter2","percm2","Per Unit Surface",(1./cm2)); 182 new G4UnitDefinition("percentimeter2", "perc << 166 new G4UnitDefinition("permillimeter2","permm2","Per Unit Surface",(1./mm2)); 183 (1. / cm2)); << 167 new G4UnitDefinition("permeter2","perm2","Per Unit Surface",(1./m2)); 184 new G4UnitDefinition("permillimeter2", "perm << 185 (1. / mm2)); << 186 new G4UnitDefinition("permeter2", "perm2", " << 187 } 168 } 188 169 189 G4double G4PSPassageCellFlux::ComputeVolume(G4 << 170 190 { << 171 G4double G4PSPassageCellFlux::ComputeVolume(G4Step* aStep, G4int idx){ 191 return ComputeSolid(aStep, idx)->GetCubicVol << 172 >> 173 G4VPhysicalVolume* physVol = aStep->GetPreStepPoint()->GetPhysicalVolume(); >> 174 G4VPVParameterisation* physParam = physVol->GetParameterisation(); >> 175 G4VSolid* solid = 0; >> 176 if(physParam) >> 177 { // for parameterized volume >> 178 if(idx<0) >> 179 { >> 180 G4ExceptionDescription ED; >> 181 ED << "Incorrect replica number --- GetReplicaNumber : " << idx << G4endl; >> 182 G4Exception("G4PSPassageCellFlux::ComputeVolume","DetPS0013",JustWarning,ED); >> 183 } >> 184 solid = physParam->ComputeSolid(idx, physVol); >> 185 solid->ComputeDimensions(physParam,idx,physVol); >> 186 } >> 187 else >> 188 { // for ordinary volume >> 189 solid = physVol->GetLogicalVolume()->GetSolid(); >> 190 } >> 191 >> 192 return solid->GetCubicVolume(); 192 } 193 } 193 194