Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/digits_hits/utils/src/G4ScoringBox.cc

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

Diff markup

Differences between /digits_hits/utils/src/G4ScoringBox.cc (Version 11.3.0) and /digits_hits/utils/src/G4ScoringBox.cc (Version 8.2)


  1 //                                                  1 
  2 // *******************************************    
  3 // * License and Disclaimer                       
  4 // *                                              
  5 // * The  Geant4 software  is  copyright of th    
  6 // * the Geant4 Collaboration.  It is provided    
  7 // * conditions of the Geant4 Software License    
  8 // * LICENSE and available at  http://cern.ch/    
  9 // * include a list of copyright holders.         
 10 // *                                              
 11 // * Neither the authors of this software syst    
 12 // * institutes,nor the agencies providing fin    
 13 // * work  make  any representation or  warran    
 14 // * regarding  this  software system or assum    
 15 // * use.  Please see the license in the file     
 16 // * for the full disclaimer and the limitatio    
 17 // *                                              
 18 // * This  code  implementation is the result     
 19 // * technical work of the GEANT4 collaboratio    
 20 // * By using,  copying,  modifying or  distri    
 21 // * any work based  on the software)  you  ag    
 22 // * use  in  resulting  scientific  publicati    
 23 // * acceptance of all terms of the Geant4 Sof    
 24 // *******************************************    
 25 //                                                
 26 // G4ScoringBox                                   
 27 // -------------------------------------------    
 28                                                   
 29 #include "G4ScoringBox.hh"                        
 30                                                   
 31 #include "G4Box.hh"                               
 32 #include "G4LogicalVolume.hh"                     
 33 #include "G4VPhysicalVolume.hh"                   
 34 #include "G4PVPlacement.hh"                       
 35 #include "G4PVReplica.hh"                         
 36 #include "G4PVDivision.hh"                        
 37 #include "G4VisAttributes.hh"                     
 38 #include "G4VVisManager.hh"                       
 39 #include "G4VScoreColorMap.hh"                    
 40                                                   
 41 #include "G4MultiFunctionalDetector.hh"           
 42 #include "G4VPrimitiveScorer.hh"                  
 43 #include "G4Polyhedron.hh"                        
 44                                                   
 45 #include "G4ScoringManager.hh"                    
 46 #include "G4StatDouble.hh"                        
 47                                                   
 48 #include "G4SystemOfUnits.hh"                     
 49                                                   
 50 #include <map>                                    
 51 #include <fstream>                                
 52                                                   
 53 G4ScoringBox::G4ScoringBox(const G4String& wNa    
 54   : G4VScoringMesh(wName)                         
 55   , fSegmentDirection(-1)                         
 56 {                                                 
 57   fShape                = MeshShape::box;         
 58   fDivisionAxisNames[0] = "X";                    
 59   fDivisionAxisNames[1] = "Y";                    
 60   fDivisionAxisNames[2] = "Z";                    
 61 }                                                 
 62                                                   
 63 void G4ScoringBox::SetupGeometry(G4VPhysicalVo    
 64 {                                                 
 65   if(verboseLevel > 9)                            
 66     G4cout << "G4ScoringBox::SetupGeometry() .    
 67                                                   
 68   // World                                        
 69   G4VPhysicalVolume* scoringWorld = fWorldPhys    
 70   G4LogicalVolume* worldLogical   = scoringWor    
 71                                                   
 72   // Scoring Mesh                                 
 73   if(verboseLevel > 9)                            
 74     G4cout << fWorldName << G4endl;               
 75   G4String boxName = fWorldName;                  
 76                                                   
 77   if(verboseLevel > 9)                            
 78     G4cout << fSize[0] << ", " << fSize[1] <<     
 79   G4VSolid* boxSolid = new G4Box(boxName + "0"    
 80   auto  boxLogical =                              
 81     new G4LogicalVolume(boxSolid, nullptr, box    
 82   new G4PVPlacement(fRotationMatrix, fCenterPo    
 83                     worldLogical, false, 0);      
 84                                                   
 85   G4String layerName[2] = { boxName + "_1", bo    
 86   G4VSolid* layerSolid[2];                        
 87   G4LogicalVolume* layerLogical[2];               
 88                                                   
 89   //-- fisrt nested layer (replicated to x dir    
 90   if(verboseLevel > 9)                            
 91     G4cout << "layer 1 :" << G4endl;              
 92   layerSolid[0] =                                 
 93     new G4Box(layerName[0], fSize[0] / fNSegme    
 94   layerLogical[0] = new G4LogicalVolume(layerS    
 95   if(fNSegment[0] > 1)                            
 96   {                                               
 97     if(verboseLevel > 9)                          
 98       G4cout << "G4ScoringBox::Construct() : R    
 99              << G4endl;                           
100     if(G4ScoringManager::GetReplicaLevel() > 0    
101     {                                             
102       new G4PVReplica(layerName[0], layerLogic    
103                       fNSegment[0], fSize[0] /    
104     }                                             
105     else                                          
106     {                                             
107       new G4PVDivision(layerName[0], layerLogi    
108                        fNSegment[0], 0.);         
109     }                                             
110   }                                               
111   else if(fNSegment[0] == 1)                      
112   {                                               
113     if(verboseLevel > 9)                          
114       G4cout << "G4ScoringBox::Construct() : P    
115     new G4PVPlacement(nullptr, G4ThreeVector(0    
116                       layerName[0], boxLogical    
117   }                                               
118   else                                            
119     G4cerr << "ERROR : G4ScoringBox::SetupGeom    
120            << fNSegment[0] << ") "                
121            << "in placement of the first neste    
122                                                   
123   if(verboseLevel > 9)                            
124   {                                               
125     G4cout << fSize[0] / fNSegment[0] << ", "     
126            << G4endl;                             
127     G4cout << layerName[0] << ": kXAxis, " <<     
128            << 2. * fSize[0] / fNSegment[0] <<     
129   }                                               
130                                                   
131   // second nested layer (replicated to y dire    
132   if(verboseLevel > 9)                            
133     G4cout << "layer 2 :" << G4endl;              
134   layerSolid[1]   = new G4Box(layerName[1], fS    
135                             fSize[1] / fNSegme    
136   layerLogical[1] = new G4LogicalVolume(layerS    
137   if(fNSegment[1] > 1)                            
138   {                                               
139     if(verboseLevel > 9)                          
140       G4cout << "G4ScoringBox::Construct() : R    
141              << G4endl;                           
142     if(G4ScoringManager::GetReplicaLevel() > 1    
143     {                                             
144       new G4PVReplica(layerName[1], layerLogic    
145                       fNSegment[1], fSize[1] /    
146     }                                             
147     else                                          
148     {                                             
149       new G4PVDivision(layerName[1], layerLogi    
150                        fNSegment[1], 0.);         
151     }                                             
152   }                                               
153   else if(fNSegment[1] == 1)                      
154   {                                               
155     if(verboseLevel > 9)                          
156       G4cout << "G4ScoringBox::Construct() : P    
157     new G4PVPlacement(nullptr, G4ThreeVector(0    
158                       layerName[1], layerLogic    
159   }                                               
160   else                                            
161     G4cerr << "ERROR : G4ScoringBox::SetupGeom    
162            << fNSegment[1] << ") "                
163            << "in placement of the second nest    
164                                                   
165   if(verboseLevel > 9)                            
166   {                                               
167     G4cout << fSize[0] / fNSegment[0] << ", "     
168            << fSize[2] << G4endl;                 
169     G4cout << layerName[1] << ": kYAxis, " <<     
170            << 2. * fSize[1] / fNSegment[1] <<     
171   }                                               
172                                                   
173   // mesh elements (replicated to z direction)    
174   if(verboseLevel > 9)                            
175     G4cout << "mesh elements :" << G4endl;        
176   G4String elementName = boxName + "_3";          
177   G4VSolid* elementSolid =                        
178     new G4Box(elementName, fSize[0] / fNSegmen    
179               fSize[2] / fNSegment[2]);           
180   fMeshElementLogical = new G4LogicalVolume(el    
181   if(fNSegment[2] > 1)                            
182   {                                               
183     if(verboseLevel > 9)                          
184       G4cout << "G4ScoringBox::Construct() : R    
185              << G4endl;                           
186                                                   
187     if(G4ScoringManager::GetReplicaLevel() > 2    
188     {                                             
189       new G4PVReplica(elementName, fMeshElemen    
190                       fNSegment[2], 2. * fSize    
191     }                                             
192     else                                          
193     {                                             
194       new G4PVDivision(elementName, fMeshEleme    
195                        kZAxis, fNSegment[2], 0    
196     }                                             
197   }                                               
198   else if(fNSegment[2] == 1)                      
199   {                                               
200     if(verboseLevel > 9)                          
201       G4cout << "G4ScoringBox::Construct() : P    
202     new G4PVPlacement(nullptr, G4ThreeVector(0    
203                       elementName, layerLogica    
204   }                                               
205   else                                            
206     G4cerr << "ERROR : G4ScoringBox::SetupGeom    
207            << "invalid parameter (" << fNSegme    
208            << "in mesh element placement." <<     
209                                                   
210   if(verboseLevel > 9)                            
211   {                                               
212     G4cout << fSize[0] / fNSegment[0] << ", "     
213            << fSize[2] / fNSegment[2] << G4end    
214     G4cout << elementName << ": kZAxis, " << f    
215            << 2. * fSize[2] / fNSegment[2] <<     
216   }                                               
217                                                   
218   // set the sensitive detector                   
219   fMeshElementLogical->SetSensitiveDetector(fM    
220                                                   
221   // vis. attributes                              
222   auto  visatt = new G4VisAttributes(G4Colour(    
223   visatt->SetVisibility(false);                   
224   layerLogical[0]->SetVisAttributes(visatt);      
225   layerLogical[1]->SetVisAttributes(visatt);      
226   visatt->SetVisibility(true);                    
227   fMeshElementLogical->SetVisAttributes(visatt    
228 }                                                 
229                                                   
230 void G4ScoringBox::List() const                   
231 {                                                 
232   G4cout << "G4ScoringBox : " << fWorldName <<    
233   G4cout << " Size (x, y, z): (" << fSize[0] /    
234          << ", " << fSize[2] / cm << ") [cm]"     
235                                                   
236   G4VScoringMesh::List();                         
237 }                                                 
238                                                   
239 void G4ScoringBox::Draw(RunScore* map, G4VScor    
240 {                                                 
241   G4VVisManager* pVisManager = G4VVisManager::    
242   if(pVisManager != nullptr)                      
243   {                                               
244     // cell vectors                               
245     std::vector<std::vector<std::vector<G4doub    
246     std::vector<G4double> ez;                     
247     for(G4int z = 0; z < fNSegment[2]; z++)       
248       ez.push_back(0.);                           
249     std::vector<std::vector<G4double>> eyz;       
250     for(G4int y = 0; y < fNSegment[1]; y++)       
251       eyz.push_back(ez);                          
252     for(G4int x = 0; x < fNSegment[0]; x++)       
253       cell.push_back(eyz);                        
254                                                   
255     std::vector<std::vector<G4double>> xycell;    
256     std::vector<G4double> ey;                     
257     for(G4int y = 0; y < fNSegment[1]; y++)       
258       ey.push_back(0.);                           
259     for(G4int x = 0; x < fNSegment[0]; x++)       
260       xycell.push_back(ey);                       
261                                                   
262     std::vector<std::vector<G4double>> yzcell;    
263     for(G4int y = 0; y < fNSegment[1]; y++)       
264       yzcell.push_back(ez);                       
265                                                   
266     std::vector<std::vector<G4double>> xzcell;    
267     for(G4int x = 0; x < fNSegment[0]; x++)       
268       xzcell.push_back(ez);                       
269                                                   
270     // projections                                
271     G4int q[3];                                   
272     auto itr = map->GetMap()->begin();            
273     for(; itr != map->GetMap()->end(); itr++)     
274     {                                             
275       GetXYZ(itr->first, q);                      
276                                                   
277       xycell[q[0]][q[1]] += (itr->second->sum_    
278       yzcell[q[1]][q[2]] += (itr->second->sum_    
279       xzcell[q[0]][q[2]] += (itr->second->sum_    
280     }                                             
281                                                   
282     // search max. & min. values in each slice    
283     G4double xymin = DBL_MAX, yzmin = DBL_MAX,    
284     G4double xymax = 0., yzmax = 0., xzmax = 0    
285     for(G4int x = 0; x < fNSegment[0]; x++)       
286     {                                             
287       for(G4int y = 0; y < fNSegment[1]; y++)     
288       {                                           
289         if(xymin > xycell[x][y])                  
290           xymin = xycell[x][y];                   
291         if(xymax < xycell[x][y])                  
292           xymax = xycell[x][y];                   
293       }                                           
294       for(G4int z = 0; z < fNSegment[2]; z++)     
295       {                                           
296         if(xzmin > xzcell[x][z])                  
297           xzmin = xzcell[x][z];                   
298         if(xzmax < xzcell[x][z])                  
299           xzmax = xzcell[x][z];                   
300       }                                           
301     }                                             
302     for(G4int y = 0; y < fNSegment[1]; y++)       
303     {                                             
304       for(G4int z = 0; z < fNSegment[2]; z++)     
305       {                                           
306         if(yzmin > yzcell[y][z])                  
307           yzmin = yzcell[y][z];                   
308         if(yzmax < yzcell[y][z])                  
309           yzmax = yzcell[y][z];                   
310       }                                           
311     }                                             
312                                                   
313     G4VisAttributes att;                          
314     att.SetForceSolid(true);                      
315     att.SetForceAuxEdgeVisible(true);             
316     G4double thick = 0.01;                        
317                                                   
318     G4Scale3D scale;                              
319     if(axflg / 100 == 1)                          
320     {                                             
321       pVisManager->BeginDraw();                   
322                                                   
323       // xy plane                                 
324       if(colorMap->IfFloatMinMax())               
325       {                                           
326         colorMap->SetMinMax(xymin, xymax);        
327       }                                           
328       G4ThreeVector zhalf(0., 0., fSize[2] / f    
329       for(G4int x = 0; x < fNSegment[0]; x++)     
330       {                                           
331         for(G4int y = 0; y < fNSegment[1]; y++    
332         {                                         
333           G4ThreeVector pos(GetReplicaPosition    
334           G4ThreeVector pos2(GetReplicaPositio    
335                              zhalf);              
336           G4Transform3D trans, trans2;            
337           if(fRotationMatrix != nullptr)          
338           {                                       
339             trans = G4Rotate3D(*fRotationMatri    
340             trans = G4Translate3D(fCenterPosit    
341             trans2 =                              
342               G4Rotate3D(*fRotationMatrix).inv    
343             trans2 = G4Translate3D(fCenterPosi    
344           }                                       
345           else                                    
346           {                                       
347             trans  = G4Translate3D(pos) * G4Tr    
348             trans2 = G4Translate3D(pos2) * G4T    
349           }                                       
350           G4double c[4];                          
351           colorMap->GetMapColor(xycell[x][y],     
352           att.SetColour(c[0], c[1], c[2]);  //    
353                                                   
354           G4Box xyplate("xy", fSize[0] / fNSeg    
355                         thick);                   
356           G4Polyhedron* poly = xyplate.GetPoly    
357           poly->Transform(trans);                 
358           poly->SetVisAttributes(&att);           
359           pVisManager->Draw(*poly);               
360                                                   
361           G4Box xyplate2      = xyplate;          
362           G4Polyhedron* poly2 = xyplate2.GetPo    
363           poly2->Transform(trans2);               
364           poly2->SetVisAttributes(&att);          
365           pVisManager->Draw(*poly2);              
366         }                                         
367       }                                           
368       pVisManager->EndDraw();                     
369     }                                             
370     axflg = axflg % 100;                          
371     if(axflg / 10 == 1)                           
372     {                                             
373       pVisManager->BeginDraw();                   
374                                                   
375       // yz plane                                 
376       if(colorMap->IfFloatMinMax())               
377       {                                           
378         colorMap->SetMinMax(yzmin, yzmax);        
379       }                                           
380       G4ThreeVector xhalf(fSize[0] / fNSegment    
381       for(G4int y = 0; y < fNSegment[1]; y++)     
382       {                                           
383         for(G4int z = 0; z < fNSegment[2]; z++    
384         {                                         
385           G4ThreeVector pos(GetReplicaPosition    
386           G4ThreeVector pos2(GetReplicaPositio    
387                              xhalf);              
388           G4Transform3D trans, trans2;            
389           if(fRotationMatrix != nullptr)          
390           {                                       
391             trans = G4Rotate3D(*fRotationMatri    
392             trans = G4Translate3D(fCenterPosit    
393             trans2 =                              
394               G4Rotate3D(*fRotationMatrix).inv    
395             trans2 = G4Translate3D(fCenterPosi    
396           }                                       
397           else                                    
398           {                                       
399             trans  = G4Translate3D(pos) * G4Tr    
400             trans2 = G4Translate3D(pos2) * G4T    
401           }                                       
402           G4double c[4];                          
403           colorMap->GetMapColor(yzcell[y][z],     
404           att.SetColour(c[0], c[1], c[2]);  //    
405                                                   
406           G4Box yzplate("yz", thick,  // fSize    
407                         fSize[1] / fNSegment[1    
408           G4Polyhedron* poly = yzplate.GetPoly    
409           poly->Transform(trans);                 
410           poly->SetVisAttributes(&att);           
411           pVisManager->Draw(*poly);               
412                                                   
413           G4Box yzplate2      = yzplate;          
414           G4Polyhedron* poly2 = yzplate2.GetPo    
415           poly2->Transform(trans2);               
416           poly2->SetVisAttributes(&att);          
417           pVisManager->Draw(*poly2);              
418         }                                         
419       }                                           
420       pVisManager->EndDraw();                     
421     }                                             
422     axflg = axflg % 10;                           
423     if(axflg == 1)                                
424     {                                             
425       pVisManager->BeginDraw();                   
426                                                   
427       // xz plane                                 
428       if(colorMap->IfFloatMinMax())               
429       {                                           
430         colorMap->SetMinMax(xzmin, xzmax);        
431       }                                           
432       G4ThreeVector yhalf(0., fSize[1] / fNSeg    
433       for(G4int x = 0; x < fNSegment[0]; x++)     
434       {                                           
435         for(G4int z = 0; z < fNSegment[2]; z++    
436         {                                         
437           G4ThreeVector pos(GetReplicaPosition    
438           G4ThreeVector pos2(GetReplicaPositio    
439                              yhalf);              
440           G4Transform3D trans, trans2;            
441           if(fRotationMatrix != nullptr)          
442           {                                       
443             trans = G4Rotate3D(*fRotationMatri    
444             trans = G4Translate3D(fCenterPosit    
445             trans2 =                              
446               G4Rotate3D(*fRotationMatrix).inv    
447             trans2 = G4Translate3D(fCenterPosi    
448           }                                       
449           else                                    
450           {                                       
451             trans  = G4Translate3D(pos) * G4Tr    
452             trans2 = G4Translate3D(pos2) * G4T    
453           }                                       
454           G4double c[4];                          
455           colorMap->GetMapColor(xzcell[x][z],     
456           att.SetColour(c[0], c[1], c[2]);  //    
457                                                   
458           G4Box xzplate("xz", fSize[0] / fNSeg    
459                         thick,  // fSize[1]/fN    
460                         fSize[2] / fNSegment[2    
461           G4Polyhedron* poly = xzplate.GetPoly    
462           poly->Transform(trans);                 
463           poly->SetVisAttributes(&att);           
464           pVisManager->Draw(*poly);               
465                                                   
466           G4Box xzplate2      = xzplate;          
467           G4Polyhedron* poly2 = xzplate2.GetPo    
468           poly2->Transform(trans2);               
469           poly2->SetVisAttributes(&att);          
470           pVisManager->Draw(*poly2);              
471         }                                         
472       }                                           
473       pVisManager->EndDraw();                     
474     }                                             
475   }                                               
476   colorMap->SetPSUnit(fDrawUnit);                 
477   colorMap->SetPSName(fDrawPSName);               
478   colorMap->DrawColorChart();                     
479 }                                                 
480                                                   
481 G4ThreeVector G4ScoringBox::GetReplicaPosition    
482 {                                                 
483   G4ThreeVector width(fSize[0] / fNSegment[0],    
484                       fSize[2] / fNSegment[2])    
485   G4ThreeVector pos(-fSize[0] + 2 * (x + 0.5)     
486                     -fSize[1] + 2 * (y + 0.5)     
487                     -fSize[2] + 2 * (z + 0.5)     
488                                                   
489   return pos;                                     
490 }                                                 
491                                                   
492 void G4ScoringBox::GetXYZ(G4int index, G4int q    
493 {                                                 
494   q[0] = index / (fNSegment[2] * fNSegment[1])    
495   q[1] = (index - q[0] * fNSegment[2] * fNSegm    
496   q[2] = index - q[1] * fNSegment[2] - q[0] *     
497 }                                                 
498                                                   
499 G4int G4ScoringBox::GetIndex(G4int x, G4int y,    
500 {                                                 
501   return x + y * fNSegment[0] + z * fNSegment[    
502 }                                                 
503                                                   
504 void G4ScoringBox::DrawColumn(RunScore* map, G    
505                               G4int idxProj, G    
506 {                                                 
507   G4int iColumn[3] = { 2, 0, 1 };                 
508   if(idxColumn < 0 || idxColumn >= fNSegment[i    
509   {                                               
510     G4cerr << "ERROR : Column number " << idxC    
511            << " is out of scoring mesh [0," <<    
512            << "]. Method ignored." << G4endl;     
513     return;                                       
514   }                                               
515   G4VVisManager* pVisManager = G4VVisManager::    
516   if(pVisManager != nullptr)                      
517   {                                               
518     pVisManager->BeginDraw();                     
519                                                   
520     // cell vectors                               
521     std::vector<std::vector<std::vector<G4doub    
522     std::vector<G4double> ez;                     
523     for(G4int z = 0; z < fNSegment[2]; z++)       
524       ez.push_back(0.);                           
525     std::vector<std::vector<G4double>> eyz;       
526     for(G4int y = 0; y < fNSegment[1]; y++)       
527       eyz.push_back(ez);                          
528     for(G4int x = 0; x < fNSegment[0]; x++)       
529       cell.push_back(eyz);                        
530                                                   
531     std::vector<std::vector<G4double>> xycell;    
532     std::vector<G4double> ey;                     
533     for(G4int y = 0; y < fNSegment[1]; y++)       
534       ey.push_back(0.);                           
535     for(G4int x = 0; x < fNSegment[0]; x++)       
536       xycell.push_back(ey);                       
537                                                   
538     std::vector<std::vector<G4double>> yzcell;    
539     for(G4int y = 0; y < fNSegment[1]; y++)       
540       yzcell.push_back(ez);                       
541                                                   
542     std::vector<std::vector<G4double>> xzcell;    
543     for(G4int x = 0; x < fNSegment[0]; x++)       
544       xzcell.push_back(ez);                       
545                                                   
546     // projections                                
547     G4int q[3];                                   
548     auto itr = map->GetMap()->begin();            
549     for(; itr != map->GetMap()->end(); itr++)     
550     {                                             
551       GetXYZ(itr->first, q);                      
552                                                   
553       if(idxProj == 0 && q[2] == idxColumn)       
554       {  // xy plane                              
555         xycell[q[0]][q[1]] += (itr->second->su    
556       }                                           
557       if(idxProj == 1 && q[0] == idxColumn)       
558       {  // yz plane                              
559         yzcell[q[1]][q[2]] += (itr->second->su    
560       }                                           
561       if(idxProj == 2 && q[1] == idxColumn)       
562       {  // zx plane                              
563         xzcell[q[0]][q[2]] += (itr->second->su    
564       }                                           
565     }                                             
566                                                   
567     // search max. & min. values in each slice    
568     G4double xymin = DBL_MAX, yzmin = DBL_MAX,    
569     G4double xymax = 0., yzmax = 0., xzmax = 0    
570     for(G4int x = 0; x < fNSegment[0]; x++)       
571     {                                             
572       for(G4int y = 0; y < fNSegment[1]; y++)     
573       {                                           
574         if(xymin > xycell[x][y])                  
575           xymin = xycell[x][y];                   
576         if(xymax < xycell[x][y])                  
577           xymax = xycell[x][y];                   
578       }                                           
579       for(G4int z = 0; z < fNSegment[2]; z++)     
580       {                                           
581         if(xzmin > xzcell[x][z])                  
582           xzmin = xzcell[x][z];                   
583         if(xzmax < xzcell[x][z])                  
584           xzmax = xzcell[x][z];                   
585       }                                           
586     }                                             
587     for(G4int y = 0; y < fNSegment[1]; y++)       
588     {                                             
589       for(G4int z = 0; z < fNSegment[2]; z++)     
590       {                                           
591         if(yzmin > yzcell[y][z])                  
592           yzmin = yzcell[y][z];                   
593         if(yzmax < yzcell[y][z])                  
594           yzmax = yzcell[y][z];                   
595       }                                           
596     }                                             
597                                                   
598     G4VisAttributes att;                          
599     att.SetForceSolid(true);                      
600     att.SetForceAuxEdgeVisible(true);             
601                                                   
602     G4Scale3D scale;                              
603     // xy plane                                   
604     if(idxProj == 0)                              
605     {                                             
606       if(colorMap->IfFloatMinMax())               
607       {                                           
608         colorMap->SetMinMax(xymin, xymax);        
609       }                                           
610       for(G4int x = 0; x < fNSegment[0]; x++)     
611       {                                           
612         for(G4int y = 0; y < fNSegment[1]; y++    
613         {                                         
614           G4Box xyplate("xy", fSize[0] / fNSeg    
615                         fSize[2] / fNSegment[2    
616                                                   
617           G4ThreeVector pos(GetReplicaPosition    
618           G4Transform3D trans;                    
619           if(fRotationMatrix != nullptr)          
620           {                                       
621             trans = G4Rotate3D(*fRotationMatri    
622             trans = G4Translate3D(fCenterPosit    
623           }                                       
624           else                                    
625           {                                       
626             trans = G4Translate3D(pos) * G4Tra    
627           }                                       
628           G4double c[4];                          
629           colorMap->GetMapColor(xycell[x][y],     
630           att.SetColour(c[0], c[1], c[2]);        
631                                                   
632           G4Polyhedron* poly = xyplate.GetPoly    
633           poly->Transform(trans);                 
634           poly->SetVisAttributes(att);            
635           pVisManager->Draw(*poly);               
636         }                                         
637       }                                           
638     }                                             
639     else                                          
640       // yz plane                                 
641       if(idxProj == 1)                            
642     {                                             
643       if(colorMap->IfFloatMinMax())               
644       {                                           
645         colorMap->SetMinMax(yzmin, yzmax);        
646       }                                           
647       for(G4int y = 0; y < fNSegment[1]; y++)     
648       {                                           
649         for(G4int z = 0; z < fNSegment[2]; z++    
650         {                                         
651           G4Box yzplate("yz", fSize[0] / fNSeg    
652                         fSize[2] / fNSegment[2    
653                                                   
654           G4ThreeVector pos(GetReplicaPosition    
655           G4Transform3D trans;                    
656           if(fRotationMatrix != nullptr)          
657           {                                       
658             trans = G4Rotate3D(*fRotationMatri    
659             trans = G4Translate3D(fCenterPosit    
660           }                                       
661           else                                    
662           {                                       
663             trans = G4Translate3D(pos) * G4Tra    
664           }                                       
665           G4double c[4];                          
666           colorMap->GetMapColor(yzcell[y][z],     
667           att.SetColour(c[0], c[1], c[2]);  //    
668                                                   
669           G4Polyhedron* poly = yzplate.GetPoly    
670           poly->Transform(trans);                 
671           poly->SetVisAttributes(att);            
672           pVisManager->Draw(*poly);               
673         }                                         
674       }                                           
675     }                                             
676     else                                          
677       // xz plane                                 
678       if(idxProj == 2)                            
679     {                                             
680       if(colorMap->IfFloatMinMax())               
681       {                                           
682         colorMap->SetMinMax(xzmin, xzmax);        
683       }                                           
684       for(G4int x = 0; x < fNSegment[0]; x++)     
685       {                                           
686         for(G4int z = 0; z < fNSegment[2]; z++    
687         {                                         
688           G4Box xzplate("xz", fSize[0] / fNSeg    
689                         fSize[2] / fNSegment[2    
690                                                   
691           G4ThreeVector pos(GetReplicaPosition    
692           G4Transform3D trans;                    
693           if(fRotationMatrix != nullptr)          
694           {                                       
695             trans = G4Rotate3D(*fRotationMatri    
696             trans = G4Translate3D(fCenterPosit    
697           }                                       
698           else                                    
699           {                                       
700             trans = G4Translate3D(pos) * G4Tra    
701           }                                       
702           G4double c[4];                          
703           colorMap->GetMapColor(xzcell[x][z],     
704           att.SetColour(c[0], c[1], c[2]);  //    
705                                                   
706           G4Polyhedron* poly = xzplate.GetPoly    
707           poly->Transform(trans);                 
708           poly->SetVisAttributes(att);            
709           pVisManager->Draw(*poly);               
710         }                                         
711       }                                           
712     }                                             
713     pVisManager->EndDraw();                       
714   }                                               
715                                                   
716   colorMap->SetPSUnit(fDrawUnit);                 
717   colorMap->SetPSName(fDrawPSName);               
718   colorMap->DrawColorChart();                     
719 }                                                 
720