Geant4 Cross Reference |
1 // 2 // ******************************************************************** 3 // * License and Disclaimer * 4 // * * 5 // * The Geant4 software is copyright of the Copyright Holders of * 6 // * the Geant4 Collaboration. It is provided under the terms and * 7 // * conditions of the Geant4 Software License, included in the file * 8 // * LICENSE and available at http://cern.ch/geant4/license . These * 9 // * include a list of copyright holders. * 10 // * * 11 // * Neither the authors of this software system, nor their employing * 12 // * institutes,nor the agencies providing financial support for this * 13 // * work make any representation or warranty, express or implied, * 14 // * regarding this software system or assume any liability for its * 15 // * use. Please see the license in the file LICENSE and URL above * 16 // * for the full disclaimer and the limitation of liability. * 17 // * * 18 // * This code implementation is the result of the scientific and * 19 // * technical work of the GEANT4 collaboration. * 20 // * By using, copying, modifying or distributing the software (or * 21 // * any work based on the software) you agree to acknowledge its * 22 // * use in resulting scientific publications, and indicate your * 23 // * acceptance of all terms of the Geant4 Software license. * 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& wName) 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(G4VPhysicalVolume* fWorldPhys) 64 { 65 if(verboseLevel > 9) 66 G4cout << "G4ScoringBox::SetupGeometry() ..." << G4endl; 67 68 // World 69 G4VPhysicalVolume* scoringWorld = fWorldPhys; 70 G4LogicalVolume* worldLogical = scoringWorld->GetLogicalVolume(); 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] << ", " << fSize[2] << G4endl; 79 G4VSolid* boxSolid = new G4Box(boxName + "0", fSize[0], fSize[1], fSize[2]); 80 auto boxLogical = 81 new G4LogicalVolume(boxSolid, nullptr, boxName + "_0"); 82 new G4PVPlacement(fRotationMatrix, fCenterPosition, boxLogical, boxName + "0", 83 worldLogical, false, 0); 84 85 G4String layerName[2] = { boxName + "_1", boxName + "_2" }; 86 G4VSolid* layerSolid[2]; 87 G4LogicalVolume* layerLogical[2]; 88 89 //-- fisrt nested layer (replicated to x direction) 90 if(verboseLevel > 9) 91 G4cout << "layer 1 :" << G4endl; 92 layerSolid[0] = 93 new G4Box(layerName[0], fSize[0] / fNSegment[0], fSize[1], fSize[2]); 94 layerLogical[0] = new G4LogicalVolume(layerSolid[0], nullptr, layerName[0]); 95 if(fNSegment[0] > 1) 96 { 97 if(verboseLevel > 9) 98 G4cout << "G4ScoringBox::Construct() : Replicate to x direction" 99 << G4endl; 100 if(G4ScoringManager::GetReplicaLevel() > 0) 101 { 102 new G4PVReplica(layerName[0], layerLogical[0], boxLogical, kXAxis, 103 fNSegment[0], fSize[0] / fNSegment[0] * 2.); 104 } 105 else 106 { 107 new G4PVDivision(layerName[0], layerLogical[0], boxLogical, kXAxis, 108 fNSegment[0], 0.); 109 } 110 } 111 else if(fNSegment[0] == 1) 112 { 113 if(verboseLevel > 9) 114 G4cout << "G4ScoringBox::Construct() : Placement" << G4endl; 115 new G4PVPlacement(nullptr, G4ThreeVector(0., 0., 0.), layerLogical[0], 116 layerName[0], boxLogical, false, 0); 117 } 118 else 119 G4cerr << "ERROR : G4ScoringBox::SetupGeometry() : invalid parameter (" 120 << fNSegment[0] << ") " 121 << "in placement of the first nested layer." << G4endl; 122 123 if(verboseLevel > 9) 124 { 125 G4cout << fSize[0] / fNSegment[0] << ", " << fSize[1] << ", " << fSize[2] 126 << G4endl; 127 G4cout << layerName[0] << ": kXAxis, " << fNSegment[0] << ", " 128 << 2. * fSize[0] / fNSegment[0] << G4endl; 129 } 130 131 // second nested layer (replicated to y direction) 132 if(verboseLevel > 9) 133 G4cout << "layer 2 :" << G4endl; 134 layerSolid[1] = new G4Box(layerName[1], fSize[0] / fNSegment[0], 135 fSize[1] / fNSegment[1], fSize[2]); 136 layerLogical[1] = new G4LogicalVolume(layerSolid[1], nullptr, layerName[1]); 137 if(fNSegment[1] > 1) 138 { 139 if(verboseLevel > 9) 140 G4cout << "G4ScoringBox::Construct() : Replicate to y direction" 141 << G4endl; 142 if(G4ScoringManager::GetReplicaLevel() > 1) 143 { 144 new G4PVReplica(layerName[1], layerLogical[1], layerLogical[0], kYAxis, 145 fNSegment[1], fSize[1] / fNSegment[1] * 2.); 146 } 147 else 148 { 149 new G4PVDivision(layerName[1], layerLogical[1], layerLogical[0], kYAxis, 150 fNSegment[1], 0.); 151 } 152 } 153 else if(fNSegment[1] == 1) 154 { 155 if(verboseLevel > 9) 156 G4cout << "G4ScoringBox::Construct() : Placement" << G4endl; 157 new G4PVPlacement(nullptr, G4ThreeVector(0., 0., 0.), layerLogical[1], 158 layerName[1], layerLogical[0], false, 0); 159 } 160 else 161 G4cerr << "ERROR : G4ScoringBox::SetupGeometry() : invalid parameter (" 162 << fNSegment[1] << ") " 163 << "in placement of the second nested layer." << G4endl; 164 165 if(verboseLevel > 9) 166 { 167 G4cout << fSize[0] / fNSegment[0] << ", " << fSize[1] / fNSegment[1] << ", " 168 << fSize[2] << G4endl; 169 G4cout << layerName[1] << ": kYAxis, " << fNSegment[1] << ", " 170 << 2. * fSize[1] / fNSegment[1] << G4endl; 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] / fNSegment[0], fSize[1] / fNSegment[1], 179 fSize[2] / fNSegment[2]); 180 fMeshElementLogical = new G4LogicalVolume(elementSolid, nullptr, elementName); 181 if(fNSegment[2] > 1) 182 { 183 if(verboseLevel > 9) 184 G4cout << "G4ScoringBox::Construct() : Replicate to z direction" 185 << G4endl; 186 187 if(G4ScoringManager::GetReplicaLevel() > 2) 188 { 189 new G4PVReplica(elementName, fMeshElementLogical, layerLogical[1], kZAxis, 190 fNSegment[2], 2. * fSize[2] / fNSegment[2]); 191 } 192 else 193 { 194 new G4PVDivision(elementName, fMeshElementLogical, layerLogical[1], 195 kZAxis, fNSegment[2], 0.); 196 } 197 } 198 else if(fNSegment[2] == 1) 199 { 200 if(verboseLevel > 9) 201 G4cout << "G4ScoringBox::Construct() : Placement" << G4endl; 202 new G4PVPlacement(nullptr, G4ThreeVector(0., 0., 0.), fMeshElementLogical, 203 elementName, layerLogical[1], false, 0); 204 } 205 else 206 G4cerr << "ERROR : G4ScoringBox::SetupGeometry() : " 207 << "invalid parameter (" << fNSegment[2] << ") " 208 << "in mesh element placement." << G4endl; 209 210 if(verboseLevel > 9) 211 { 212 G4cout << fSize[0] / fNSegment[0] << ", " << fSize[1] / fNSegment[1] << ", " 213 << fSize[2] / fNSegment[2] << G4endl; 214 G4cout << elementName << ": kZAxis, " << fNSegment[2] << ", " 215 << 2. * fSize[2] / fNSegment[2] << G4endl; 216 } 217 218 // set the sensitive detector 219 fMeshElementLogical->SetSensitiveDetector(fMFD); 220 221 // vis. attributes 222 auto visatt = new G4VisAttributes(G4Colour(.5, .5, .5)); 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 << " --- Shape: Box mesh" << G4endl; 233 G4cout << " Size (x, y, z): (" << fSize[0] / cm << ", " << fSize[1] / cm 234 << ", " << fSize[2] / cm << ") [cm]" << G4endl; 235 236 G4VScoringMesh::List(); 237 } 238 239 void G4ScoringBox::Draw(RunScore* map, G4VScoreColorMap* colorMap, G4int axflg) 240 { 241 G4VVisManager* pVisManager = G4VVisManager::GetConcreteInstance(); 242 if(pVisManager != nullptr) 243 { 244 // cell vectors 245 std::vector<std::vector<std::vector<G4double>>> cell; // cell[X][Y][Z] 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; // xycell[X][Y] 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; // yzcell[Y][Z] 263 for(G4int y = 0; y < fNSegment[1]; y++) 264 yzcell.push_back(ez); 265 266 std::vector<std::vector<G4double>> xzcell; // xzcell[X][Z] 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_wx()) / fDrawUnitValue; 278 yzcell[q[1]][q[2]] += (itr->second->sum_wx()) / fDrawUnitValue; 279 xzcell[q[0]][q[2]] += (itr->second->sum_wx()) / fDrawUnitValue; 280 } 281 282 // search max. & min. values in each slice 283 G4double xymin = DBL_MAX, yzmin = DBL_MAX, xzmin = 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] / fNSegment[2] - thick); 329 for(G4int x = 0; x < fNSegment[0]; x++) 330 { 331 for(G4int y = 0; y < fNSegment[1]; y++) 332 { 333 G4ThreeVector pos(GetReplicaPosition(x, y, 0) - zhalf); 334 G4ThreeVector pos2(GetReplicaPosition(x, y, fNSegment[2] - 1) + 335 zhalf); 336 G4Transform3D trans, trans2; 337 if(fRotationMatrix != nullptr) 338 { 339 trans = G4Rotate3D(*fRotationMatrix).inverse() * G4Translate3D(pos); 340 trans = G4Translate3D(fCenterPosition) * trans; 341 trans2 = 342 G4Rotate3D(*fRotationMatrix).inverse() * G4Translate3D(pos2); 343 trans2 = G4Translate3D(fCenterPosition) * trans2; 344 } 345 else 346 { 347 trans = G4Translate3D(pos) * G4Translate3D(fCenterPosition); 348 trans2 = G4Translate3D(pos2) * G4Translate3D(fCenterPosition); 349 } 350 G4double c[4]; 351 colorMap->GetMapColor(xycell[x][y], c); 352 att.SetColour(c[0], c[1], c[2]); //, c[3]); 353 354 G4Box xyplate("xy", fSize[0] / fNSegment[0], fSize[1] / fNSegment[1], 355 thick); 356 G4Polyhedron* poly = xyplate.GetPolyhedron(); 357 poly->Transform(trans); 358 poly->SetVisAttributes(&att); 359 pVisManager->Draw(*poly); 360 361 G4Box xyplate2 = xyplate; 362 G4Polyhedron* poly2 = xyplate2.GetPolyhedron(); 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[0] - thick, 0., 0.); 381 for(G4int y = 0; y < fNSegment[1]; y++) 382 { 383 for(G4int z = 0; z < fNSegment[2]; z++) 384 { 385 G4ThreeVector pos(GetReplicaPosition(0, y, z) - xhalf); 386 G4ThreeVector pos2(GetReplicaPosition(fNSegment[0] - 1, y, z) + 387 xhalf); 388 G4Transform3D trans, trans2; 389 if(fRotationMatrix != nullptr) 390 { 391 trans = G4Rotate3D(*fRotationMatrix).inverse() * G4Translate3D(pos); 392 trans = G4Translate3D(fCenterPosition) * trans; 393 trans2 = 394 G4Rotate3D(*fRotationMatrix).inverse() * G4Translate3D(pos2); 395 trans2 = G4Translate3D(fCenterPosition) * trans2; 396 } 397 else 398 { 399 trans = G4Translate3D(pos) * G4Translate3D(fCenterPosition); 400 trans2 = G4Translate3D(pos2) * G4Translate3D(fCenterPosition); 401 } 402 G4double c[4]; 403 colorMap->GetMapColor(yzcell[y][z], c); 404 att.SetColour(c[0], c[1], c[2]); //, c[3]); 405 406 G4Box yzplate("yz", thick, // fSize[0]/fNSegment[0]*0.001, 407 fSize[1] / fNSegment[1], fSize[2] / fNSegment[2]); 408 G4Polyhedron* poly = yzplate.GetPolyhedron(); 409 poly->Transform(trans); 410 poly->SetVisAttributes(&att); 411 pVisManager->Draw(*poly); 412 413 G4Box yzplate2 = yzplate; 414 G4Polyhedron* poly2 = yzplate2.GetPolyhedron(); 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] / fNSegment[1] - thick, 0.); 433 for(G4int x = 0; x < fNSegment[0]; x++) 434 { 435 for(G4int z = 0; z < fNSegment[2]; z++) 436 { 437 G4ThreeVector pos(GetReplicaPosition(x, 0, z) - yhalf); 438 G4ThreeVector pos2(GetReplicaPosition(x, fNSegment[1] - 1, z) + 439 yhalf); 440 G4Transform3D trans, trans2; 441 if(fRotationMatrix != nullptr) 442 { 443 trans = G4Rotate3D(*fRotationMatrix).inverse() * G4Translate3D(pos); 444 trans = G4Translate3D(fCenterPosition) * trans; 445 trans2 = 446 G4Rotate3D(*fRotationMatrix).inverse() * G4Translate3D(pos2); 447 trans2 = G4Translate3D(fCenterPosition) * trans2; 448 } 449 else 450 { 451 trans = G4Translate3D(pos) * G4Translate3D(fCenterPosition); 452 trans2 = G4Translate3D(pos2) * G4Translate3D(fCenterPosition); 453 } 454 G4double c[4]; 455 colorMap->GetMapColor(xzcell[x][z], c); 456 att.SetColour(c[0], c[1], c[2]); //, c[3]); 457 458 G4Box xzplate("xz", fSize[0] / fNSegment[0], 459 thick, // fSize[1]/fNSegment[1]*0.001, 460 fSize[2] / fNSegment[2]); 461 G4Polyhedron* poly = xzplate.GetPolyhedron(); 462 poly->Transform(trans); 463 poly->SetVisAttributes(&att); 464 pVisManager->Draw(*poly); 465 466 G4Box xzplate2 = xzplate; 467 G4Polyhedron* poly2 = xzplate2.GetPolyhedron(); 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(G4int x, G4int y, G4int z) 482 { 483 G4ThreeVector width(fSize[0] / fNSegment[0], fSize[1] / fNSegment[1], 484 fSize[2] / fNSegment[2]); 485 G4ThreeVector pos(-fSize[0] + 2 * (x + 0.5) * width.x(), 486 -fSize[1] + 2 * (y + 0.5) * width.y(), 487 -fSize[2] + 2 * (z + 0.5) * width.z()); 488 489 return pos; 490 } 491 492 void G4ScoringBox::GetXYZ(G4int index, G4int q[3]) const 493 { 494 q[0] = index / (fNSegment[2] * fNSegment[1]); 495 q[1] = (index - q[0] * fNSegment[2] * fNSegment[1]) / fNSegment[2]; 496 q[2] = index - q[1] * fNSegment[2] - q[0] * fNSegment[2] * fNSegment[1]; 497 } 498 499 G4int G4ScoringBox::GetIndex(G4int x, G4int y, G4int z) const 500 { 501 return x + y * fNSegment[0] + z * fNSegment[0] * fNSegment[1]; 502 } 503 504 void G4ScoringBox::DrawColumn(RunScore* map, G4VScoreColorMap* colorMap, 505 G4int idxProj, G4int idxColumn) 506 { 507 G4int iColumn[3] = { 2, 0, 1 }; 508 if(idxColumn < 0 || idxColumn >= fNSegment[iColumn[idxProj]]) 509 { 510 G4cerr << "ERROR : Column number " << idxColumn 511 << " is out of scoring mesh [0," << fNSegment[iColumn[idxProj]] - 1 512 << "]. Method ignored." << G4endl; 513 return; 514 } 515 G4VVisManager* pVisManager = G4VVisManager::GetConcreteInstance(); 516 if(pVisManager != nullptr) 517 { 518 pVisManager->BeginDraw(); 519 520 // cell vectors 521 std::vector<std::vector<std::vector<G4double>>> cell; // cell[X][Y][Z] 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; // xycell[X][Y] 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; // yzcell[Y][Z] 539 for(G4int y = 0; y < fNSegment[1]; y++) 540 yzcell.push_back(ez); 541 542 std::vector<std::vector<G4double>> xzcell; // xzcell[X][Z] 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->sum_wx()) / fDrawUnitValue; 556 } 557 if(idxProj == 1 && q[0] == idxColumn) 558 { // yz plane 559 yzcell[q[1]][q[2]] += (itr->second->sum_wx()) / fDrawUnitValue; 560 } 561 if(idxProj == 2 && q[1] == idxColumn) 562 { // zx plane 563 xzcell[q[0]][q[2]] += (itr->second->sum_wx()) / fDrawUnitValue; 564 } 565 } 566 567 // search max. & min. values in each slice 568 G4double xymin = DBL_MAX, yzmin = DBL_MAX, xzmin = 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] / fNSegment[0], fSize[1] / fNSegment[1], 615 fSize[2] / fNSegment[2]); 616 617 G4ThreeVector pos(GetReplicaPosition(x, y, idxColumn)); 618 G4Transform3D trans; 619 if(fRotationMatrix != nullptr) 620 { 621 trans = G4Rotate3D(*fRotationMatrix).inverse() * G4Translate3D(pos); 622 trans = G4Translate3D(fCenterPosition) * trans; 623 } 624 else 625 { 626 trans = G4Translate3D(pos) * G4Translate3D(fCenterPosition); 627 } 628 G4double c[4]; 629 colorMap->GetMapColor(xycell[x][y], c); 630 att.SetColour(c[0], c[1], c[2]); 631 632 G4Polyhedron* poly = xyplate.GetPolyhedron(); 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] / fNSegment[0], fSize[1] / fNSegment[1], 652 fSize[2] / fNSegment[2]); 653 654 G4ThreeVector pos(GetReplicaPosition(idxColumn, y, z)); 655 G4Transform3D trans; 656 if(fRotationMatrix != nullptr) 657 { 658 trans = G4Rotate3D(*fRotationMatrix).inverse() * G4Translate3D(pos); 659 trans = G4Translate3D(fCenterPosition) * trans; 660 } 661 else 662 { 663 trans = G4Translate3D(pos) * G4Translate3D(fCenterPosition); 664 } 665 G4double c[4]; 666 colorMap->GetMapColor(yzcell[y][z], c); 667 att.SetColour(c[0], c[1], c[2]); //, c[3]); 668 669 G4Polyhedron* poly = yzplate.GetPolyhedron(); 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] / fNSegment[0], fSize[1] / fNSegment[1], 689 fSize[2] / fNSegment[2]); 690 691 G4ThreeVector pos(GetReplicaPosition(x, idxColumn, z)); 692 G4Transform3D trans; 693 if(fRotationMatrix != nullptr) 694 { 695 trans = G4Rotate3D(*fRotationMatrix).inverse() * G4Translate3D(pos); 696 trans = G4Translate3D(fCenterPosition) * trans; 697 } 698 else 699 { 700 trans = G4Translate3D(pos) * G4Translate3D(fCenterPosition); 701 } 702 G4double c[4]; 703 colorMap->GetMapColor(xzcell[x][z], c); 704 att.SetColour(c[0], c[1], c[2]); //, c[3]); 705 706 G4Polyhedron* poly = xzplate.GetPolyhedron(); 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