#include <G4PSFlatSurfaceFlux.hh>
Inheritance diagram for G4PSFlatSurfaceFlux:
Public Member Functions | |
G4PSFlatSurfaceFlux (G4String name, G4int direction, G4int depth=0) | |
G4PSFlatSurfaceFlux (G4String name, G4int direction, const G4String &unit, G4int depth=0) | |
virtual | ~G4PSFlatSurfaceFlux () |
void | Weighted (G4bool flg=true) |
void | DivideByArea (G4bool flg=true) |
virtual void | Initialize (G4HCofThisEvent *) |
virtual void | EndOfEvent (G4HCofThisEvent *) |
virtual void | clear () |
virtual void | DrawAll () |
virtual void | PrintAll () |
virtual void | SetUnit (const G4String &unit) |
Protected Member Functions | |
virtual G4bool | ProcessHits (G4Step *, G4TouchableHistory *) |
G4int | IsSelectedSurface (G4Step *, G4Box *) |
virtual void | DefineUnitAndCategory () |
Definition at line 62 of file G4PSFlatSurfaceFlux.hh.
Definition at line 60 of file G4PSFlatSurfaceFlux.cc.
References DefineUnitAndCategory(), and SetUnit().
00062 : G4VPrimitiveScorer(name,depth),HCID(-1),fDirection(direction), 00063 weighted(true),divideByArea(true) 00064 { 00065 DefineUnitAndCategory(); 00066 SetUnit("percm2"); 00067 }
G4PSFlatSurfaceFlux::G4PSFlatSurfaceFlux | ( | G4String | name, | |
G4int | direction, | |||
const G4String & | unit, | |||
G4int | depth = 0 | |||
) |
Definition at line 69 of file G4PSFlatSurfaceFlux.cc.
References DefineUnitAndCategory(), and SetUnit().
00073 : G4VPrimitiveScorer(name,depth),HCID(-1),fDirection(direction), 00074 weighted(true),divideByArea(true) 00075 { 00076 DefineUnitAndCategory(); 00077 SetUnit(unit); 00078 }
G4PSFlatSurfaceFlux::~G4PSFlatSurfaceFlux | ( | ) | [virtual] |
void G4PSFlatSurfaceFlux::clear | ( | ) | [virtual] |
void G4PSFlatSurfaceFlux::DefineUnitAndCategory | ( | ) | [protected, virtual] |
Definition at line 219 of file G4PSFlatSurfaceFlux.cc.
Referenced by G4PSFlatSurfaceFlux().
00219 { 00220 // Per Unit Surface 00221 new G4UnitDefinition("percentimeter2","percm2","Per Unit Surface",(1./cm2)); 00222 new G4UnitDefinition("permillimeter2","permm2","Per Unit Surface",(1./mm2)); 00223 new G4UnitDefinition("permeter2","perm2","Per Unit Surface",(1./m2)); 00224 }
void G4PSFlatSurfaceFlux::DivideByArea | ( | G4bool | flg = true |
) | [inline] |
Definition at line 73 of file G4PSFlatSurfaceFlux.hh.
Referenced by G4ScoreQuantityMessenger::SetNewValue().
void G4PSFlatSurfaceFlux::DrawAll | ( | ) | [virtual] |
void G4PSFlatSurfaceFlux::EndOfEvent | ( | G4HCofThisEvent * | ) | [virtual] |
void G4PSFlatSurfaceFlux::Initialize | ( | G4HCofThisEvent * | ) | [virtual] |
Reimplemented from G4VPrimitiveScorer.
Definition at line 172 of file G4PSFlatSurfaceFlux.cc.
References G4HCofThisEvent::AddHitsCollection(), G4VPrimitiveScorer::GetCollectionID(), G4VPrimitiveScorer::GetMultiFunctionalDetector(), G4VPrimitiveScorer::GetName(), and G4VSensitiveDetector::GetName().
00173 { 00174 EvtMap = new G4THitsMap<G4double>(GetMultiFunctionalDetector()->GetName(), 00175 GetName()); 00176 if ( HCID < 0 ) HCID = GetCollectionID(0); 00177 HCE->AddHitsCollection(HCID, (G4VHitsCollection*)EvtMap); 00178 }
Definition at line 143 of file G4PSFlatSurfaceFlux.cc.
References fFlux_In, fFlux_Out, fGeomBoundary, G4GeometryTolerance::GetInstance(), G4StepPoint::GetPosition(), G4Step::GetPostStepPoint(), G4Step::GetPreStepPoint(), G4StepPoint::GetStepStatus(), G4GeometryTolerance::GetSurfaceTolerance(), G4StepPoint::GetTouchableHandle(), and G4Box::GetZHalfLength().
Referenced by ProcessHits().
00143 { 00144 00145 G4TouchableHandle theTouchable = 00146 aStep->GetPreStepPoint()->GetTouchableHandle(); 00147 G4double kCarTolerance=G4GeometryTolerance::GetInstance()->GetSurfaceTolerance(); 00148 00149 if (aStep->GetPreStepPoint()->GetStepStatus() == fGeomBoundary ){ 00150 // Entering Geometry 00151 G4ThreeVector stppos1= aStep->GetPreStepPoint()->GetPosition(); 00152 G4ThreeVector localpos1 = 00153 theTouchable->GetHistory()->GetTopTransform().TransformPoint(stppos1); 00154 if(std::fabs( localpos1.z() + boxSolid->GetZHalfLength())<kCarTolerance ){ 00155 return fFlux_In; 00156 } 00157 } 00158 00159 if (aStep->GetPostStepPoint()->GetStepStatus() == fGeomBoundary ){ 00160 // Exiting Geometry 00161 G4ThreeVector stppos2= aStep->GetPostStepPoint()->GetPosition(); 00162 G4ThreeVector localpos2 = 00163 theTouchable->GetHistory()->GetTopTransform().TransformPoint(stppos2); 00164 if(std::fabs( localpos2.z() + boxSolid->GetZHalfLength())<kCarTolerance ){ 00165 return fFlux_Out; 00166 } 00167 } 00168 00169 return -1; 00170 }
void G4PSFlatSurfaceFlux::PrintAll | ( | ) | [virtual] |
Reimplemented from G4VPrimitiveScorer.
Definition at line 190 of file G4PSFlatSurfaceFlux.cc.
References G4VPrimitiveScorer::detector, G4cout, G4endl, G4VPrimitiveScorer::GetName(), G4VSensitiveDetector::GetName(), G4VPrimitiveScorer::GetUnit(), and G4VPrimitiveScorer::GetUnitValue().
00191 { 00192 G4cout << " MultiFunctionalDet " << detector->GetName() << G4endl; 00193 G4cout << " PrimitiveScorer" << GetName() <<G4endl; 00194 G4cout << " Number of entries " << EvtMap->entries() << G4endl; 00195 std::map<G4int,G4double*>::iterator itr = EvtMap->GetMap()->begin(); 00196 for(; itr != EvtMap->GetMap()->end(); itr++) { 00197 G4cout << " copy no.: " << itr->first 00198 << " flux : " << *(itr->second)/GetUnitValue() 00199 << " [" << GetUnit() <<"]" 00200 << G4endl; 00201 } 00202 }
G4bool G4PSFlatSurfaceFlux::ProcessHits | ( | G4Step * | , | |
G4TouchableHistory * | ||||
) | [protected, virtual] |
Implements G4VPrimitiveScorer.
Definition at line 83 of file G4PSFlatSurfaceFlux.cc.
References G4VSolid::ComputeDimensions(), G4VPVParameterisation::ComputeSolid(), FALSE, fFlux_In, fFlux_InOut, fFlux_Out, G4cout, G4endl, G4VPrimitiveScorer::GetIndex(), G4VPhysicalVolume::GetLogicalVolume(), G4StepPoint::GetMomentumDirection(), G4VPhysicalVolume::GetParameterisation(), G4StepPoint::GetPhysicalVolume(), G4Step::GetPostStepPoint(), G4Step::GetPreStepPoint(), G4LogicalVolume::GetSolid(), G4StepPoint::GetTouchable(), G4StepPoint::GetTouchableHandle(), G4StepPoint::GetWeight(), G4Box::GetXHalfLength(), G4Box::GetYHalfLength(), G4VPrimitiveScorer::indexDepth, IsSelectedSurface(), and TRUE.
00084 { 00085 G4StepPoint* preStep = aStep->GetPreStepPoint(); 00086 G4VPhysicalVolume* physVol = preStep->GetPhysicalVolume(); 00087 G4VPVParameterisation* physParam = physVol->GetParameterisation(); 00088 G4VSolid * solid = 0; 00089 if(physParam) 00090 { // for parameterized volume 00091 G4int idx = ((G4TouchableHistory*)(aStep->GetPreStepPoint()->GetTouchable())) 00092 ->GetReplicaNumber(indexDepth); 00093 solid = physParam->ComputeSolid(idx, physVol); 00094 solid->ComputeDimensions(physParam,idx,physVol); 00095 } 00096 else 00097 { // for ordinary volume 00098 solid = physVol->GetLogicalVolume()->GetSolid(); 00099 } 00100 00101 G4Box* boxSolid = (G4Box*)(solid); 00102 00103 G4int dirFlag =IsSelectedSurface(aStep,boxSolid); 00104 if ( dirFlag > 0 ) { 00105 if ( fDirection == fFlux_InOut || fDirection == dirFlag ){ 00106 00107 G4StepPoint* thisStep=0; 00108 if ( dirFlag == fFlux_In ){ 00109 thisStep = preStep; 00110 }else if ( dirFlag == fFlux_Out ){ 00111 thisStep = aStep->GetPostStepPoint(); 00112 }else{ 00113 return FALSE; 00114 } 00115 00116 G4TouchableHandle theTouchable = thisStep->GetTouchableHandle(); 00117 G4ThreeVector pdirection = thisStep->GetMomentumDirection(); 00118 G4ThreeVector localdir = 00119 theTouchable->GetHistory()->GetTopTransform().TransformAxis(pdirection); 00120 // 00121 G4double angleFactor = localdir.z(); 00122 if ( angleFactor < 0 ) angleFactor *= -1.; 00123 G4double flux = 1.0; 00124 if ( weighted ) flux *=preStep->GetWeight(); // Current (Particle Weight) 00125 // 00126 G4double square = 4.*boxSolid->GetXHalfLength()*boxSolid->GetYHalfLength(); 00127 // 00128 flux = flux/angleFactor; // Flux with angle. 00129 if ( divideByArea ) flux /= square; 00130 // 00131 G4int index = GetIndex(aStep); 00132 EvtMap->add(index,flux); 00133 } 00134 } 00135 #ifdef debug 00136 G4cout << " PASSED vol " 00137 << index << " trk "<<trkid<<" len " << fFlatSurfaceFlux<<G4endl; 00138 #endif 00139 00140 return TRUE; 00141 }
void G4PSFlatSurfaceFlux::SetUnit | ( | const G4String & | unit | ) | [virtual] |
Reimplemented from G4VPrimitiveScorer.
Definition at line 204 of file G4PSFlatSurfaceFlux.cc.
References G4VPrimitiveScorer::CheckAndSetUnit(), G4Exception(), G4VPrimitiveScorer::GetName(), G4VPrimitiveScorer::GetUnit(), JustWarning, G4VPrimitiveScorer::unitName, and G4VPrimitiveScorer::unitValue.
Referenced by G4PSFlatSurfaceFlux(), G4PSFlatSurfaceFlux3D::G4PSFlatSurfaceFlux3D(), and G4ScoreQuantityMessenger::SetNewValue().
00205 { 00206 if ( divideByArea ) { 00207 CheckAndSetUnit(unit,"Per Unit Surface"); 00208 } else { 00209 if (unit == "" ){ 00210 unitName = unit; 00211 unitValue = 1.0; 00212 }else{ 00213 G4String msg = "Invalid unit ["+unit+"] (Current unit is [" +GetUnit()+"] ) for " + GetName(); 00214 G4Exception("G4PSFlatSurfaceFlux::SetUnit","DetPS0008",JustWarning,msg); 00215 } 00216 } 00217 }
void G4PSFlatSurfaceFlux::Weighted | ( | G4bool | flg = true |
) | [inline] |
Definition at line 70 of file G4PSFlatSurfaceFlux.hh.
Referenced by G4ScoreQuantityMessenger::SetNewValue().