G4ScoringCylinder Class Reference

#include <G4ScoringCylinder.hh>

Inheritance diagram for G4ScoringCylinder:

G4VScoringMesh

Public Member Functions

 G4ScoringCylinder (G4String wName)
 ~G4ScoringCylinder ()
virtual void Construct (G4VPhysicalVolume *fWorldPhys)
virtual void List () const
virtual void Draw (std::map< G4int, G4double * > *map, G4VScoreColorMap *colorMap, G4int axflg=111)
virtual void DrawColumn (std::map< G4int, G4double * > *map, G4VScoreColorMap *colorMap, G4int idxProj, G4int idxColumn)
void SetRMax (G4double rMax)
void SetZSize (G4double zSize)
void RegisterPrimitives (std::vector< G4VPrimitiveScorer * > &vps)
void GetRZPhi (G4int index, G4int q[3]) const

Detailed Description

Definition at line 41 of file G4ScoringCylinder.hh.


Constructor & Destructor Documentation

G4ScoringCylinder::G4ScoringCylinder ( G4String  wName  ) 

Definition at line 55 of file G4ScoringCylinder.cc.

References cylinderMesh, G4VScoringMesh::fDivisionAxisNames, and G4VScoringMesh::fShape.

00056   :G4VScoringMesh(wName), fMeshElementLogical(0)
00057 {
00058   fShape = cylinderMesh;
00059 
00060   fDivisionAxisNames[0] = "Z";
00061   fDivisionAxisNames[1] = "PHI";
00062   fDivisionAxisNames[2] = "R";
00063 }

G4ScoringCylinder::~G4ScoringCylinder (  ) 

Definition at line 65 of file G4ScoringCylinder.cc.

00066 {;}


Member Function Documentation

void G4ScoringCylinder::Construct ( G4VPhysicalVolume fWorldPhys  )  [virtual]

Implements G4VScoringMesh.

Definition at line 68 of file G4ScoringCylinder.cc.

References G4VScoringMesh::fConstructed, G4cout, G4endl, G4VPhysicalVolume::GetName(), G4VScoringMesh::ResetScore(), and G4VScoringMesh::verboseLevel.

00069 {
00070   if(fConstructed) {
00071 
00072     if(verboseLevel > 0) 
00073       G4cout << fWorldPhys->GetName() << " --- All quantities are reset." << G4endl;
00074     ResetScore();
00075 
00076   } else {
00077     fConstructed = true;
00078     SetupGeometry(fWorldPhys);
00079   }
00080 }

void G4ScoringCylinder::Draw ( std::map< G4int, G4double * > *  map,
G4VScoreColorMap colorMap,
G4int  axflg = 111 
) [virtual]

Implements G4VScoringMesh.

Definition at line 225 of file G4ScoringCylinder.cc.

References DBL_MAX, G4VVisManager::Draw(), G4VScoreColorMap::DrawColorChart(), G4VScoringMesh::fCenterPosition, G4VScoringMesh::fDrawPSName, G4VScoringMesh::fDrawUnit, G4VScoringMesh::fDrawUnitValue, G4VScoringMesh::fNSegment, G4VScoringMesh::fRotationMatrix, G4VScoringMesh::fSize, G4cout, G4VVisManager::GetConcreteInstance(), G4VScoreColorMap::GetMapColor(), GetRZPhi(), G4VScoreColorMap::IfFloatMinMax(), G4VisAttributes::SetColour(), G4VisAttributes::SetForceAuxEdgeVisible(), G4VisAttributes::SetForceSolid(), G4VScoreColorMap::SetMinMax(), G4VScoreColorMap::SetPSName(), and G4VScoreColorMap::SetPSUnit().

00225                                                                                                     {
00226 
00227   G4VVisManager * pVisManager = G4VVisManager::GetConcreteInstance();
00228   if(pVisManager) {
00229     
00230     // cell vectors
00231     std::vector<double> ephi;
00232     for(int phi = 0; phi < fNSegment[IPHI]; phi++) ephi.push_back(0.);
00233     //-
00234     std::vector<std::vector<double> > zphicell; // zphicell[Z][PHI]
00235     for(int z = 0; z < fNSegment[IZ]; z++) zphicell.push_back(ephi);
00236     //-
00237     std::vector<std::vector<double> > rphicell; // rphicell[R][PHI]
00238     for(int r = 0; r < fNSegment[IR]; r++) rphicell.push_back(ephi);
00239 
00240     // projections
00241     G4int q[3];
00242     std::map<G4int, G4double*>::iterator itr = map->begin();
00243     for(; itr != map->end(); itr++) {
00244       if(itr->first < 0) {
00245         G4cout << itr->first << G4endl;
00246         continue;
00247       }
00248       GetRZPhi(itr->first, q);
00249 
00250       zphicell[q[IZ]][q[IPHI]] += *(itr->second)/fDrawUnitValue;
00251       rphicell[q[IR]][q[IPHI]] += *(itr->second)/fDrawUnitValue;
00252     }  
00253     
00254     // search min./max. values
00255     G4double zphimin = DBL_MAX, rphimin = DBL_MAX;
00256     G4double zphimax = 0., rphimax = 0.;
00257     for(int iphi = 0; iphi < fNSegment[IPHI]; iphi++) {
00258       for(int iz = 0; iz < fNSegment[IZ]; iz++) {
00259         if(zphimin > zphicell[iz][iphi]) zphimin = zphicell[iz][iphi];
00260         if(zphimax < zphicell[iz][iphi]) zphimax = zphicell[iz][iphi];
00261       }
00262       for(int ir = 0; ir < fNSegment[IR]; ir++) {
00263         if(rphimin > rphicell[ir][iphi]) rphimin = rphicell[ir][iphi];
00264         if(rphimax < rphicell[ir][iphi]) rphimax = rphicell[ir][iphi];
00265       }
00266     }
00267 
00268     G4VisAttributes att;
00269     att.SetForceSolid(true);
00270     att.SetForceAuxEdgeVisible(true);
00271 
00272 
00273     G4Scale3D scale;
00274     if(axflg/100==1) {
00275       // rz plane
00276     }
00277     axflg = axflg%100;
00278     if(axflg/10==1) {
00279       // z-phi plane
00280       if(colorMap->IfFloatMinMax()) { colorMap->SetMinMax(zphimin, zphimax); }
00281 
00282       G4double zhalf = fSize[1]/fNSegment[IZ];
00283       for(int phi = 0; phi < fNSegment[IPHI]; phi++) {
00284         for(int z = 0; z < fNSegment[IZ]; z++) {
00285           //-
00286           G4double angle = twopi/fNSegment[IPHI]*phi;
00287           G4double dphi = twopi/fNSegment[IPHI];
00288           G4Tubs cylinder("z-phi",                    // name
00289                           fSize[0]*0.99, fSize[0],  // inner radius, outer radius
00290                           zhalf,                      // half length in z
00291                           angle, dphi*0.99999);       // starting phi angle, delta angle
00292           //-
00293           G4ThreeVector zpos(0., 0., -fSize[1] + fSize[1]/fNSegment[IZ]*(1 + 2.*z));
00294           G4Transform3D trans;
00295           if(fRotationMatrix) {
00296             trans = G4Rotate3D(*fRotationMatrix).inverse()*G4Translate3D(zpos);
00297             trans = G4Translate3D(fCenterPosition)*trans;
00298           } else {
00299             trans = G4Translate3D(zpos)*G4Translate3D(fCenterPosition);
00300           }
00301           G4double c[4];
00302           colorMap->GetMapColor(zphicell[z][phi], c);
00303           att.SetColour(c[0], c[1], c[2]);//, c[3]);
00304           //-
00305           pVisManager->Draw(cylinder, att, trans);
00306         }
00307       }
00308     }
00309     axflg = axflg%10;
00310     if(axflg==1) {
00311       // r-phi plane
00312       if(colorMap->IfFloatMinMax()) { colorMap->SetMinMax(rphimin, rphimax); }
00313 
00314       G4double rsize = fSize[0]/fNSegment[IR];
00315       for(int phi = 0; phi < fNSegment[IPHI]; phi++) {
00316         for(int r = 0; r < fNSegment[IR]; r++) {
00317 
00318           G4double rs[2] = {rsize*r, rsize*(r+1)};
00319           G4double angle = twopi/fNSegment[IPHI]*phi;
00320           G4double dphi = twopi/fNSegment[IPHI];
00321           G4Tubs cylinder("z-phi", rs[0], rs[1], 0.001,
00322                           angle, dphi*0.99999);
00323           /*
00324           G4cout << ">>>> "
00325                  << rs[0] << " - " << rs[1] << " : "
00326                  << angle << " - " << angle + dphi
00327                  << G4endl;
00328           */
00329 
00330           G4ThreeVector zposn(0., 0., -fSize[1]);
00331           G4ThreeVector zposp(0., 0.,  fSize[1]);
00332           G4Transform3D transn, transp;
00333           if(fRotationMatrix) {
00334             transn = G4Rotate3D(*fRotationMatrix).inverse()*G4Translate3D(zposn);
00335             transn = G4Translate3D(fCenterPosition)*transn;
00336             transp = G4Rotate3D(*fRotationMatrix).inverse()*G4Translate3D(zposp);
00337             transp = G4Translate3D(fCenterPosition)*transp;
00338           } else {
00339             transn = G4Translate3D(zposn)*G4Translate3D(fCenterPosition);
00340             transp = G4Translate3D(zposp)*G4Translate3D(fCenterPosition);
00341           }
00342           G4double c[4];
00343           colorMap->GetMapColor(rphicell[r][phi], c);
00344           att.SetColour(c[0], c[1], c[2]);//, c[3]);
00345           /*
00346           G4cout << "   " << c[0] << ", "
00347                  << c[1] << ", " << c[2] << G4endl;
00348           */
00349           pVisManager->Draw(cylinder, att, transn);
00350           pVisManager->Draw(cylinder, att, transp);
00351         }
00352       }
00353     }
00354 
00355     colorMap->SetPSUnit(fDrawUnit);
00356     colorMap->SetPSName(fDrawPSName);
00357     colorMap->DrawColorChart();
00358 
00359   }
00360 }

void G4ScoringCylinder::DrawColumn ( std::map< G4int, G4double * > *  map,
G4VScoreColorMap colorMap,
G4int  idxProj,
G4int  idxColumn 
) [virtual]

Implements G4VScoringMesh.

Definition at line 362 of file G4ScoringCylinder.cc.

References DBL_MAX, G4VVisManager::Draw(), G4VScoreColorMap::DrawColorChart(), G4VScoringMesh::fCenterPosition, G4VScoringMesh::fDrawPSName, G4VScoringMesh::fDrawUnit, G4VScoringMesh::fDrawUnitValue, G4VScoringMesh::fNSegment, G4VScoringMesh::fRotationMatrix, G4VScoringMesh::fSize, G4cerr, G4cout, G4VVisManager::GetConcreteInstance(), G4VScoreColorMap::GetMapColor(), GetRZPhi(), G4VScoreColorMap::IfFloatMinMax(), G4VisAttributes::SetColour(), G4VisAttributes::SetForceAuxEdgeVisible(), G4VisAttributes::SetForceSolid(), G4VScoreColorMap::SetMinMax(), G4VScoreColorMap::SetPSName(), and G4VScoreColorMap::SetPSUnit().

00364 {
00365   G4int projAxis = 0;
00366   switch(idxProj) {
00367   case 0: 
00368     projAxis = IR;
00369     break;
00370   case 1:
00371     projAxis = IZ;
00372     break;
00373   case 2:
00374     projAxis = IPHI;
00375     break;
00376   }
00377 
00378   if(idxColumn<0 || idxColumn>=fNSegment[projAxis])
00379   {
00380     G4cerr << "Warning : Column number " << idxColumn << " is out of scoring mesh [0," << fNSegment[projAxis]-1 <<
00381     "]. Method ignored." << G4endl;
00382     return;
00383   }
00384   G4VVisManager * pVisManager = G4VVisManager::GetConcreteInstance();
00385   if(pVisManager) {
00386 
00387     // cell vectors
00388     std::vector<std::vector<std::vector<double> > > cell; // cell[R][Z][PHI]
00389     std::vector<double> ephi;
00390     for(int phi = 0; phi < fNSegment[IPHI]; phi++) ephi.push_back(0.);
00391     std::vector<std::vector<double> > ezphi;
00392     for(int z = 0; z < fNSegment[IZ]; z++) ezphi.push_back(ephi);
00393     for(int r = 0; r < fNSegment[IR]; r++) cell.push_back(ezphi);
00394 
00395     std::vector<std::vector<double> > rzcell; // rzcell[R][Z]
00396     std::vector<double> ez;
00397     for(int z = 0; z < fNSegment[IZ]; z++) ez.push_back(0.);
00398     for(int r = 0; r < fNSegment[IR]; r++) rzcell.push_back(ez);
00399 
00400     std::vector<std::vector<double> > zphicell; // zphicell[Z][PHI]
00401     for(int z = 0; z < fNSegment[IZ]; z++) zphicell.push_back(ephi);
00402 
00403     std::vector<std::vector<double> > rphicell; // rphicell[R][PHI]
00404     for(int r = 0; r < fNSegment[IR]; r++) rphicell.push_back(ephi);
00405 
00406     // projections
00407     G4int q[3];
00408     std::map<G4int, G4double*>::iterator itr = map->begin();
00409     for(; itr != map->end(); itr++) {
00410       if(itr->first < 0) {
00411         G4cout << itr->first << G4endl;
00412         continue;
00413       }
00414       GetRZPhi(itr->first, q);
00415 
00416       if(projAxis == IR && q[IR] == idxColumn) { // zphi plane
00417         zphicell[q[IZ]][q[IPHI]] += *(itr->second)/fDrawUnitValue;
00418       }
00419       if(projAxis == IZ && q[IZ] == idxColumn) { // rphi plane
00420         rphicell[q[IR]][q[IPHI]] += *(itr->second)/fDrawUnitValue;
00421       }
00422       if(projAxis == IPHI && q[IPHI] == idxColumn) { // rz plane
00423         rzcell[q[IR]][q[IZ]] += *(itr->second)/fDrawUnitValue; 
00424       }
00425     }  
00426 
00427     // search min./max. values
00428     G4double rzmin = DBL_MAX, zphimin = DBL_MAX, rphimin = DBL_MAX;
00429     G4double rzmax = 0., zphimax = 0., rphimax = 0.;
00430     for(int r = 0; r < fNSegment[IR]; r++) {
00431       for(int phi = 0; phi < fNSegment[IPHI]; phi++) {
00432         if(rphimin > rphicell[r][phi]) rphimin = rphicell[r][phi];
00433         if(rphimax < rphicell[r][phi]) rphimax = rphicell[r][phi];
00434       }
00435       for(int z = 0; z < fNSegment[IZ]; z++) {
00436         if(rzmin > rzcell[r][z]) rzmin = rzcell[r][z];
00437         if(rzmax < rzcell[r][z]) rzmax = rzcell[r][z];
00438       }
00439     }
00440     for(int z = 0; z < fNSegment[IZ]; z++) {
00441       for(int phi = 0; phi < fNSegment[IPHI]; phi++) {
00442         if(zphimin > zphicell[z][phi]) zphimin = zphicell[z][phi];
00443         if(zphimax < zphicell[z][phi]) zphimax = zphicell[z][phi];
00444       }
00445     }
00446 
00447 
00448     G4VisAttributes att;
00449     att.SetForceSolid(true);
00450     att.SetForceAuxEdgeVisible(true);
00451 
00452 
00453     G4Scale3D scale;
00454     // z-phi plane
00455     if(projAxis == IR) {
00456       if(colorMap->IfFloatMinMax()) { colorMap->SetMinMax(zphimin,zphimax); }
00457 
00458       G4double zhalf = fSize[1]/fNSegment[IZ];
00459       G4double rsize[2] = {fSize[0]/fNSegment[IR]*idxColumn,
00460                             fSize[0]/fNSegment[IR]*(idxColumn+1)};
00461       for(int phi = 0; phi < fNSegment[IPHI]; phi++) {
00462         for(int z = 0; z < fNSegment[IZ]; z++) {
00463 
00464           G4double angle = twopi/fNSegment[IPHI]*phi*radian;
00465           G4double dphi = twopi/fNSegment[IPHI]*radian;
00466           G4Tubs cylinder("z-phi", rsize[0], rsize[1], zhalf,
00467                           angle, dphi*0.99999);
00468 
00469           G4ThreeVector zpos(0., 0., -fSize[1] + fSize[1]/fNSegment[IZ]*(1 + 2.*z));
00470           G4Transform3D trans;
00471           if(fRotationMatrix) {
00472             trans = G4Rotate3D(*fRotationMatrix).inverse()*G4Translate3D(zpos);
00473             trans = G4Translate3D(fCenterPosition)*trans;
00474           } else {
00475             trans = G4Translate3D(zpos)*G4Translate3D(fCenterPosition);
00476           }
00477           G4double c[4];
00478           colorMap->GetMapColor(zphicell[z][phi], c);
00479           att.SetColour(c[0], c[1], c[2]);//, c[3]);
00480           pVisManager->Draw(cylinder, att, trans);
00481         }
00482       }
00483 
00484       // r-phi plane
00485     } else if(projAxis == IZ) {
00486       if(colorMap->IfFloatMinMax()) { colorMap->SetMinMax(rphimin,rphimax); }
00487 
00488       G4double rsize = fSize[0]/fNSegment[IR];
00489       for(int phi = 0; phi < fNSegment[IPHI]; phi++) {
00490         for(int r = 0; r < fNSegment[IR]; r++) {
00491 
00492           G4double rs[2] = {rsize*r, rsize*(r+1)};
00493           G4double angle = twopi/fNSegment[IPHI]*phi*radian;
00494           G4double dz = fSize[1]/fNSegment[IZ];
00495           G4double dphi = twopi/fNSegment[IPHI]*radian;
00496           G4Tubs cylinder("r-phi", rs[0], rs[1], dz,
00497                           angle, dphi*0.99999);
00498           G4ThreeVector zpos(0., 0.,
00499                              -fSize[1]+fSize[1]/fNSegment[IZ]*(idxColumn*2+1));
00500           G4Transform3D trans;
00501           if(fRotationMatrix) {
00502             trans = G4Rotate3D(*fRotationMatrix).inverse()*G4Translate3D(zpos);
00503             trans = G4Translate3D(fCenterPosition)*trans;
00504           } else {
00505             trans = G4Translate3D(zpos)*G4Translate3D(fCenterPosition);
00506           }
00507           G4double c[4];
00508           colorMap->GetMapColor(rphicell[r][phi], c);
00509           att.SetColour(c[0], c[1], c[2]);//, c[3]);
00510           pVisManager->Draw(cylinder, att, trans);
00511         }
00512       }
00513 
00514       // r-z plane
00515     } else if(projAxis == IPHI) {
00516       if(colorMap->IfFloatMinMax()) { colorMap->SetMinMax(rzmin,rzmax); }
00517 
00518       G4double rsize = fSize[0]/fNSegment[IR];
00519       G4double zhalf = fSize[1]/fNSegment[IZ];
00520       G4double angle = twopi/fNSegment[IPHI]*idxColumn*radian;
00521       G4double dphi = twopi/fNSegment[IPHI]*radian;
00522       for(int z = 0; z < fNSegment[IZ]; z++) {
00523         for(int r = 0; r < fNSegment[IR]; r++) {
00524 
00525           G4double rs[2] = {rsize*r, rsize*(r+1)};
00526           G4Tubs cylinder("z-phi", rs[0], rs[1], zhalf,
00527                           angle, dphi);
00528 
00529           G4ThreeVector zpos(0., 0.,
00530                              -fSize[1]+fSize[1]/fNSegment[IZ]*(2.*z+1));
00531           G4Transform3D trans;
00532           if(fRotationMatrix) {
00533             trans = G4Rotate3D(*fRotationMatrix).inverse()*G4Translate3D(zpos);
00534             trans = G4Translate3D(fCenterPosition)*trans;
00535           } else {
00536             trans = G4Translate3D(zpos)*G4Translate3D(fCenterPosition);
00537           }
00538           G4double c[4];
00539           colorMap->GetMapColor(rzcell[r][z], c);
00540           att.SetColour(c[0], c[1], c[2]);//, c[3]);
00541           pVisManager->Draw(cylinder, att, trans);
00542         }
00543       }
00544     }
00545   }
00546 
00547   colorMap->SetPSUnit(fDrawUnit);
00548   colorMap->SetPSName(fDrawPSName);
00549   colorMap->DrawColorChart();
00550 
00551 }

void G4ScoringCylinder::GetRZPhi ( G4int  index,
G4int  q[3] 
) const

Definition at line 553 of file G4ScoringCylinder.cc.

Referenced by Draw(), and DrawColumn().

00553                                                               {
00554   // index = k + j * k-size + i * jk-plane-size
00555 
00556   // nested : z -> phi -> r
00557   G4int i = IZ;
00558   G4int j = IPHI;
00559   G4int k = IR;
00560   G4int jk = fNSegment[j]*fNSegment[k];
00561   q[i] = index/jk;
00562   q[j] = (index - q[i]*jk)/fNSegment[k];
00563   q[k] = index - q[j]*fNSegment[k] - q[i]*jk;
00564 }

void G4ScoringCylinder::List (  )  const [virtual]

Reimplemented from G4VScoringMesh.

Definition at line 213 of file G4ScoringCylinder.cc.

References G4VScoringMesh::fSize, G4VScoringMesh::fWorldName, G4cout, and G4VScoringMesh::List().

00213                                    {
00214   G4cout << "G4ScoringCylinder : " << fWorldName << " --- Shape: Cylindrical mesh" << G4endl;
00215 
00216   G4cout << " Size (R, Dz): ("
00217          << fSize[0]/cm << ", "
00218          << fSize[1]/cm << ") [cm]"
00219          << G4endl;
00220 
00221   G4VScoringMesh::List();
00222 }

void G4ScoringCylinder::RegisterPrimitives ( std::vector< G4VPrimitiveScorer * > &  vps  ) 

void G4ScoringCylinder::SetRMax ( G4double  rMax  )  [inline]

Definition at line 54 of file G4ScoringCylinder.hh.

References G4VScoringMesh::fSize.

00054 {fSize[0] = rMax;}

void G4ScoringCylinder::SetZSize ( G4double  zSize  )  [inline]

Definition at line 55 of file G4ScoringCylinder.hh.

References G4VScoringMesh::fSize.

00055 {fSize[1] = zSize;} // half height


The documentation for this class was generated from the following files:
Generated on Mon May 27 17:53:21 2013 for Geant4 by  doxygen 1.4.7