#include <G4PSCylinderSurfaceFlux.hh>
Inheritance diagram for G4PSCylinderSurfaceFlux:
Public Member Functions | |
G4PSCylinderSurfaceFlux (G4String name, G4int direction, G4int depth=0) | |
G4PSCylinderSurfaceFlux (G4String name, G4int direction, const G4String &unit, G4int depth=0) | |
virtual | ~G4PSCylinderSurfaceFlux () |
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 *, G4Tubs *) |
virtual void | DefineUnitAndCategory () |
Definition at line 59 of file G4PSCylinderSurfaceFlux.hh.
G4PSCylinderSurfaceFlux::G4PSCylinderSurfaceFlux | ( | G4String | name, | |
G4int | direction, | |||
G4int | depth = 0 | |||
) |
Definition at line 58 of file G4PSCylinderSurfaceFlux.cc.
References DefineUnitAndCategory(), and SetUnit().
00060 : G4VPrimitiveScorer(name,depth),HCID(-1),fDirection(direction), 00061 weighted(true),divideByArea(true) 00062 { 00063 DefineUnitAndCategory(); 00064 SetUnit("percm2"); 00065 }
G4PSCylinderSurfaceFlux::G4PSCylinderSurfaceFlux | ( | G4String | name, | |
G4int | direction, | |||
const G4String & | unit, | |||
G4int | depth = 0 | |||
) |
Definition at line 67 of file G4PSCylinderSurfaceFlux.cc.
References DefineUnitAndCategory(), and SetUnit().
00071 : G4VPrimitiveScorer(name,depth),HCID(-1),fDirection(direction), 00072 weighted(true),divideByArea(true) 00073 { 00074 DefineUnitAndCategory(); 00075 SetUnit(unit); 00076 }
G4PSCylinderSurfaceFlux::~G4PSCylinderSurfaceFlux | ( | ) | [virtual] |
void G4PSCylinderSurfaceFlux::clear | ( | ) | [virtual] |
void G4PSCylinderSurfaceFlux::DefineUnitAndCategory | ( | ) | [protected, virtual] |
Definition at line 240 of file G4PSCylinderSurfaceFlux.cc.
Referenced by G4PSCylinderSurfaceFlux().
00240 { 00241 // Per Unit Surface 00242 new G4UnitDefinition("percentimeter2","percm2","Per Unit Surface",(1./cm2)); 00243 new G4UnitDefinition("permillimeter2","permm2","Per Unit Surface",(1./mm2)); 00244 new G4UnitDefinition("permeter2","perm2","Per Unit Surface",(1./m2)); 00245 }
void G4PSCylinderSurfaceFlux::DivideByArea | ( | G4bool | flg = true |
) | [inline] |
void G4PSCylinderSurfaceFlux::DrawAll | ( | ) | [virtual] |
void G4PSCylinderSurfaceFlux::EndOfEvent | ( | G4HCofThisEvent * | ) | [virtual] |
void G4PSCylinderSurfaceFlux::Initialize | ( | G4HCofThisEvent * | ) | [virtual] |
Reimplemented from G4VPrimitiveScorer.
Definition at line 193 of file G4PSCylinderSurfaceFlux.cc.
References G4HCofThisEvent::AddHitsCollection(), G4VPrimitiveScorer::GetCollectionID(), G4VPrimitiveScorer::GetMultiFunctionalDetector(), G4VPrimitiveScorer::GetName(), and G4VSensitiveDetector::GetName().
00194 { 00195 EvtMap = new G4THitsMap<G4double>(GetMultiFunctionalDetector()->GetName(), 00196 GetName()); 00197 if ( HCID < 0 ) HCID = GetCollectionID(0); 00198 HCE->AddHitsCollection(HCID, (G4VHitsCollection*)EvtMap); 00199 }
Definition at line 150 of file G4PSCylinderSurfaceFlux.cc.
References fFlux_In, fFlux_Out, fGeomBoundary, G4Tubs::GetInnerRadius(), G4GeometryTolerance::GetInstance(), G4StepPoint::GetPosition(), G4Step::GetPostStepPoint(), G4Step::GetPreStepPoint(), G4StepPoint::GetStepStatus(), G4GeometryTolerance::GetSurfaceTolerance(), G4StepPoint::GetTouchableHandle(), and G4Tubs::GetZHalfLength().
Referenced by ProcessHits().
00150 { 00151 00152 G4TouchableHandle theTouchable = 00153 aStep->GetPreStepPoint()->GetTouchableHandle(); 00154 G4double kCarTolerance = G4GeometryTolerance::GetInstance()->GetSurfaceTolerance(); 00155 00156 if (aStep->GetPreStepPoint()->GetStepStatus() == fGeomBoundary ){ 00157 // Entering Geometry 00158 G4ThreeVector stppos1= aStep->GetPreStepPoint()->GetPosition(); 00159 G4ThreeVector localpos1 = 00160 theTouchable->GetHistory()->GetTopTransform().TransformPoint(stppos1); 00161 if ( std::fabs(localpos1.z()) > tubsSolid->GetZHalfLength() ) return -1; 00162 //if(std::fabs( localpos1.x()*localpos1.x()+localpos1.y()*localpos1.y() 00163 // - (tubsSolid->GetInnerRadius()*tubsSolid->GetInnerRadius())) 00164 // <kCarTolerance ){ 00165 G4double localR2 = localpos1.x()*localpos1.x()+localpos1.y()*localpos1.y(); 00166 G4double InsideRadius = tubsSolid->GetInnerRadius(); 00167 if (localR2 > (InsideRadius-kCarTolerance)*(InsideRadius-kCarTolerance) 00168 &&localR2 < (InsideRadius+kCarTolerance)*(InsideRadius+kCarTolerance)){ 00169 return fFlux_In; 00170 } 00171 } 00172 00173 if (aStep->GetPostStepPoint()->GetStepStatus() == fGeomBoundary ){ 00174 // Exiting Geometry 00175 G4ThreeVector stppos2= aStep->GetPostStepPoint()->GetPosition(); 00176 G4ThreeVector localpos2 = 00177 theTouchable->GetHistory()->GetTopTransform().TransformPoint(stppos2); 00178 if ( std::fabs(localpos2.z()) > tubsSolid->GetZHalfLength() ) return -1; 00179 //if(std::fabs( localpos2.x()*localpos2.x()+localpos2.y()*localpos2.y() 00180 // - (tubsSolid->GetInnerRadius()*tubsSolid->GetInnerRadius())) 00181 // <kCarTolerance ){ 00182 G4double localR2 = localpos2.x()*localpos2.x()+localpos2.y()*localpos2.y(); 00183 G4double InsideRadius = tubsSolid->GetInnerRadius(); 00184 if (localR2 > (InsideRadius-kCarTolerance)*(InsideRadius-kCarTolerance) 00185 &&localR2 < (InsideRadius+kCarTolerance)*(InsideRadius+kCarTolerance)){ 00186 return fFlux_Out; 00187 } 00188 } 00189 00190 return -1; 00191 }
void G4PSCylinderSurfaceFlux::PrintAll | ( | ) | [virtual] |
Reimplemented from G4VPrimitiveScorer.
Definition at line 211 of file G4PSCylinderSurfaceFlux.cc.
References G4VPrimitiveScorer::detector, G4cout, G4endl, G4VPrimitiveScorer::GetName(), G4VSensitiveDetector::GetName(), G4VPrimitiveScorer::GetUnit(), and G4VPrimitiveScorer::GetUnitValue().
00212 { 00213 G4cout << " MultiFunctionalDet " << detector->GetName() << G4endl; 00214 G4cout << " PrimitiveScorer" << GetName() <<G4endl; 00215 G4cout << " Number of entries " << EvtMap->entries() << G4endl; 00216 std::map<G4int,G4double*>::iterator itr = EvtMap->GetMap()->begin(); 00217 for(; itr != EvtMap->GetMap()->end(); itr++) { 00218 G4cout << " copy no.: " << itr->first 00219 << " flux : " << *(itr->second)/GetUnitValue() 00220 << " ["<<GetUnit()<<"]" 00221 << G4endl; 00222 } 00223 }
G4bool G4PSCylinderSurfaceFlux::ProcessHits | ( | G4Step * | , | |
G4TouchableHistory * | ||||
) | [protected, virtual] |
Implements G4VPrimitiveScorer.
Definition at line 81 of file G4PSCylinderSurfaceFlux.cc.
References G4VSolid::ComputeDimensions(), G4VPVParameterisation::ComputeSolid(), FALSE, fFlux_In, fFlux_InOut, fFlux_Out, G4Tubs::GetDeltaPhiAngle(), G4VPrimitiveScorer::GetIndex(), G4Tubs::GetInnerRadius(), G4VPhysicalVolume::GetLogicalVolume(), G4StepPoint::GetMomentumDirection(), G4VPhysicalVolume::GetParameterisation(), G4StepPoint::GetPhysicalVolume(), G4StepPoint::GetPosition(), G4Step::GetPostStepPoint(), G4Step::GetPreStepPoint(), G4LogicalVolume::GetSolid(), G4StepPoint::GetTouchable(), G4StepPoint::GetTouchableHandle(), G4StepPoint::GetWeight(), G4Tubs::GetZHalfLength(), G4VPrimitiveScorer::indexDepth, IsSelectedSurface(), position, and TRUE.
00082 { 00083 G4StepPoint* preStep = aStep->GetPreStepPoint(); 00084 00085 G4VPhysicalVolume* physVol = preStep->GetPhysicalVolume(); 00086 G4VPVParameterisation* physParam = physVol->GetParameterisation(); 00087 G4VSolid * solid = 0; 00088 if(physParam) 00089 { // for parameterized volume 00090 G4int idx = ((G4TouchableHistory*)(aStep->GetPreStepPoint()->GetTouchable())) 00091 ->GetReplicaNumber(indexDepth); 00092 solid = physParam->ComputeSolid(idx, physVol); 00093 solid->ComputeDimensions(physParam,idx,physVol); 00094 } 00095 else 00096 { // for ordinary volume 00097 solid = physVol->GetLogicalVolume()->GetSolid(); 00098 } 00099 00100 G4Tubs* tubsSolid = (G4Tubs*)(solid); 00101 00102 G4int dirFlag =IsSelectedSurface(aStep,tubsSolid); 00103 00104 if ( dirFlag > 0 ){ 00105 if (fDirection == fFlux_InOut || dirFlag == fDirection ){ 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 G4ThreeVector position = thisStep->GetPosition(); 00121 G4ThreeVector localpos = 00122 theTouchable->GetHistory()->GetTopTransform().TransformAxis(position); 00123 G4double angleFactor = (localdir.x()*localpos.x()+localdir.y()*localpos.y()) 00124 /std::sqrt(localdir.x()*localdir.x() 00125 +localdir.y()*localdir.y()+localdir.z()*localdir.z()) 00126 /std::sqrt(localpos.x()*localpos.x()+localpos.y()*localpos.y()); 00127 00128 if ( angleFactor < 0 ) angleFactor *= -1.; 00129 G4double square = 2.*tubsSolid->GetZHalfLength() 00130 *tubsSolid->GetInnerRadius()* tubsSolid->GetDeltaPhiAngle()/radian; 00131 00132 G4double flux = 1.0; 00133 if ( weighted ) flux *=preStep->GetWeight(); 00134 // Current (Particle Weight) 00135 00136 flux = flux/angleFactor; 00137 if ( divideByArea ) flux /= square; 00138 //Flux with angle. 00139 G4int index = GetIndex(aStep); 00140 EvtMap->add(index,flux); 00141 return TRUE; 00142 }else{ 00143 return FALSE; 00144 } 00145 }else{ 00146 return FALSE; 00147 } 00148 }
void G4PSCylinderSurfaceFlux::SetUnit | ( | const G4String & | unit | ) | [virtual] |
Reimplemented from G4VPrimitiveScorer.
Definition at line 225 of file G4PSCylinderSurfaceFlux.cc.
References G4VPrimitiveScorer::CheckAndSetUnit(), G4Exception(), G4VPrimitiveScorer::GetName(), G4VPrimitiveScorer::GetUnit(), JustWarning, G4VPrimitiveScorer::unitName, and G4VPrimitiveScorer::unitValue.
Referenced by G4PSCylinderSurfaceFlux(), and G4PSCylinderSurfaceFlux3D::G4PSCylinderSurfaceFlux3D().
00226 { 00227 if ( divideByArea ) { 00228 CheckAndSetUnit(unit,"Per Unit Surface"); 00229 } else { 00230 if (unit == "" ){ 00231 unitName = unit; 00232 unitValue = 1.0; 00233 }else{ 00234 G4String msg = "Invalid unit ["+unit+"] (Current unit is [" +GetUnit()+"] ) for " + GetName(); 00235 G4Exception("G4PSCylinderSurfaceFlux::SetUnit","DetPS0003",JustWarning,msg); 00236 } 00237 } 00238 }
void G4PSCylinderSurfaceFlux::Weighted | ( | G4bool | flg = true |
) | [inline] |