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
00032
00033
00034
00035
00036
00037 #include "G4VScoringMesh.hh"
00038
00039 #include "G4SystemOfUnits.hh"
00040 #include "G4VPhysicalVolume.hh"
00041 #include "G4MultiFunctionalDetector.hh"
00042 #include "G4VPrimitiveScorer.hh"
00043 #include "G4VSDFilter.hh"
00044 #include "G4SDManager.hh"
00045
00046 G4VScoringMesh::G4VScoringMesh(const G4String& wName)
00047 : fWorldName(wName),fCurrentPS(0),fConstructed(false),fActive(true),
00048 fRotationMatrix(0), fMFD(new G4MultiFunctionalDetector(wName)),
00049 verboseLevel(0),sizeIsSet(false),nMeshIsSet(false),
00050 fDrawUnit(""), fDrawUnitValue(1.)
00051 {
00052 G4SDManager::GetSDMpointer()->AddNewDetector(fMFD);
00053
00054 fSize[0] = fSize[1] = fSize[2] = 0.;
00055 fNSegment[0] = fNSegment[1] = fNSegment[2] = 1;
00056 fDivisionAxisNames[0] = fDivisionAxisNames[1] = fDivisionAxisNames[2] = "";
00057 }
00058
00059 G4VScoringMesh::~G4VScoringMesh()
00060 {
00061 ;
00062 }
00063
00064 void G4VScoringMesh::ResetScore() {
00065 if(verboseLevel > 9) G4cout << "G4VScoringMesh::ResetScore() is called." << G4endl;
00066 std::map<G4String, G4THitsMap<G4double>* >::iterator itr = fMap.begin();
00067 for(; itr != fMap.end(); itr++) {
00068 if(verboseLevel > 9) G4cout << "G4VScoringMesh::ResetScore()" << itr->first << G4endl;
00069 itr->second->clear();
00070 }
00071 }
00072 void G4VScoringMesh::SetSize(G4double size[3]) {
00073 if ( !sizeIsSet ){
00074 for(int i = 0; i < 3; i++) fSize[i] = size[i];
00075 sizeIsSet = true;
00076 }else{
00077 G4String message = " The size of scoring mesh can not be changed.";
00078 G4Exception("G4VScoringMesh::SetSize()",
00079 "DigiHitsUtilsScoreVScoringMesh000", JustWarning,
00080 message);
00081 }
00082 }
00083 G4ThreeVector G4VScoringMesh::GetSize() const {
00084 if(sizeIsSet)
00085 return G4ThreeVector(fSize[0], fSize[1], fSize[2]);
00086 else
00087 return G4ThreeVector(0., 0., 0.);
00088 }
00089 void G4VScoringMesh::SetCenterPosition(G4double centerPosition[3]) {
00090 fCenterPosition = G4ThreeVector(centerPosition[0], centerPosition[1], centerPosition[2]);
00091 }
00092 void G4VScoringMesh::SetNumberOfSegments(G4int nSegment[3]) {
00093 if ( !nMeshIsSet ){
00094 for(int i = 0; i < 3; i++) fNSegment[i] = nSegment[i];
00095 nMeshIsSet = true;
00096 } else {
00097 G4String message = " The size of scoring segments can not be changed.";
00098 G4Exception("G4VScoringMesh::SetNumberOfSegments()",
00099 "DigiHitsUtilsScoreVScoringMesh000", JustWarning,
00100 message);
00101 }
00102 }
00103 void G4VScoringMesh::GetNumberOfSegments(G4int nSegment[3]) {
00104 for(int i = 0; i < 3; i++) nSegment[i] = fNSegment[i];
00105 }
00106 void G4VScoringMesh::RotateX(G4double delta) {
00107 if(fRotationMatrix == 0) fRotationMatrix = new G4RotationMatrix();
00108 fRotationMatrix->rotateX(delta);
00109 }
00110
00111 void G4VScoringMesh::RotateY(G4double delta) {
00112 if(fRotationMatrix == 0) fRotationMatrix = new G4RotationMatrix();
00113 fRotationMatrix->rotateY(delta);
00114 }
00115
00116 void G4VScoringMesh::RotateZ(G4double delta) {
00117 if(fRotationMatrix == 0) fRotationMatrix = new G4RotationMatrix();
00118 fRotationMatrix->rotateZ(delta);
00119 }
00120
00121 void G4VScoringMesh::SetPrimitiveScorer(G4VPrimitiveScorer * ps) {
00122
00123 if(!ReadyForQuantity())
00124 {
00125 G4cerr << "ERROR : G4VScoringMesh::SetPrimitiveScorer() : " << ps->GetName()
00126 << " does not yet have mesh size or number of bins. Set them first." << G4endl
00127 << "This Method is ignored." << G4endl;
00128 return;
00129 }
00130 if(verboseLevel > 0) G4cout << "G4VScoringMesh::SetPrimitiveScorer() : "
00131 << ps->GetName() << " is registered."
00132 << " 3D size: ("
00133 << fNSegment[0] << ", "
00134 << fNSegment[1] << ", "
00135 << fNSegment[2] << ")" << G4endl;
00136
00137 ps->SetNijk(fNSegment[0], fNSegment[1], fNSegment[2]);
00138 fCurrentPS = ps;
00139 fMFD->RegisterPrimitive(ps);
00140 G4THitsMap<G4double> * map = new G4THitsMap<G4double>(fWorldName, ps->GetName());
00141 fMap[ps->GetName()] = map;
00142 }
00143
00144 void G4VScoringMesh::SetFilter(G4VSDFilter * filter) {
00145
00146 if(fCurrentPS == 0) {
00147 G4cerr << "ERROR : G4VScoringMesh::SetSDFilter() : a quantity must be defined first. This method is ignored." << G4endl;
00148 return;
00149 }
00150 if(verboseLevel > 0) G4cout << "G4VScoringMesh::SetFilter() : "
00151 << filter->GetName()
00152 << " is set to "
00153 << fCurrentPS->GetName() << G4endl;
00154
00155 G4VSDFilter* oldFilter = fCurrentPS->GetFilter();
00156 if(oldFilter)
00157 {
00158 G4cout << "WARNING : G4VScoringMesh::SetFilter() : " << oldFilter->GetName()
00159 << " is overwritten by " << filter->GetName() << G4endl;
00160 }
00161 fCurrentPS->SetFilter(filter);
00162 }
00163
00164 void G4VScoringMesh::SetCurrentPrimitiveScorer(const G4String & name) {
00165 fCurrentPS = GetPrimitiveScorer(name);
00166 if(fCurrentPS == 0) {
00167 G4cerr << "ERROR : G4VScoringMesh::SetCurrentPrimitiveScorer() : The primitive scorer <"
00168 << name << "> does not found." << G4endl;
00169 }
00170 }
00171
00172 G4bool G4VScoringMesh::FindPrimitiveScorer(const G4String & psname) {
00173 std::map<G4String, G4THitsMap<G4double>* >::iterator itr = fMap.find(psname);;
00174 if(itr == fMap.end()) return false;
00175 return true;
00176 }
00177
00178 G4String G4VScoringMesh::GetPSUnit(const G4String & psname) {
00179 std::map<G4String, G4THitsMap<G4double>* >::iterator itr = fMap.find(psname);;
00180 if(itr == fMap.end()) {
00181 return G4String("");
00182 } else {
00183 return GetPrimitiveScorer(psname)->GetUnit();
00184 }
00185 }
00186
00187 G4String G4VScoringMesh::GetCurrentPSUnit(){
00188 G4String unit = "";
00189 if(fCurrentPS == 0) {
00190 G4String msg = "ERROR : G4VScoringMesh::GetCurrentPSUnit() : ";
00191 msg += " Current primitive scorer is null.";
00192 G4cerr << msg << G4endl;
00193 }else{
00194 unit = fCurrentPS->GetUnit();
00195 }
00196 return unit;
00197 }
00198
00199 void G4VScoringMesh::SetCurrentPSUnit(const G4String& unit){
00200 if(fCurrentPS == 0) {
00201 G4String msg = "ERROR : G4VScoringMesh::GetCurrentPSUnit() : ";
00202 msg += " Current primitive scorer is null.";
00203 G4cerr << msg << G4endl;
00204 }else{
00205 fCurrentPS->SetUnit(unit);
00206 }
00207 }
00208
00209 G4double G4VScoringMesh::GetPSUnitValue(const G4String & psname) {
00210 std::map<G4String, G4THitsMap<G4double>* >::iterator itr = fMap.find(psname);;
00211 if(itr == fMap.end()) {
00212 return 1.;
00213 } else {
00214 return GetPrimitiveScorer(psname)->GetUnitValue();
00215 }
00216 }
00217
00218 void G4VScoringMesh::GetDivisionAxisNames(G4String divisionAxisNames[3]) {
00219 for(int i = 0; i < 3; i++) divisionAxisNames[i] = fDivisionAxisNames[i];
00220 }
00221
00222 G4VPrimitiveScorer * G4VScoringMesh::GetPrimitiveScorer(const G4String & name) {
00223 if(fMFD == 0) return 0;
00224
00225 G4int nps = fMFD->GetNumberOfPrimitives();
00226 for(G4int i = 0; i < nps; i++) {
00227 G4VPrimitiveScorer * ps = fMFD->GetPrimitive(i);
00228 if(name == ps->GetName()) return ps;
00229 }
00230
00231 return 0;
00232 }
00233 void G4VScoringMesh::List() const {
00234
00235 G4cout << " # of segments: ("
00236 << fNSegment[0] << ", "
00237 << fNSegment[1] << ", "
00238 << fNSegment[2] << ")"
00239 << G4endl;
00240 G4cout << " displacement: ("
00241 << fCenterPosition.x()/cm << ", "
00242 << fCenterPosition.y()/cm << ", "
00243 << fCenterPosition.z()/cm << ") [cm]"
00244 << G4endl;
00245 if(fRotationMatrix != 0) {
00246 G4cout << " rotation matrix: "
00247 << fRotationMatrix->xx() << " "
00248 << fRotationMatrix->xy() << " "
00249 << fRotationMatrix->xz() << G4endl
00250 << " "
00251 << fRotationMatrix->yx() << " "
00252 << fRotationMatrix->yy() << " "
00253 << fRotationMatrix->yz() << G4endl
00254 << " "
00255 << fRotationMatrix->zx() << " "
00256 << fRotationMatrix->zy() << " "
00257 << fRotationMatrix->zz() << G4endl;
00258 }
00259
00260
00261 G4cout << " registered primitve scorers : " << G4endl;
00262 G4int nps = fMFD->GetNumberOfPrimitives();
00263 G4VPrimitiveScorer * ps;
00264 for(int i = 0; i < nps; i++) {
00265 ps = fMFD->GetPrimitive(i);
00266 G4cout << " " << i << " " << ps->GetName();
00267 if(ps->GetFilter() != 0) G4cout << " with " << ps->GetFilter()->GetName();
00268 G4cout << G4endl;
00269 }
00270
00271
00272 }
00273
00274 void G4VScoringMesh::Dump() {
00275 G4cout << "scoring mesh name: " << fWorldName << G4endl;
00276
00277 G4cout << "# of G4THitsMap : " << fMap.size() << G4endl;
00278 std::map<G4String, G4THitsMap<G4double>* >::iterator itr = fMap.begin();
00279 for(; itr != fMap.end(); itr++) {
00280 G4cout << "[" << itr->first << "]" << G4endl;
00281 itr->second->PrintAllHits();
00282 }
00283 G4cout << G4endl;
00284
00285 }
00286
00287
00288 void G4VScoringMesh::DrawMesh(const G4String& psName,G4VScoreColorMap* colorMap,G4int axflg)
00289 {
00290 fDrawPSName = psName;
00291 std::map<G4String, G4THitsMap<G4double>* >::const_iterator fMapItr = fMap.find(psName);
00292 if(fMapItr!=fMap.end()) {
00293 fDrawUnit = GetPSUnit(psName);
00294 fDrawUnitValue = GetPSUnitValue(psName);
00295 Draw(fMapItr->second->GetMap(), colorMap,axflg);
00296 } else {
00297 G4cerr << "Scorer <" << psName << "> is not defined. Method ignored." << G4endl;
00298 }
00299 }
00300
00301 void G4VScoringMesh::DrawMesh(const G4String& psName,G4int idxPlane,G4int iColumn,G4VScoreColorMap* colorMap)
00302 {
00303 fDrawPSName = psName;
00304 std::map<G4String, G4THitsMap<G4double>* >::const_iterator fMapItr = fMap.find(psName);
00305 if(fMapItr!=fMap.end()) {
00306 fDrawUnit = GetPSUnit(psName);
00307 fDrawUnitValue = GetPSUnitValue(psName);
00308 DrawColumn(fMapItr->second->GetMap(),colorMap,idxPlane,iColumn);
00309 } else {
00310 G4cerr << "Scorer <" << psName << "> is not defined. Method ignored." << G4endl;
00311 }
00312 }
00313