00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031 #include "G4PSCellFlux.hh"
00032
00033 #include "G4SystemOfUnits.hh"
00034 #include "G4Track.hh"
00035 #include "G4VSolid.hh"
00036 #include "G4VPhysicalVolume.hh"
00037 #include "G4VPVParameterisation.hh"
00038 #include "G4UnitsTable.hh"
00040
00041
00042
00043
00044
00045
00046
00047
00048
00049
00050
00051
00052
00053
00055
00056 G4PSCellFlux::G4PSCellFlux(G4String name, G4int depth)
00057 :G4VPrimitiveScorer(name,depth),HCID(-1),weighted(true)
00058 {
00059 DefineUnitAndCategory();
00060 SetUnit("percm2");
00061
00062 }
00063
00064 G4PSCellFlux::G4PSCellFlux(G4String name, const G4String& unit, G4int depth)
00065 :G4VPrimitiveScorer(name,depth),HCID(-1),weighted(true)
00066 {
00067 DefineUnitAndCategory();
00068 SetUnit(unit);
00069 }
00070
00071 G4PSCellFlux::~G4PSCellFlux()
00072 {;}
00073
00074 G4bool G4PSCellFlux::ProcessHits(G4Step* aStep,G4TouchableHistory*)
00075 {
00076 G4double stepLength = aStep->GetStepLength();
00077 if ( stepLength == 0. ) return FALSE;
00078
00079 G4int idx = ((G4TouchableHistory*)
00080 (aStep->GetPreStepPoint()->GetTouchable()))
00081 ->GetReplicaNumber(indexDepth);
00082 G4double cubicVolume = ComputeVolume(aStep, idx);
00083
00084 G4double CellFlux = stepLength / cubicVolume;
00085 if (weighted) CellFlux *= aStep->GetPreStepPoint()->GetWeight();
00086 G4int index = GetIndex(aStep);
00087 EvtMap->add(index,CellFlux);
00088
00089 return TRUE;
00090 }
00091
00092 void G4PSCellFlux::Initialize(G4HCofThisEvent* HCE)
00093 {
00094 EvtMap = new G4THitsMap<G4double>(detector->GetName(),
00095 GetName());
00096 if ( HCID < 0 ) HCID = GetCollectionID(0);
00097 HCE->AddHitsCollection(HCID,EvtMap);
00098 }
00099
00100 void G4PSCellFlux::EndOfEvent(G4HCofThisEvent*)
00101 {;}
00102
00103 void G4PSCellFlux::clear(){
00104 EvtMap->clear();
00105 }
00106
00107 void G4PSCellFlux::DrawAll()
00108 {;}
00109
00110 void G4PSCellFlux::PrintAll()
00111 {
00112 G4cout << " MultiFunctionalDet " << detector->GetName() << G4endl;
00113 G4cout << " PrimitiveScorer " << GetName() <<G4endl;
00114 G4cout << " Number of entries " << EvtMap->entries() << G4endl;
00115 std::map<G4int,G4double*>::iterator itr = EvtMap->GetMap()->begin();
00116 for(; itr != EvtMap->GetMap()->end(); itr++) {
00117 G4cout << " copy no.: " << itr->first
00118 << " cell flux : " << *(itr->second)/GetUnitValue()
00119 << " [" << GetUnit() << "]"
00120 << G4endl;
00121 }
00122 }
00123
00124 void G4PSCellFlux::SetUnit(const G4String& unit)
00125 {
00126 CheckAndSetUnit(unit,"Per Unit Surface");
00127 }
00128
00129 void G4PSCellFlux::DefineUnitAndCategory(){
00130
00131 new G4UnitDefinition("percentimeter2","percm2","Per Unit Surface",(1./cm2));
00132 new G4UnitDefinition("permillimeter2","permm2","Per Unit Surface",(1./mm2));
00133 new G4UnitDefinition("permeter2","perm2","Per Unit Surface",(1./m2));
00134 }
00135
00136 G4double G4PSCellFlux::ComputeVolume(G4Step* aStep, G4int idx){
00137
00138 G4VPhysicalVolume* physVol = aStep->GetPreStepPoint()->GetPhysicalVolume();
00139 G4VPVParameterisation* physParam = physVol->GetParameterisation();
00140 G4VSolid* solid = 0;
00141 if(physParam)
00142 {
00143 if(idx<0)
00144 {
00145 G4ExceptionDescription ED;
00146 ED << "Incorrect replica number --- GetReplicaNumber : " << idx << G4endl;
00147 G4Exception("G4PSCellFlux::ComputeVolume","DetPS0001",JustWarning,ED);
00148 }
00149 solid = physParam->ComputeSolid(idx, physVol);
00150 solid->ComputeDimensions(physParam,idx,physVol);
00151 }
00152 else
00153 {
00154 solid = physVol->GetLogicalVolume()->GetSolid();
00155 }
00156
00157 return solid->GetCubicVolume();
00158 }